File parsing using NASA data

AIM-

1. 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 file data.

2. Calculate the molecular weight of each species and display it in the command window.

3. Plot the Cp, Enthalpy and Entropy for the local temperature range (low temperature : high temperature) specific for each species.

Introduction-

File parsing is reading data froma file and manipulating or using this data for our convienience.

Here we read data from a THEMO .dat file issued by nasa which contains data about few elements.

Data file-

https://drive.google.com/file/d/1AikCHDlz9s_fu8whE81qee_GD4AyYpCP/view.

(The first seven coeffiecints for each gas are for low temperature and the next seven are for high temperature)

The formulas we have are-

Where R is universal gas constant

Main program code-

The file is opened using fopen command in'read 'format to take in the data

clear all
close all
clc
f1 = fopen('THERMO.dat','r');

now we need to use the fgetl command to read each line of the file .the first two lines are read.

h1 =fgetl(f1);
h2 =fgetl(f1);

now we need to skip the next three lines as we can see in the file.so we read them and supress them.

for i = 1: 3
fgetl(f1);
end

Now we need to start another loop to read data till the end of the file.

Now we read line 6 contains element name ,as well as the temperatures we split the line wherever we encounter a space using strsplit command and store the first element into sp which is the element name.

To read and seperate our temperature values,now we use command findtsr(l6,'.') so that wherever in the array we encounter a . its position is recorded.

As we can see the first temp value is places before the first dot and last temp value is 23 places after the first dot.we read them into temprange.now inorder to spererate them into three different values we need to identify the spaces betwwen them using strsplit command and read the individual temperatures.

for i = 1: 53
l6 = fgetl(f1);
s = strsplit(l6,' ');
sp = s(1);
d = findstr(l6 ,'.');
temprange = l6((d(1)-3):(d(1)+23));
temp = strsplit (temprange,' ') ;
locallowtemp = str2num (temp {1});
templ = locallowtemp;
localhightemp = str2num(temp{2});
temp2 = localhightemp;
localmidtemp = str2num(temp{3});
temp3=localhightemp;

Now as we can see in the data file the line 7 ,line 8,line9 contains the coefficients for the specific gas. we read three of them using the fgetl command .now inorder to seperate the each coefficient we use strfind to find E so that using the position of E we can split the coefficients.

now using str2double we convert the string into integers so that we van make our futute caluclations.

coeff_line1=fgetl(f1);
coeff_line2=fgetl(f1);
coeff_line3=fgetl(f1);
A=strfind(coeff_line1,'E');
B=strfind(coeff_line2,'E');
C=strfind(coeff_line3,'E');
a1=coeff_line1(1:A(1)+3);
a1=str2double(a1);
a2=coeff_line1(A(1)+4:A(2)+3);
a2=str2double(a2);
a3=coeff_line1(A(2)+4:A(3)+3);
a3=str2double(a3);
a4=coeff_line1(A(3)+4:A(4)+3);
a4=str2double(a4);
a5=coeff_line1(A(4)+4:A(5)+3);
a5=str2double(a5);
a6=coeff_line2(1:B(1)+3);
a6=str2double(a6);
a7=coeff_line2(B(1)+4:B(2)+3);
a7=str2double(a7);
b1=coeff_line2(B(2)+4:B(3)+3);
b1=str2double(b1);
b2=coeff_line2(B(3)+4:B(4)+3);
b2=str2double(b2);
b3=coeff_line3(B(4)+4:B(5)+3);
b3=str2double(b3);
b4=coeff_line3(1:C(1)+3);
b4=str2double(b4);
b5=coeff_line3(C(1)+4:C(2)+3);
b5=str2double(b5);
b6=coeff_line3(C(2)+4:C(3)+3);
b6=str2double(b6);
b7=coeff_line3(C(3)+4:C(4)+3);
b7=str2double(b7);
molecule_name=char(sp) ;

Here a1 to a6 are for low temperature coefficients and b1 to b7 are high temperature coefficients.

now we have read our required data from file.

now we need to further caluclate the

Enthalpy

Entropy

Specific heat

molecular weight for each of the gases.

so for each specific molecule let us create a seperate file so that we can store its own data .we can do this by using mkdir command.

mkdir(['D:fileparsemolecule',molecule_name]);

now let us caluclate the cp,s,h values.

As we have already read our temperature ranges of each element we create a array taking the temperature ranges and considering 20 values between them so that we can get a smooth graph.

We have the constant value of universal gas constant R=8.314.

now let us run a loop for various various values of temperatures.

M = linspace (templ, temp2, 20);
T = int16(M(1:20));
R = 8.314;
k =1: length (T) ;

Also we have different formulas for cp,h,s for different temperature ranges. so we compare their temperature with global medium temperature in order to choose which temperature coefficients .if temperature is less tha the global average then we use the low temperature coeffiecients otherwise the high temperature coefficients .

