Parsing NASA's Thermodynamic Data file using Matlab program

Parsing NASA’s Thermodynamic Data file

Aim: To parse the NASA’s thermodynamics data file to extract the co-efficient and use it to calculate thermodynamic properties using Matlab program.

Theory:

File Parsing basically refers to how information can be read from particular file or text file and how that information can be used in a computer program. By Parsing a data file we are able to read data from external file or external database and we can store the data or include the data into our program as required.

In this project, we extract the co-efficient of polynomials that are stored in NASA’s Thermodynamics data file by performing file parsing and use this data in our program to calculate the thermodynamic properties like Specific heat, Enthalpy and Entropy using the following Equations:-

`(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*ln(T) + a2*T + a3*T^2/2 + a4*T^3/3 + a5*T^4/4 + a7`

Where Cp = Specific heat , H = Enthalpy, S = Entropy, T = Temperature, R = Universal Gas Constant,

a1,a2,a3,a4,a5,a6,a7 are co-efficient of polynomial which we need to extract from the data file.

Program Description:

(1). In this Program first we define the input value for Universal Gas constant which will be later used in the program. After that we start the file parsing by first Opening the ‘THERMO.dat’ file for reading the file.

(2). Then, we use fgetl command to ignore the comments and store the Temperatures in form of string and split the string to store the minimum, maxm and mid temperatures.

(3). After that, to extract the co-efficient for each species, we use for loop. Then we read the first three lines using fgetl command and store the data in form of string. Then to identify each co-efficient we trace the location of ‘E’ then extract and convert the co-efficient from string to numbers using str2double command.

(4). We distinguish the first 7 co-efficient as High temperature co-eff and next 7 as Low Temperature co-efficient.

(5). We create separate functions for calculation of Specific Heat, Enthalpy, Entropy and using the above co-efficient.

(6). After that we plot the Specific Heat, Enthalpy, Entropy vs Temperature.

(7). Now, we save the plots of each species in individual species folder. For that we make new subfolder using the ‘mkdir’ cmd by name of each species and save the jpg file of plots in it

(8). Now, we create separate function for calculation of Molecular Wt and return the values in main program to display the molecular wt in command window.

(9). Finally we close the file.

Matlab program (main file):

%Program for Parsing NASA's Thermodynamics Data file
clear all
close all
clc

%Defining input for universal gas constant
R = 8.314; % J/mol K
 
%Opening the thermo.dat file for reading data
 f1 = fopen('THERMO.dat','r');     
 l1 = fgetl(f1);                  %using fgetl cmd to ignore line
 l2 = fgetl(f1);                  %using fgetl cmd to store the temp data
 t = strsplit(l2);            %separating the string by splitting it into 3
 global_min = str2double(t{2});  %assigining the values of min,mid and
 global_mid = str2double(t{3});  %max temp after converting from string
 global_max = str2double(t{4});  % to numbers
 
 temp = linspace(global_min,global_max,400); %dividing the temp in 400parts
 
 % we can skip the next 3 header lines by using fgetl in for loop
 for i = 1:3
     fgetl(f1);
 end
 
 for i= 1:53
     %first reading the line with info of species
     line_1 = fgetl(f1);
     A = strsplit(line_1, ' ');
     Species_name = (A{1});
      % now extracting the 14 High and low temperature coefficients
    % Reading the three lines from which the co-efficents are extracted
     tline_2 = fgetl(f1);
     tline_3 = fgetl(f1);
     tline_4 = fgetl(f1);
     %finding the location of E to extract the coefficient as string
     a = strfind(tline_2,'E');    
     b = strfind(tline_3,'E');
     c = strfind(tline_4,'E');
     
     %extraction of 7 high temperature coefficients a1,a2,a3,a4,a5,a6,a7
     a1 = str2double(tline_2(1:a(1)+3));    
     a2 = str2double(tline_2(a(1)+4:a(2)+3));
     a3 = str2double(tline_2(a(2)+4:a(3)+3));     %extracting coeff by
     a4 = str2double(tline_2(a(3)+4:a(4)+3));     %converting into number
     a5 = str2double(tline_2(a(4)+4:a(5)+3));     %using str2double
     a6 = str2double(tline_3(1:b(1)+3));
     a7 = str2double(tline_3(b(1)+4: b(2)+3));
     
     %extraction of 7 low temperature coefficients b1,b2,b3,b4,b5,b6,b7 
     b1 = str2double(tline_3(b(2)+4: b(3)+3));
     b2 = str2double(tline_3(b(3)+4: b(4)+3));
     b3 = str2double(tline_3(b(4)+4: b(5)+3));
     b4 = str2double(tline_4(1:c(1)+3));
     b5 = str2double(tline_4(c(1)+4: c(2)+3));
     b6 = str2double(tline_4(c(2)+4: c(3)+3));
     b7 = str2double(tline_4(c(3)+4: c(4)+3));
     
 % Making new subfolder species to save the plots of individual species
 mkdir(['C:\Users\dell\Documents\MATLAB\Parsing_New\species\',Species_name]);
     
 %Calculation of Specific Heat by calling fn    
     C_p = specific_heat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,R,temp,global_mid);
 %Calculation of Enthalpy by calling fn    
       h = Enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,R,temp,global_mid);
 %Calvulation of Entropy by calling fn      
       S = entropy(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,R,temp,global_mid);

 %Plotting the values of Specific heat,enthalpy, entropy vs Temperature       
       figure(1)
       %Specific Heat vs Temperature plot
       subplot(3,1,1)
       plot(temp, C_p, 'linewidth', 4,'color', 'b');%plot based on the observed data
       title(Species_name)
       xlabel("Temperature")
       ylabel("Specific Heat")
       
       %Enthalpy vs Temperature plot
       subplot(3,1,2);
       plot(temp, h, 'linewidth', 4,'color', 'r');
       title(Species_name)
       xlabel("Temperature")
       ylabel("Enthalpy")
       
       %Entropy vs Temperature plot
       subplot(3,1,3);
       plot(temp, S, 'linewidth', 4,'color', 'g');
       title(Species_name)
       xlabel("Temperature")
       ylabel("Entropy")
       
