Parsing NASA Thermodynamic Data File & Using for Coefficient Extraction and Calculation using MATLAB

Objective: To parse the NASA Thermodynamic data file & use it to calculate various thermodynamic properties using MATLAB Scripting.

Theory:

What is File Parsing?

File Parsing is all about reading the information from the text file & how that information can be used in a computer language. By using file parsing we can decode the large file which is in another format. We can read the file and write down the file using the file parsing method. File Parsing can be used to do the calculations required for the specific purpose of the application. Parsing helps us make something understandable by analyzing its parts & bits. This is done by partially analyzing the data understanding its underlying structure (by making some assumptions based on what you're expecting to see), and then explicitly representing that structure in the code.

Thermo.dat file illustration:

NASA came up with polynomials that can be used to evaluate thermodynamic properties such as Cp, H, and S. They have also documented the coefficients that are required to evaluate these polynomials for 1000+ species.

Nasa Thermodynamic File can be found through below link

https://drive.google.com/file/d/1AikCHDlz9s_fu8whE81qee_GD4AyYpCP/view?usp=sharing

The NASA polynomials have the form:

    Cp/R = a1 + a2 T + a3 T^2 + a4 T^3 + a5 T^4
    
    H/RT = a1 + a2 T /2 + a3 T^2 /3 + a4 T^3 /4 + a5 T^4 /5 + a6/T
    
    S/R  = a1 lnT + a2 T + a3 T^2 /2 + a4 T^3 /3 + a5 T^4 /4 + a7

where a1, a2, a3, a4, a5, a6, and a7 are the numerical coefficients supplied in the file. The first 7 numbers starting on the second line of each species entry (five of the second line and the first two of the third line) are the seven coefficients (a1 through a7, respectively) for the high-temperature range (above 1000 K, the upper boundary is specified on the first line of the species entry). The following seven numbers are the coefficients (a1 through a7, respectively) for the low-temperature range (below 1000 K, the lower boundary is specified on the first line of the species entry).

Function Coding for extracting the values of 14 coefficients and local temperatures:

• The values of local temperature and theie respective coefficients are extracted by the below function. In this function LT, HT & MT is local low,high & mid temperature for specific species and co_1, co_2, co_3, co_4, co_5, co_6, co_7, co_8, co_9, co_10, co_11, co_12, co_13 & co_14 are coefficients for each species.
• This is extracted by using 2 commands first is findstr which is used for finding the specific word/letter from the text and second is str2num which converts the string form of the digit to number. We can also use str2double to extract the number instead of str2num.
• Here co_1, co_2, co_3, co_4, co_5, co_6, co_7, co_8 are High Coefficients & co_7, co_8, co_9, co_10, co_11, co_12, co_13 & co_14 are Low Coeffiecients.
• This function is used in the main program to plot the graph for species.

function [elements,LT,HT,MT,co_1,co_2,co_3,co_4,co_5,co_6,co_7,co_8,co_9,co_10,co_11,co_12,co_13,co_14] = Coefficients(Nasa,line_1);

%Local temperatures extraction from q line of each element
elements = line_1(1);
LT = str2num(line_1{1,(length(line_1))-3}); %Local Low Temp
HT = str2num(line_1{1,(length(line_1))-2}); %Local High Temp
MT = str2num(line_1{1,(length(line_1))-1}); %Local Medium Temp

%5 coefficient extraction from line no.2
Line_2 = fgetl(Nasa);
N= findstr(Line_2, 'E');
co_1 = Line_2(1:(N(1))+3);
co_1 = str2num(co_1);
co_2 = Line_2((N(1)+4:(N(2))+3));
co_2 = str2num(co_2);
co_3 = Line_2((N(2)+4:(N(3))+3));
co_3 = str2num(co_3);
co_4 = Line_2((N(3)+4:(N(4))+3));
co_4 = str2num(co_4);
co_5 = Line_2((N(4)+4:(N(5))+3));
co_5 = str2num(co_5);

%5 coefficient extraction from line no.3
Line_3 = fgetl(Nasa);
O = findstr(Line_3, 'E');
co_6 = Line_3((1:(O(1))+3));
co_6 = str2num(co_6);
co_7 = Line_3((O(1)+4:(O(2))+3));
co_7 = str2num(co_7);
co_8 = Line_3((O(2)+4:(O(3))+3));
co_8 = str2num(co_8);
co_9 = Line_3((O(3)+4:(O(4))+3));
co_9 = str2num(co_9);
co_10 = Line_3((O(4)+4:(O(5))+3));
co_10 = str2num(co_10);

%4 coefficient extraction from line no.4
Line_4 = fgetl(Nasa);
P = findstr(Line_4,'E');
co_11 = Line_4((1:(P(1))+3));
co_11 = str2num(co_11);
co_12 = Line_4((P(1)+4:(P(2))+3));
co_12 = str2num(co_12);
co_13 = Line_4((P(2)+4:(P(3))+3));
co_13 = str2num(co_13);
co_14 = Line_4((P(3)+4:(P(4))+3));
co_14 = str2num(co_14);

end

Function Coding for finding the values of Specific Heat, Enthalpy & Entropy:

• This function is used to find the values for Specific Heat(CP), Entropy(S) & Enthalpy of all the species present in the NASA Thermodynamic file using the below formulae.

`(CP)/R = a1 + a2 T + a3 T^2 + a4 T^3 + a5 T^4`

`(H)/(RT) = a1 + (a2 T) /2 + (a3 T^2) /3 + (a4 T^3) /4 + (a5 T^4) /5 + (a6)/T`

`(S)/R = a1 lnT + a2 T +( a3 T^2) /2 + (a4 T^3) /3 +( a5 T^4) /4 + a7`

Here a1,a2,a3,a4,a5,a6,a7 are replaced by coefficients as co_1, co_2, co_3, co_4, co_5, co_6, co_7, co_8, co_9, co_10, co_11, co_12, co_13 & co_14.

• R is Universal Gas Constant and its value is 8.314.
• We have given conditions that if T is greater than Local Mid Temperature then the function should use the High Coefficients and if T is lesser than local Mid Temperature then the function should use the Low Coefficients values.
• This function is used in the main program to plot the graph of Entropy, Specific Heat and Enthalpy for species.

function [CP,H,S] = Thermal_Reln(T,MT,co_1,co_2,co_3,co_4,co_5,co_6,co_7,co_8,co_9,co_10,co_11,co_12,co_13,co_14);

%Value for Universal Gas Constant
R = 8.314;

%Use first 7 coefficients if T > Medium Temp, else use the next 7 coefficients
if T >= MT
CP = R*(co_1 + (co_2*T) + (co_3*(T^2)) + (co_4*(T^3)) + (co_5*(T^4)));
H = R*T*(co_1 + ((co_2*T)/2) + ((co_3*(T^2))/3) + ((co_4*(T^3))/4)+ ((co_5*(T^4))/5) + (co_6/T));
S= R*((co_1*(log(T))) + co_2*T + ((co_3*(T^2))/2) +((co_4*(T^3))/3) + ((co_5*(T^4))/4) + (co_7));
else
CP = R*(co_8 + co_9*T + (co_10*(T^2)) + (co_11*(T^3)) + (co_12*(T^4)));
H = R*T*(co_8 + ((co_9*T)/2) + ((co_10*(T^2))/3)+ ((co_11*(T^3))/4) + ((co_12*(T^4))/5) + (co_13/T));
S = R*( (co_8*(log(T))) + (co_9*T) + ((co_10*(T^2))/2) +((co_11*(T^3))/3) + ((co_12*(T^4))/4) + (co_14));
end
end

Steps of Main Program Coding:

1) First, we use fopen to open our Thermo.dat file with permission to read the file with 'r'.