for j=1: length (k)
if T(k) <= 1000
Cp = R*(b1 + (b2 * T(k)) + (b3*T (k).^2)+(b4*T (k).^3) + (b5*T(k).^4));
H= (R*T (k)).*( b1+(b2.*T(k)./2)+(b3.*T (k).^2 /3)+(b4.*T (k).^3/4) + (b5.*T(k).^4/5)+ (b6 /T(k))) ;
S = R*(b1.*log(double(T (k)) +(b2 .*T (k))+(b3.*T(k).^2/2)+(b4.* T (k).^3/3) + (b5.*T (k).^4/4) +b7));
else
Cp = R*(a1+(a2* T(k))+(a3*T(k).^2) +(a4 *T(k).^3)+(a5* T(k).^4));
H= (R* T (k)).*( a1+(a2*T(k))./ 2)+((a3*T(k).^2)./ 3) + (a4*T (k).^3./4) + ((a5*T(k).^4)./5) +(a6./T(k));
S = R *(a1.*log(double(T(k))+double((a2 .*T(k))+(a3.*T(k).^2./2)+(a4.*T(k).^3./3) +(a5.*T(k).^4./4) + a7)));
end
end

Now we have completed caluclation of cp,h,s and we need to plot them for each element

figure (1)
plot (M ,H,'-*');
title (molecule_name);
xlabel ('temperature')
ylabel (' enthalpy' )
figure (2)
plot (M, Cp, 'o')
title (molecule_name);
xlabel (' temperature ')
ylabel (' specific heat ')
figure (3)
plot (M, S, '-')
title (molecule_name);
xlabel (' temperature ')
ylabel (' entropy ')

Now we have three graphs for each element we need to save them in their respective folder which we have created for them.

so now we open the directory where our folder is created using the cd command and suing the saveas function we store the three plots in the coressponding molecule file

directory =(['D:fileparsemolecule', molecule_name]);
cd (directory );
saveas (1, 'Enthalpy.jpg');
saveas (2, 'Specificheat .jpg');
saveas (3,' Entropy.jpg' );
cd ('D:fileparse' );

NOW inorder to caluclate the molecular weight od each element we write a function molecular _mass.

we initialise mm=0;

as we can see that all the atoms in the molecule are one of H,C,N,O,A,S we make a array of their names as atoms and their coressponding atomic weights as an aw array..

now we need to run a for loop to check each element of the molecule , if the element matches then its coresspondig aw is added to the mm.

but for elements like O2,N2 where we need to conisder the mass of the element twice we use a if command. we convert the character into a double and check if it is >1 .

If the condition is true we add the aw of that specific element

function [mm]=molecular_mass(molecule_name)
fprintf('molecular weight of')
fprintf(molecule_name)
fprintf(':')
mm=0;
atoms=['H','C','N','O','A','S'];
aw=[1,12,14,16,40,32];
for i=1:length(molecule_name)
for j=1:length(atoms)
if ( strcmp(atoms(j),molecule_name(i))==1)
mm=mm+aw(j);
k=j;

end

end
factor=str2double(molecule_name(i));
if (factor>1)
mm=mm+((factor*aw(k))-aw(k)) ;

end

end
fprintf('%d',mm)
disp(' ');

Now we need to call this function into the main program and terminate the whole loop.

weight (i) = molecular_mass (molecule_name);
end

Program code-

clear all
close all
clc
f1 = fopen('THERMO.dat','r');
% header of f1
h1 =fgetl(f1);
h2 =fgetl(f1);
% removing lines from 3 to 5, since  the next line after header is 3rd line,
% indexing from 1 to 3 will give lines 3 to 5
for i = 1: 3 
    fgetl(f1); 
end
% 1:  53 loop since there are 53 species 
for i = 1: 53 
    l6 = fgetl(f1);
    s = strsplit(l6,' ');
   sp = s(1);
    d = findstr(l6 ,'.');
   temprange = l6((d(1)-3):(d(1)+23));
   temp = strsplit (temprange,' ') ;
   locallowtemp = str2num (temp {1}); 
 templ = locallowtemp;
 localhightemp = str2num(temp{2}); 
 temp2 = localhightemp; 
 localmidtemp = str2num(temp{3});
 temp3=localhightemp; 
 %coefficients from line 7 
 coeff_line1=fgetl(f1);