%Saving the jpg file of plots for each species in individual species folder
  directory = (['C:\Users\dell\Documents\MATLAB\Parsing_New\species\',Species_name]);
  old_folder = cd(directory);    %saving old path to old_folder 
  saveas(1,'output1.jpg');       %saving the plot in species folder
  %changing current directory back to the original folder     
  cd(old_folder);
  %displaying values of molecular weight in the command window
  MW = molecular_weight(Species_name); 
      
 end
 %closing the file 
 fclose(f1)

Matlab program for Specific Heat (function):

%separate Function for Specific Heat Calculation
function [C_p] = specific_heat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,R,temp,global_mid)
 if temp <= global_mid
C_p = R*(b1 + (b2*temp)+ (b3*(temp).^2)+(b4 * (temp).^3) + (b5 *(temp).^4));
else
    C_p =  R*(a1 + (a2*temp)+ (a3*(temp).^2)+(a4 * (temp).^3) + (a5 *(temp).^4));
end

Matlab program for Enthalpy (function):

% Separate function for Enthalpy calculation
function [h] = Enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,R,temp, global_mid)
if temp <= global_mid
h = R.*temp.*(b1 + (b2.*(temp./2))+ (b3.*(temp./3))+(b4.*(temp./4)) + (b5.*(temp./5)) + (b6./temp));
 else 
    h =  R.*temp.*(a1 + (a2.*(temp./2))+ (a3.*(temp./3))+(a4.*(temp./4)) + (a5.*(temp./5)) + (a6./temp));
end

Matlab program for Entropy (function):

%Separate Funtion for Entropy Calculation
function [S] = entropy(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,R,temp, global_mid)
if temp <= global_mid
S = R.*(b1.*log(temp) + (b2.*temp)+ (b3.*(((temp).^2)./2))+ (b4.*(((temp).^3)./3)) + (b5.*(((temp).^4)./4))+ b7);
 else 
    S =  R.*(a1.*log(temp) + (a2.*temp)+ (a3.*(((temp).^2)./2))+ (a4.*(((temp).^3)./3)) + (a5.*(((temp).^4)./4))+ a7);
end

Matlab program for Molecular Weight (function):

%Separate Function for Molecular Weight Calculation
function [MW] = molecular_weight(Species_name)

element= {'H' 'C' 'N' 'O' 'Ar'};
atomic_weight = [1 12 14 16 40];

MW = 0;

for i = 1:length(Species_name)
    for j = 1:length(atomic_weight);
        if strcmpi(Species_name(i), element(j));
            MW = MW + atomic_weight(j);
            
            q = j;
        end
    end
    
    m = str2double(Species_name(i));
    if m>1
        MW = MW + (atomic_weight(q)*(m-1));
    end
    
end

Errors:

• This error occurred because the input for Universal Gas Constant ‘R’ was defined or missing. The error was rectified by defining the constant ‘R’ = 8.314.

Results:

• Plots of Specific heat, Enthalpy, Entropy vs Temperature for O2, N2, CO2 are as follows:-

For O2:

For N2:

For CO2:

• Screenshot of Folders of all species with plot saved as images is shown below:

• After execution of program we get the values of Molecular Weight of respective species in command section as follows:

Conclusion:

We have successfully extracted the co-efficient from the NASA’s thermodynamics data and used in our program to calculate the required properties and plots and saved the figures in individual Species folders using this program.