Parsing NASA thermodynamic data using MATLAB

Aim: NASA thermodynamic data file parcing using MATLAB and calculate the thermodynamic properties of various gas species.

Theory:

File Parsing:

It is the splitting the file or other input into a piece of data that can be easily stored or manipulated. By doing this we are going to split a string into parts then recognizing the parts to convert it into something simpler.

Project we are going to split the NASA thermodynamic data which is having 53 gas species values in 213 lines. To perform the calculation for thermodynamic properties, we need some sort of co-efficient values and temperature ranges which is available in the NASA file. So, we are going to parse the values which are all needed for our calculation from file.

Formula calculating the Cp, H, S:

here formule,

(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,

Cp = Specific Heat (J/K kg)

H = Enthalpy (J)

S = Entropy (J/K)

R = Universal gas constant = 8.314 (J/mol.K)

T = Temperature (K)

Code:

clear all
close all
clc

f1 = fopen('THERMO.dat','r')

line = fgetl(f1)
globel_temp = fgetl(f1)

for i= 1:3
    line = fgetl(f1)
    
end
for i = 1:53
line = fgetl(f1);
line_1 = strsplit(line,' ')

specie_name = line_1{1};
m_t = str2double(line_1(length(line_1)-1));
h_t = str2double(line_1(length(line_1)-2));
l_t = str2double(line_1(length(line_1)-3));

line_1 =fgetl(f1);
pos = strfind(line_1,'E');

a1  = str2double(line_1(1:pos(1)+3))
a2 = str2double(line_1(pos(1)+4:pos(2)+3))
a3 = str2double(line_1(pos(2)+4:pos(3)+3))
a4 = str2double(line_1(pos(3)+4:pos(4)+3))
a5 = str2double(line_1(pos(4)+4:pos(5)+3))

line_2 =fgetl(f1);
pos = strfind(line_2,'E');

a6  = str2double(line_2(1:pos(1)+3))
a7 = str2double(line_2(pos(1)+4:pos(2)+3))
a8 = str2double(line_2(pos(2)+4:pos(3)+3))
a9 = str2double(line_2(pos(3)+4:pos(4)+3))
a10 = str2double(line_2(pos(4)+4:pos(5)+3))

line_3 =fgetl(f1);
pos = strfind(line_3,'E');

a11  = str2double(line_3(1:pos(1)+3))
a12= str2double(line_3(pos(1)+4:pos(2)+3))
a13= str2double(line_3(pos(2)+4:pos(3)+3))
a14= str2double(line_3(pos(3)+4:pos(4)+3))
%%

T = linspace(l_t,h_t,1000)

for i = 1:length(T)
 [Cp(i),H(i),S(i)] = coeff(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T(i),m_t);
end



%%
elements = ['H','C','O','N','A','S'];                           % Declaring species elements
atomic_weight = [1.0079 12.0107 15.999 14.0067 39.948 32.065];  % Declaring the atomic weight value
mol_mass = 0;                                                   % Declaring molecular weight to be '0' in initial stage 

for l = 1:length(specie_name)
    for m = 1:length(atomic_weight)
        if strcmp(specie_name(l), elements(m))
            mol_mass = mol_mass + atomic_weight(m);
            position = m;
        end
    end
    n = str2double(specie_name(l));
    if n>1
        mol_mass = mol_mass + atomic_weight(position)*(n-1);
    end
end

% Printing the molecular weight of the gas species in the command window

sprintf('Molecular Weight of %s is %f u',specie_name,mol_mass) 


%%

mkdir(specie_name)
figure(1)
plot(T,Cp)

figure(2)
plot(T,H)

figure(3)
plot(T,S)
saveas(figure(1),'cp.jpg')
saveas(figure(2),'H.jpg')
saveas(figure(3),'S.jpg')
cd('D:matlabwork')


end
fclose(f1)

fuction for coifficients:

function [Cp,H,S] = coeff(c1,c2,c3,c4,c5,c6,c7,c8,c9,c10,c11,c12,c13,c14,T,m_t)

R = 8.314

if T>m_t
%inputs formule
Cp= (c1 + c2.*T + c3.*T.^2 + c4.*T.^3 +c5.*T.^4).*R;
H = (c1 + (c2.*T)./2 + (c3.*T.^2)./3 + (c4.*T.^3)./4 +(c5.*T.^4)./5 + c6./T).*T.*R;
S = (c1.*log(T) + c2.*T + (c3.*T.^2)./2 + (c4.*T.^3)./3 +(c5.*T.^4)./4  + c7).*R;

else
    
Cp= (c8 + c9.*T + c10.*T.^2 + c11.*T.^3 +c12.*T.^4).*R;
H = (c8 + (c9.*T)./2 + (c10.*T.^2)./3 + (c11.*T.^3)./4 +(c12.*T.^4)./5 + c13./T).*T.*R;
S = (c8.*log(T) + c9.*T + (c10.*T.^2)./2 + (c11*T.^3)./3 +(c12.*T.^4)./4  + c14).*R;    

end
end

Code explanation:

1. read the data from the external NASA thermodynamic data file using fopen command which is used to open the file and f command which fgetl is used to read the content in the file line by line.
2. After we use the fopen command we will be assigning to the variable called f. Thus, the fgetl(f) is the command acts like a line cursor, which will retrieve the next time while running this command for the second time and so on.
3. Now, we will be skipping the first 5 lines to enter the pattern of values of the respective gas species. We will be using the for loop command to process the same for all the gas species.
4. The name of the gas species we will be found out using findstr command and for finding out the species we are going to locate the string ‘ ’ which is present in after each species name.
5. And then we are going to use A(1:z(1)-1) which tells us the exact location of the species name, using which the name of the gas species is stored in the variable species_name.
6. We are going to extract these values similarly like we did for finding out the local temperatures by using the findstr and str2num commands and store it in the variables a1 to a14. For extracting we are going to locate the string ‘E’ which is present in the each co-efficient value and so it will be followed for all the co-efficients values in the same line and then we can proceed for the next lines.
7. The values from a1 to a7 represent the high temperature co-efficients and the values from a8 to a14 represent the low temperature co-efficients.
8. Temperature range values from low temperature to high temperature is assigned using the linspace command.
9. The for loop command is used to calculate the values of specific heat, enthalpy and entropy and the if-else command is also used to find out the required values if the temperature below the mid temperature or above the mid temperature.
10. In order to calculate the molecular weight of the species we declare the species elements in the form of cell array. Their atomic mass is declared in the form of array. Then initially molecular weight is declared as zero.
11. Now, for loop command is used to detect if there are more than one species in the variable species_name and the respective values of the gas species values are calculated and found out.
12. Now, we are going to create plots for the obtain specific heat, enthalpy and entropy values and save it in respective species_name folder automatically.
13. In order to create a new folder with the species_name we are going to use the mkdir command so that a new folder is created in the main folder. Now, we have to change the current directory from the main folder to the species_name folder so the cd command is used to do the same.
14. The plots for specific heat vs temperature, enthalpy vs temperature and entropy vs temperature is plotted and the labelling of axis and the naming of the title is also done along with the value of molecular weight.
15. To save the plots in the species_name folder, saveas command is used to save the plots in the current folder. After which to change the current directory to the main folder.

Outputs:

O2:

Species in folders:

Molecular weight:

Errors found:

Error in line number can be overcome by making corretion in line number 83 as mkir converted into mkdir.