Parsing NASA thermodynamic data
Introduction:
Parsing the data or information means breaking down the data into component parts so that its syntax can be analyzed, categorized and understood. Parsing is complementary to templating, which produces formatted output. These may be applied to different domains, but often appear together, such as the scanf/printf pair, or the input and output stages.
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. There are 53 species in the Thermodynamic data file and we require three plots per species.
Aim:
To parse the NASA Thermodynamic data and plot the graphs of for the local temperature range.
Objective:
- To write a code that extracts the 14 co-efficients and calculate Cp, H and S for all the species in the data file(THERMO.dat).
- To plot the graphs of Cp, H and S for the local temperature range for each species.
- To save the plots as images in separate folder for each species.
- To Calculate the Molecular weight of each species.
Formulae:
CODE:
%Program to extract the 14 co-efficients and calculate Cp, H and S for all the species in the data file
clear all
close all
clc
%Opening the file
f1 = fopen('THERMO.dat','r');
fgetl(f1);
%Read the header file data
H = fgetl(f1);
%In THERMO.dat there is space between the string,we can eliminate it
%Split the lines using a delimiter
H = strsplit(H,' ');
%converting strings and using the data in cells.
global_temp_low = str2double(H{1});
global_temp_mid = str2double(H{2});
global_temp_high = str2double(H{3});
%Next three comments are skipped after the temperature range.
for i = 1:3
fgetl(f1);
end
%Working with the species
for i = 1:53 %No. of species is 53.
header = fgetl(f1); %Skipping the header(O)
spheader = strsplit(header,' ' );
species_name = spheader{1};
en = length(spheader); %en is the last column of the local temp line.
t_high = str2double(spheader(en-2));
t_mid = str2double(spheader(en-1));
t_low = str2double(spheader(en-3));
%Extract the first five coefficients of first species
line_1 = fgetl(f1);
a = strfind(line_1,'E');
a1 = str2double(line_1(1:a(1)+3));
a2 = str2double(line_1(a(1)+4:a(2)+3));
a3 = str2double(line_1(a(2)+4:a(3)+3));
a4 = str2double(line_1(a(3)+4:a(4)+3));
a5 = str2double(line_1(a(4)+4:a(5)+3));
%Extract the Next five coefficients of first species
line_2 = fgetl(f1);
a_1 = strfind(line_2,'E');
a6 = str2double(line_2(1:a_1(1)+3));
a7 = str2double(line_2(a_1(1)+4:a_1(2)+3));
a8 = str2double(line_2(a_1(2)+4:a_1(3)+3));
a9 = str2double(line_2(a_1(3)+4:a_1(4)+3));
a10 = str2double(line_2(a_1(2)+4:a_1(5)+3));
%Extract the Next four coefficients of first species
line_3 = fgetl(f1);
a_2 = strfind(line_3,'E');
a11 = str2double(line_3(1:a_2(1)+3));
a12 = str2double(line_3(a_2(1)+4:a_1(2)+3));
a13 = str2double(line_3(a_2(2)+4:a_1(3)+3));
a14 = str2double(line_3(a_2(3)+4:a_1(4)+3));
%Creating a folder
mkdir(['D:Danish MatlabFile parsing images',species_name]);
%If T>1000 high temp coeff and if T<1000 Low temp coeff
T = linspace(t_low,t_high,100);
R = 8.314;
%Calculating Cp, H, and S
for x = 1:length(T)
%Finding high temp coefficients
if T(x) > t_mid
Cp = (a1 +( a2.*T) + (a3.*T.^2) + (a4.*T.^3) + (a5.*T.^4))*R;
H = (a1 + (a2.*T)/2 + (a3.*T.^2)/3 + (a4.*T.^3)/4 + (a5.*T.^4)/5 + (a6./T))*R.*T;
S = a1.*log(T) + (a2.*T) + (a3.*T.^2)/2 + (a4.*T.^3)/3 + (a5.*T.^4)/4 + a7;
else
%Finding Low temp coefficients
Cp = (a8 + a9.*T + a10.*T.^2 + a11.*T.^3 + a12.*T.^4)*R;
H = (a8 + (a9.*T)/2 + (a10.*T.^2)/3 + (a11.*T.^3)/4 +(a12.*T.^4)/5 + (a13./T))*R.*T;
S = a8.*(log(T)) + (a9.*T) + (a10.*T.^2)/2 + (a11.*T.^3)/3 + (a12.*T.^4)/4 + a14;
end
end
%Plotting the graphs
figure(1)
plot(T,Cp,'linewidth',3)
grid on
title('Temperature Vs Specific Heat');
xlabel('Temperature (K)');
ylabel('Specific heat');
figure(2)
plot(T,H,'linewidth',3)
grid on
title('Temperature Vs Enthalpy');
xlabel('Temperature (K)');
ylabel('Enthalpy');
figure(3)
plot(T,S,'linewidth',3)
grid on
title('Temperature Vs Entropy');
xlabel('Temperature (K)');
ylabel('Entropy');
%As we have many plots, so we save the plots in a folder.
cd(['D:Danish MatlabFile parsing images', species_name]);
saveas(1,'Temperature Vs Specific Heat','png')
saveas(2,'Temperature Vs Enthalpy','png')
saveas(3,'Temperature Vs Entropy','png')
%Calculating the moleculer weight of each species.
species = species_name;
atoms = ['H','C','O','N','A'];
atomic_wt = [1.0079, 12.0107, 15.0999, 14.0067, 39.948];
molmass = 0;
for i = 1:length(species)
for j = 1:length(atoms)
if strcmp(species(i),atoms(j)); %Comparing the strings.
molmass = molmass + atomic_wt(j);
pos = j;
end
end
end
N = str2num(species(i));
if N > 1
molmass = molmass + (atomic_wt(pos)*(N-1));
end
disp(sprintf('Molecular Weight of %s is %f ',species, molmass))
end
fclose(f1)
%End of Program
Explanation of the Code:
Initially we open the file THERMO.dat using the command ‘fopen’.
fileID = fopen(filename,permission) opens the file with the type of access specified by permission.
Then we skip the first line of the file, which is just a comment, we use ‘fgetl’ to skip a line.
tline = fgetl(fileID) returns the next line of the specified file, removing the newline characters.
To eliminate the space between the strings from the data file in order to acquire the required data, we use ‘strsplit’.
C = strsplit(str,delimiter) splits str at the delimiters specified by delimiter.
We then obtain the global temperature ranges from the data file by getting the data from a line and assigning particular columns to their respective variables.
Since there are 53 species, it is necessary for us to use loop to obtain the plots. We eliminate the space as well as the strings of the line and then we use the temperature values and assign those to respective variables(high, mid and low).
Now we extract the co-efficients from the data file. Set the index value and assign it to a variable, here we use ‘strfind’ to set it.
k = strfind(str,pattern) searches str for occurrences of pattern.
We can now find the variables with respect to the indexing variable that is already assigned. All the 14 co-efficients are found in a similar way.
Make a folder using ‘mkdir’ creates a path to save all you code data along with sapecies name.
mkdir parentFolder folderName creates folderName in parentFolder. If parentFolder does not exist, MATLAB attempts to create it.
As the middle temperature range is 1000, with reference to it we can get our higher and lower temperature co-efficients. If T>1000 high temp coeff and if T<1000 Low temp coeff.
Using the formulae we find Specific heat, enthalpy and entropy.
To plot the required graphs we use ‘plot’ command. The plots are for Temperature Vs Specific Heat, Temperature Vs enthalpy, Temperature Vs Entropy.
Save the plots using the command ‘saveas’ to the folders of respective species.
Lastly Calculate the molecular weight of each species. The atomic weights of the elements :H, C, O, N, A are assigned to variables. For this we use loop and compare the values of the two variable. If the Molecular weight is greater than 1 then the Molecular weight of the particular species is obtained.
Errors Occurred:
There was an error shown in the command window, the variable 'spheaders' was replaced by the main variable 'spheader' in order to solve the error.
Result:
The Plots for all the 53 species are obtained and the Molecular Weight of each species are shown below.
Plot for Temp. vs sp. heat (O2) :
Plot for Temp. vs Enthalpy(O2) :
Plot for Temp. vs Entropy(O2) :
Molecular weight output :
Molecular Weight of O is 15.099900
Molecular Weight of O2 is 30.199800
Molecular Weight of H is 1.007900
Molecular Weight of H2 is 2.015800
Molecular Weight of OH is 16.107800
Molecular Weight of H2O is 16.107800
Molecular Weight of HO2 is 31.207700
Molecular Weight of H2O2 is 31.207700
Molecular Weight of C is 12.010700
Molecular Weight of CH is 13.018600
Molecular Weight of CH2 is 14.026500
Molecular Weight of CH2(S) is 13.018600
Molecular Weight of CH3 is 15.034400
Molecular Weight of CH4 is 16.042300
Molecular Weight of CO is 27.110600
Molecular Weight of CO2 is 42.210500
Molecular Weight of HCO is 28.118500
Molecular Weight of CH2O is 28.118500
Molecular Weight of CH2OH is 29.126400
Molecular Weight of CH3O is 28.118500
Molecular Weight of CH3OH is 29.126400
Molecular Weight of C2H is 13.018600
Molecular Weight of C2H2 is 14.026500
Molecular Weight of C2H3 is 15.034400
Molecular Weight of C2H4 is 16.042300
Molecular Weight of C2H5 is 17.050200
Molecular Weight of C2H6 is 18.058100
Molecular Weight of CH2CO is 40.129200
Molecular Weight of HCCO is 40.129200
Molecular Weight of HCCOH is 41.137100
Molecular Weight of H2CN is 27.025300
Molecular Weight of HCN is 27.025300
Molecular Weight of HNO is 30.114500
Molecular Weight of N is 14.006700
Molecular Weight of NNH is 29.021300
Molecular Weight of N2O is 29.106600
Molecular Weight of NH is 15.014600
Molecular Weight of NH2 is 16.022500
Molecular Weight of NH3 is 17.030400
Molecular Weight of NO is 29.106600
Molecular Weight of NO2 is 44.206500
Molecular Weight of HCNO is 42.125200
Molecular Weight of HOCN is 42.125200
Molecular Weight of HNCO is 42.125200
Molecular Weight of NCO is 41.117300
Molecular Weight of CN is 26.017400
Molecular Weight of HCNN is 41.032000
Molecular Weight of N2 is 28.013400
Molecular Weight of AR is 39.948000
Molecular Weight of C3H8 is 20.073900
Molecular Weight of C3H7 is 19.066000
Molecular Weight of CH3CHO is 41.137100
Molecular Weight of CH2CHO is 41.137100
Folders for the plots:
Conclusion:
The NASA Thermodynamic file has a lot of data in it, we cannot handle it just by writing a code, file parsing is helpful in this case. File parsing allows us to handle larger data files much easily. Automation of any kind of data processing is easier using the file parsing.
Thus, the folders for the plots of each species is systematically alloted and can be viewed whenever required.