coeff_line2=fgetl(f1);
coeff_line3=fgetl(f1);
A=strfind(coeff_line1,'E');
B=strfind(coeff_line2,'E');
C=strfind(coeff_line3,'E');
a1=coeff_line1(1:A(1)+3);
a1=str2double(a1);
a2=coeff_line1(A(1)+4:A(2)+3);
a2=str2double(a2);
a3=coeff_line1(A(2)+4:A(3)+3);
a3=str2double(a3);
a4=coeff_line1(A(3)+4:A(4)+3);
a4=str2double(a4);
a5=coeff_line1(A(4)+4:A(5)+3);
a5=str2double(a5);
a6=coeff_line2(1:B(1)+3);
a6=str2double(a6);
a7=coeff_line2(B(1)+4:B(2)+3);
a7=str2double(a7);
b1=coeff_line2(B(2)+4:B(3)+3);
b1=str2double(b1);
b2=coeff_line2(B(3)+4:B(4)+3);
b2=str2double(b2);
b3=coeff_line3(B(4)+4:B(5)+3);
b3=str2double(b3);
b4=coeff_line3(1:C(1)+3);
b4=str2double(b4);
b5=coeff_line3(C(1)+4:C(2)+3);
b5=str2double(b5);
b6=coeff_line3(C(2)+4:C(3)+3);
b6=str2double(b6);
b7=coeff_line3(C(3)+4:C(4)+3);
b7=str2double(b7);
 molecule_name=char(sp) ;
  
  mkdir(['D:fileparsemolecule',molecule_name]);
  % solve for Cp, H, and S
  M = linspace (templ, temp2, 20); 
  T = int16(M(1:20));
  R = 8.314;
  k =1: length (T) ;
  for j=1: length (k)
      if T(k) <= 1000 
          Cp = R*(b1 + (b2 * T(k)) + (b3*T (k).^2)+(b4*T (k).^3) + (b5*T(k).^4));
          H= (R*T (k)).*( b1+(b2.*T(k)./2)+(b3.*T (k).^2 /3)+(b4.*T (k).^3/4) + (b5.*T(k).^4/5)+ (b6 /T(k))) ;
          S = R*(b1.*log(double(T (k)) +(b2 .*T (k))+(b3.*T(k).^2/2)+(b4.* T (k).^3/3) + (b5.*T (k).^4/4) +b7));
      else
          Cp = R*(a1+(a2* T(k))+(a3*T(k).^2) +(a4 *T(k).^3)+(a5* T(k).^4));
          H= (R* T (k)).*( a1+(a2*T(k))./ 2)+((a3*T(k).^2)./ 3) + (a4*T (k).^3./4) + ((a5*T(k).^4)./5) +(a6./T(k)); 
          S = R *(a1.*log(double(T(k))+double((a2 .*T(k))+(a3.*T(k).^2./2)+(a4.*T(k).^3./3) +(a5.*T(k).^4./4) + a7)));
      end
  end
  %plots 
  figure (1)
  plot (M ,H,'-*');
  title (molecule_name);
  xlabel ('temperature')
  ylabel (' enthalpy' ) 
  figure (2)
  plot (M, Cp, 'o')
  title (molecule_name);
  xlabel (' temperature ')
  ylabel (' specific heat ') 
  figure (3)
  plot (M, S, '-') 
  title (molecule_name);
  xlabel (' temperature ')
  ylabel (' entropy ') 
  directory =(['D:fileparsemolecule', molecule_name]);
  cd (directory );
  saveas (1, 'Enthalpy.jpg');
  % saving the plots in respective folder for each molecular species 
  saveas (2, 'Specificheat .jpg'); 
  saveas (3,' Entropy.jpg'  ); 
  cd ('D:fileparse' );
 weight (i) = molecular_mass (molecule_name);
end
 

Function code-

function [mm]=molecular_mass(molecule_name)
fprintf('molecular weight of')
fprintf(molecule_name)
fprintf(':')
mm=0;
atoms=['H','C','N','O','A','S'];
aw=[1,12,14,16,40,32];
for i=1:length(molecule_name)
    for j=1:length(atoms)
        if ( strcmp(atoms(j),molecule_name(i))==1)
            mm=mm+aw(j);
        k=j;   
             
        end
        
    end
        factor=str2double(molecule_name(i));
        if (factor>1)
           mm=mm+((factor*aw(k))-aw(k)) ;
        
end
    
end
fprintf('%d',mm)
    disp(' ');

Output-

Molecular weights of each elements-

Screenshot of files created for individual elements-

3.Cp,H,S graphs for-

1.02-

2.N2-

3.CO2-

Conclusion - Hence all the properties for various elemnts are plotted to know their variation with temperature.

Mistakes encountered and measures taken-

1. I made a error usind d(2) instead of d(1) because of which the program stopped saying index exceeded array bounds.

I have corrected the mistake by changing it to d(1)+23 so that the array bounds are satrisfioed.

2. i changed the adressv where all my files are saved.

initiallty i saved them in disk D but later i made a mistake by giving it as user disk C.

NOW I changed the directory to D and sucessfully elimintaed the error.

There are many other small errors which i have faced which include no ending my for function,improper or excess use

of paranthesis etc....