Project 1 - Parsing NASA thermodynamic data
AIM : Write a function that extracts the 14 co-efficients and calculates the enthalpy, entropy and specific heats for all the species in the data file.
NASA Thermodynamic data :
NASA came up with polynomials that can be used to evaluate thermodynamic properties such as Cp, H and S. They have also documented the co-efficients that are required to evaluate these polynomials for 1000+ species.
The formulae are shown below.
Here R is the universal gas constant, T is the temeprature
Specific heat:
Specific heat, ratio of the quantity of heat required to raise the temperature of a body one degree to that required to raise the temperature of an equal mass of water one degree.
Function to find specific heat :
for calculating specific heat 'fuction' is used and by declaring inputs paramters outside the function.
And assigng the related formula from 14 coefficient.
Now create a 'if loop' condition and make the temperature less than equal to temperature2.
Then the calculation for high temperature and low temperature co-efficient values are done respectively.
function cp = specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,temp2)
%comparing with the temperature
%calculating specific heat for higher and lower temperature respectively.
if(T<= temp2)
cp = R*(a8+(a9*T)+(a10*(T).^2)+(a11*(T).^3)+(a12*(T).^4));
else
cp = R*(a1+(a2*T)+(a3*(T).^2)+(a4*(T).^3)+(a5*(T).^4));
end
end
Enthalpy:
Enthalpy is defined as follows H = E + PV. Enthalpy is the sum of the internal energy. and pressure times volume. We cannot measure the enthalpy of a system, but we can look at changes in enthalpy.
Function to find Enthalpy: :
for calculating Enthalpy 'fuction' is used and by declaring inputs paramters outside the function.
And assigng the related formula from 14 coefficient.
Now create a 'if loop' condition and make the temperature less than equal to temperature2.
Then the calculation for high temperature and low temperature co-efficient values are done respectively.
function H = enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,temp2)
%comparing with the temperature
%calculating enthalpy for higher and lower temperature respectively.
if(T<= temp2)
H = R*T.*(a8+((a9*T)/2)+((a10*(T).^2)/3)+((a11*(T).^3)/4)+((a12*(T).^4)/5)+((a13./T)));
else
H = R*T.*(a1+((a2*T)/2)+((a3*(T).^2)/3)+((a4*(T).^3)/4)+((a5*(T).^4)/5)+((a6./T)));
end
endENTROPY :
IN thermodynamics : a measure of the unavailable energy in a closed thermodynamic system that is also usually considered to be a measure of the system's disorder, that is a property of the system's state, and that varies directly with any reversible change in heat in the system and inversely with the temperature of the system. The degree of disorder or uncertainty in a system.
Function to find ENTROPY: :
for calculating Entropy 'fuction' is used and by declaring inputs paramters outside the function.
And assigng the related formula from 14 coefficient.
Now create a 'if loop' condition and make the temperature less than equal to temperature2.
Then the calculation for high temperature and low temperature co-efficient values are done respectively.
function S = entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,temp2)
%comparing with the temperature
%calculating enthalpy for higher and lower temperature respectively.
if(T<= temp2)
S = R*(a8.*log(T)+(a9*T)+((a10*(T).^2)/2)+((a11*(T).^3)/3)+((a12*(T).^4)/4)+(a14));
else
S = R*(a1.*log(T)+(a2*T)+((a3*(T).^2)/2)+((a4*(T).^3)/3)+((a5*(T).^4)/4)+(a7));
end
endFunction to find MOLECULAR WEIGHT of different species:
To create Molecular_weight function W and assiged the species input into that.Define the main element of all 53 species in an array of atomic_name.atomic_weight = [1,12,14,16,40] from the sources. Then assign the molecular range w from 0 and then create an for loop codition.Create an 'if' condition and 'strcmp' to compare the 2 strong.now the string to the number using str2num command is used and assigned as n.now we will find if there is more elements in the species then increments atomic weight of the species if n > 1
W = W+atomic_weight(position)*(n-1).
function W = Molecular_weight(species)
%defining atomic name of the species
atomic_name = ['H','C','O','N','Ar'];
atomic_weight = [1,12,14,16,40];
W=0; %molecular weight dtartinf range
for i = 1:length(species) %loop runs to th length of species
for j = 1:length(atomic_name) %loop runs to the length of atomic names
if strcmp(species(i), atomic_name(j)) %function to compare 2 string
W = W+atomic_weight(j);
position = j;
end
end
n =str2num(species(i)); %now we will find if there is more elements in the species then increments atomic weight of the species
if n > 1
W = W+atomic_weight(position)*(n-1);
end
end
fprintf ('molecular weight of %s',species)
fprintf('%f',W);
disp('.')
end
MAIN CODE:
Step 1: NASA thermodynamic data open the file and read the data.(fopen,fgetl) command is used.
f1 = fopen('THERMO.dat','r');
get = fgetl(f1);
Open the THERMO.dat file for reading by using the fopen command and assinging to the variable as f1.reading the temperature from this file.
temp = fgetl(f1). after reading the temperature from the file we have splitted the string using strsplit command.It assigned the varriable as (temp).
splitting string from the line
s = strsplit(temp).
Then the string value will converted into number by using command 'str2double'.converting string cells into numbers
temp1 = str2double(s(2));
temp2 = str2double(s(3));
temp3 = str2double(s(4));
setting the temperature range using linspace command.T = linspace(temp1,temp3,1000). starts from temp1 and end at temp3 which will have 1000 values.
AND Then skipping the line using 'for loop' command .skipping commment lines
for i = 1:3
get = fgetl(f1).
step 2 : EXTRACTING the all species data.
using 'for loop' command. for j = 1:53
tline = fgetl(f1);
and splitting into cell array which is stored in variable A
A = strsplit(tline,' ' ).
To read the name of species.species = (A{1}). and next converting string to number by using command 'str2double'.converting string cells into numbers.
lowtemp = str2double(A{end-3});
midtemp = str2double(A{end-2});
hightemp = str2double(A{end-1});
And to getting the next line using command'fgetl'.tline1 = fgetl(f1);
To find thae string E
a = strfind(tline1,'E');
step 3 : finding first 7 higher and second 7 lower co-efficients.
finding higher co-efficient.to specify the specific heat(cp),enthalpy(H) and entropy(S) function to find these parameters using these 14 co-efficient for high and higher temperrature.And it will soecify the molecular weight W function to find the molecular wight of different species.calculate specific heat
cp = specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,temp2);
calculate enthalpy
H = enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,temp2);
calculate entropy
S = entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,temp2);
calculate molecular weight of species
W = Molecular_weight(species);
writing molecular weight file
MW_file = fopen('Molecular weight.txt','w');
fprintf(MW_file,'Molecular weight of %s is %d',species,W);
fclose(MW_file); %closing the file
step4 : ploting
ploting cp,h,s with temperature. now changing the directory using 'cd'.saving in a folder cd('C:\Users\USER\Desktop\temp pogram matlab\New folder\species') . and then creating directory the folder where as to save our files using 'mkdir'.mkdir(species). then plot specific heat,enthalpy,entropy with temperature of different species. and gives perfect title to each. then saveas command is used to save the plot of different species in the directory folder.
plot 1 cp VS t
figure(1)
plot(T,cp,'linewidth',2,'color','b');
xlabel('Temperature(K)');
ylabel('Specific heat(KJ)');
title(sprintf('specific heat VS temperature of %s',species));
grid on
saveas(figure(1),'Specific heat.jpg');
plot 2 H VS t
figure(2)
plot(T,H,'linewidth',2,'color','r');
xlabel('Temperature(K)');
ylabel('Enthalpy (KJ/mol)');
title(sprintf('Enthalpy VS temperature of %s',species));
grid on
saveas(figure(2),'Enthalpy.jpg');
plot 3 s VS t
figure(3)
plot(T,S,'linewidth',2,'color','g');
xlabel('Temperature(K)');
ylabel('Entropy (KJ/mol-K)');
title(sprintf('Entropy VS temperature of %s',species));
grid on
saveas(figure(3),'Entropy.jpg');
cd('C:\Users\USER\Desktop\temp pogram matlab\New folder') %changing the directory to original matlab path
%molecular weight
W = Molecular_weight(species).
CODE :
clear all
close all
clc
%input
R = 8.314; %gas constant
%reading data file
%open data file and reading THERMO.data file
f1 = fopen('THERMO.dat','r');
get = fgetl(f1);
%reading the temperature from this file
temp = fgetl(f1);
%splitting string from the line
s = strsplit(temp);
%converting string cells into numbers
temp1 = str2double(s(2));
temp2 = str2double(s(3));
temp3 = str2double(s(4));
%temparature range
T = linspace(temp1,temp3,1000);
%skipping commment lines
for i = 1:3
get = fgetl(f1);
end
%reading the species
for j = 1:53
tline = fgetl(f1);
A = strsplit(tline,' ' );%splitting into cell array
species = (A{1}); %to read the name of species
%converting string to number
lowtemp = str2double(A{end-3});
midtemp = str2double(A{end-2});
hightemp = str2double(A{end-1});
tline1 = fgetl(f1);% getting next line
a = strfind(tline1,'E'); %finding the string E
% finding first 7 higher and second 7 lower co-efficients
%finding higher co-efficient
a1 = tline1(1:a(1)+3);
a1 = str2double(a1); %1
a2 = tline1(a(1)+4:a(2)+3);
a2 = str2double(a2); %2
a3 = tline1(a(2)+4:a(3)+3);
a3 = str2double(a3); %3
a4 = tline1(a(3)+4:a(4)+3);
a4 = str2double(a4); %4
a5 = tline1(a(4)+4:a(5)+3);
a5 = str2double(a5); %5
tline2 = fgetl(f1);
b = strfind(tline2,'E');
a6 = tline2(1:b(1)+3);
a6 = str2double(a6); %6
a7 = tline2(b(1)+4:b(2)+3);
a7 = str2double(a7); %7
%finding lower co-efficient
a8 = tline2(b(2)+4:b(3)+3);
a8 = str2double(a8); %8
a9 = tline2(b(3)+4:b(4)+3);
a9 = str2double(a9); %9
a10 = tline2(b(4)+4:b(5)+3);
a10 = str2double(a10); %10
tline3 = fgetl(f1);
c = strfind(tline3,'E');
a11 = tline3(1:c(1)+3);
a11 = str2double(a11); %11
a12 = tline3(c(1)+4:c(2)+3);
a12 = str2double(a12); %12
a13 = tline3(c(2)+4:c(3)+3);
a13 = str2double(a13); %13
a14 = tline3(c(3)+4:c(4)+3);
a14 = str2double(a14); %14
%calculate specific heat
cp = specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,temp2);
%calculate enthalpy
H = enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,temp2);
%calculate entropy
S = entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,temp2);
%calculate molecular weight of species
W = Molecular_weight(species);
%writing molecular weight file
MW_file = fopen('Molecular weight.txt','w');
fprintf(MW_file,'Molecular weight of %s is %d',species,W);
fclose(MW_file); %closing the file
%ploting cp,h,s with temperature
%saving in afolder
cd('C:\Users\USER\Desktop\temp pogram matlab\New folder\species') %changing directory
mkdir(species); %creating directory
cd(species);
%plot 1 cp VS t
figure(1)
plot(T,cp,'linewidth',2,'color','b');
xlabel('Temperature(K)');
ylabel('Specific heat(KJ)');
title(sprintf('specific heat VS temperature of %s',species));
grid on
saveas(figure(1),'Specific heat.jpg');
%plot 2 H VS t
figure(2)
plot(T,H,'linewidth',2,'color','r');
xlabel('Temperature(K)');
ylabel('Enthalpy (KJ/mol)');
title(sprintf('Enthalpy VS temperature of %s',species));
grid on
saveas(figure(2),'Enthalpy.jpg');
%plot 3 s VS t
figure(3)
plot(T,S,'linewidth',2,'color','g');
xlabel('Temperature(K)');
ylabel('Entropy (KJ/mol-K)');
title(sprintf('Entropy VS temperature of %s',species));
grid on
saveas(figure(3),'Entropy.jpg');
cd('C:\Users\USER\Desktop\temp pogram matlab\New folder') %changing the directory to original matlab path
%molecular weight
W = Molecular_weight(species)
end
===============================x===================================x=======================
molecular weight of O14.000000.
molecular weight of O14.000000.
W =
14
molecular weight of O228.000000.
molecular weight of O228.000000.
W =
28
molecular weight of H1.000000.
molecular weight of H1.000000.
W =
1
molecular weight of H22.000000.
molecular weight of H22.000000.
W =
2
molecular weight of OH15.000000.
molecular weight of OH15.000000.
W =
15
molecular weight of H2O16.000000.
molecular weight of H2O16.000000.
W =
16
molecular weight of HO229.000000.
molecular weight of HO229.000000.
W =
29
molecular weight of H2O230.000000.
molecular weight of H2O230.000000.
W =
30
molecular weight of C12.000000.
molecular weight of C12.000000.
W =
12
molecular weight of CH13.000000.
molecular weight of CH13.000000.
W =
13
molecular weight of CH214.000000.
molecular weight of CH214.000000.
W =
14
molecular weight of CH2(S)14.000000.
molecular weight of CH2(S)14.000000.
W =
14
molecular weight of CH315.000000.
molecular weight of CH315.000000.
W =
15
molecular weight of CH416.000000.
molecular weight of CH416.000000.
W =
16
molecular weight of CO26.000000.
molecular weight of CO26.000000.
W =
26
molecular weight of CO240.000000.
molecular weight of CO240.000000.
W =
40
molecular weight of HCO27.000000.
molecular weight of HCO27.000000.
W =
27
molecular weight of CH2O28.000000.
molecular weight of CH2O28.000000.
W =
28
molecular weight of CH2OH29.000000.
molecular weight of CH2OH29.000000.
W =
29
molecular weight of CH3O29.000000.
molecular weight of CH3O29.000000.
W =
29
molecular weight of CH3OH30.000000.
molecular weight of CH3OH30.000000.
W =
30
molecular weight of C2H25.000000.
molecular weight of C2H25.000000.
W =
25
molecular weight of C2H226.000000.
molecular weight of C2H226.000000.
W =
26
molecular weight of C2H327.000000.
molecular weight of C2H327.000000.
W =
27
molecular weight of C2H428.000000.
molecular weight of C2H428.000000.
W =
28
molecular weight of C2H529.000000.
molecular weight of C2H529.000000.
W =
29
molecular weight of C2H630.000000.
molecular weight of C2H630.000000.
W =
30
molecular weight of CH2CO40.000000.
molecular weight of CH2CO40.000000.
W =
40
molecular weight of HCCO39.000000.
molecular weight of HCCO39.000000.
W =
39
molecular weight of HCCOH40.000000.
molecular weight of HCCOH40.000000.
W =
40
molecular weight of H2CN30.000000.
molecular weight of H2CN30.000000.
W =
30
molecular weight of HCN29.000000.
molecular weight of HCN29.000000.
W =
29
molecular weight of HNO31.000000.
molecular weight of HNO31.000000.
W =
31
molecular weight of N16.000000.
molecular weight of N16.000000.
W =
16
molecular weight of NNH33.000000.
molecular weight of NNH33.000000.
W =
33
molecular weight of N2O46.000000.
molecular weight of N2O46.000000.
W =
46
molecular weight of NH17.000000.
molecular weight of NH17.000000.
W =
17
molecular weight of NH218.000000.
molecular weight of NH218.000000.
W =
18
molecular weight of NH319.000000.
molecular weight of NH319.000000.
W =
19
molecular weight of NO30.000000.
molecular weight of NO30.000000.
W =
30
molecular weight of NO244.000000.
molecular weight of NO244.000000.
W =
44
molecular weight of HCNO43.000000.
molecular weight of HCNO43.000000.
W =
43
molecular weight of HOCN43.000000.
molecular weight of HOCN43.000000.
W =
43
molecular weight of HNCO43.000000.
molecular weight of HNCO43.000000.
W =
43
molecular weight of NCO42.000000.
molecular weight of NCO42.000000.
W =
42
molecular weight of CN28.000000.
molecular weight of CN28.000000.
W =
28
molecular weight of HCNN45.000000.
molecular weight of HCNN45.000000.
W =
45
molecular weight of N232.000000.
molecular weight of N232.000000.
W =
32
molecular weight of AR40.000000.
molecular weight of AR40.000000.
W =
40
molecular weight of C3H844.000000.
molecular weight of C3H844.000000.
W =
44
molecular weight of C3H743.000000.
molecular weight of C3H743.000000.
W =
43
molecular weight of CH3CHO42.000000.
molecular weight of CH3CHO42.000000.
W =
42
molecular weight of CH2CHO41.000000.
molecular weight of CH2CHO41.000000.
W =
41
>> 
RESULT :
FINAL PLOT OF O2 SPECIE:
FINAL PLOT OF N2 SPECIE:
FINAL PLOT OF Co2 SPECIE :
Screenshot of the window where all the folders are saved.
error :
get = fget1(f1). fget'1' is wrongly used. corrected as get = fgetl(f1).
CONCLUSION:
The file parsing is done using NASA Thermodynamics Data and the properties cp, H, S are calculated.