2) Next, we use fgetl to read the next line/skip the next line in the script window. Here we skip the header line.

3) Now to extract the global temperatures from line no 2 we use fscanf to scan the floating numbers and store the global temperatures data.

4) Again fgetl is used to skip the next 3 lines which were of no use for our purpose.

5) Next, we provide input to ask for the species name from a user for which we have to extract the data from the NASA file and plot the graphs.

6) As we know we have to repeat the same process to extract the values for each species we use contour in the while loop to do so and this while loop lasts till the end of the file to extract the data easily and reduce the solver time.

7) In a while loop, we use fgetl to read the first line of species and then strsplit to split the string in the first line into the cell line element.

8) Then we call the coefficient function to provide the values for coefficients and local temperatures. The local temperature values of low and high temperatures are used to create a range of 1000 values.

9) Next, if structure is used to determine the string input with the user to compare the element with their names after the function of the coefficient. If the state is true then the program will execute else the contour is increased by 1 & the while loop runs again.

10) Once the string comparison stands true then the atomic weight is calculated using for loop and the values of the weight are printed using fprintf command.

11) Further, we call the Thermal_Reln function to calculate the values of CP, H, and S which we described before.

12) As the program gets executes after the if conditions find to be true then the next part is to create new folders for storing the files of species.

13) For this, we used mkdir to create a new folder with the same name as species name for which the program is run.

14) After creating the folders the working directory is changed to the new folder and the plotted graphs are stored in it.

15) A new file with a molecular weight text file with .txt extension is opened and text is scripted in that file with the molecular weight of that species. The file is then closed and the working directory is changed again to the main directory.

16) Once the program is executed the loop is stopped using the break statement and the file is closed.

Main Program:

%Opening the file to read the data
Nasa = fopen('THERMO.dat','r');

%Skipping the header line
fgetl(Nasa);

%Extracting the data of global temperatures from line no 2
M = fscanf(Nasa,'%f',[2,inf]);
M = M(:);
gt_low = M(1); %global low temp
gt_mid_= M(2); %global mid temp
gt_high = M(3); %global high temp

%Skipping the next three lines
fgetl(Nasa);
fgetl(Nasa);
fgetl(Nasa);

%Species name input
Species = input('Name of species for which you want to find the data?','s');
ct = 0;
while ~feof(Nasa) %feof - for end of file
Line = fgetl(Nasa);%Read first line of each element
Line_1 = strsplit(Line,' '); % Splitting of first line

%Funtion calling to extract the data of 14 coefficients and their local temperatures
[elements,LT,HT,MT,co_1,co_2,co_3,co_4,co_5,co_6,co_7,co_8,co_9,co_10,co_11,co_12,co_13,co_14] = Coefficients(Nasa,Line_1);
T = linspace(LT,HT,1000);

%Cheking for a match between user input and file data
if strcmpi(Species,elements)
    
%Column array to char array conversion
elements = char(elements);
Q = ["H","N","O","C","Ar"]; %elements in data file
Atomic_weight = [1,14,16,12,40]; %respective atomic weights
weight = 0;

%Finding the molecular weight of the species
for s = 1: length(elements)
for t = 1: length(Q)
if strcmpi(elements(s),Q(t))
weight = weight + Atomic_weight(t);
a = Atomic_weight(t);
end
end
n = str2num(elements(s));
if n > 1
weight = weight + (a*(n-1));
end
end

%Printing the value of a moelcular weight
fprintf("Molecular weight of the species is %f",weight);

% Finding the values of CP,H,S using Thermal Relation Function
for i = 1:length(T)
[CP,H,S] = Thermal_Reln(T(i),MT,co_1,co_2,co_3,co_4,co_5,co_6,co_7,co_8,co_9,co_10,co_11,co_12,co_13,co_14);
CP_1(i)=CP;
H_1(i)=H;
S_1(i)=S;
end

%Converting column array to string
foldername = sprintf('%s',elements);
mkdir(foldername) %Creating folder for specific species
newfolder = cd(foldername);%Changing directory to new folder

%Creating and Saving the plot for T vs CP
figure(1)
plot(T,CP_1,'linewidth',2,'color','g');
xlabel('Temprature [K]');
ylabel('Specific Heat [KJ/Kmol.K]');
title('Plot between Temp & CP');
saveas(figure(1),'Specific_Heat.png');

%Creating and Saving the plot for T vs H
figure(2)
plot(T,H_1,'linewidth',2,'color','g');
xlabel('Temprature [K]');
ylabel('Enthalpy');
title('Plot between Temp & Enthalpy');
saveas(figure(2),'Enthalpy.png')

%Creating and Saving the plot for T vs S
figure(3)
plot(T,S_1,'linewidth',2,'color','g');
xlabel('Temprature [K]');
ylabel('Entropy');
title('Plot between Temp & Entropy');
saveas(figure(3),'Entropy.png')

%Creating a new .txt file and with the permission to write the file
%Displaying the molecular weight of the species and close a file
m_w = fopen('Molecular_Weight.txt','w');
fprintf(m_w,'Molecular weight of %s is %d',elements,weight);
fclose(m_w);

%Changing the directory back to the old folder
oldFolder = cd ('..');

%Breaking the program once the condition comes true
break
end
end
%close the data file
fclose(Nasa);

Command Window & other Outputs

Molecular Weight of the species in the command window

This is for species O. If we run the program again with the next species name we can see their molecular weight too.

Screenshot of the window where all the folders are saved

All the species files are stored in the hard disk after executing the program automatically in their respective species name folders.

Saving all the plots for each species and the text file of molecular weight in the folder

Molecular Weight Text File

Result Plots for O2, N2 & CO2 species:

1) Plots for O2

Fig - Plot of Temperature vs Specific Heat

Fig - Plot of Temperature vs Enthalpy

Fig - Plot of Temperature vs Entropy

2) Plots for N2

Fig - Plot of Temperature vs Specific Heat

Fig - Plot of Temperature vs Enthalpy

Fig - Plot of Temperature vs Entropy

3) Plots for CO2

Fig - Plot of Temperature vs Specific Heat

Fig - Plot of Temperature vs Enthalpy

Fig - Plot of Temperature vs Entropy

Conclusion:

Here we can conclude that File Parsing is the best approach to parse any file with large data and interpret the file using any Programming software like MATLAB for extracting the decoding the data available in the file. This type of File Parsing is very useful in the current era of industries where we transfer the large engineering data file from one place to another.