Project 1 - Parsing NASA thermodynamic data

AIM : To 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.

OBJECTIVE : To 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

• The FIRST 7 coefficients are HIGH-temperature coefficients and the SECOND 7 coefficients are LOW-temperature coefficients and also use the LOCAL TEMPERATURE RANGES which is unique for each individual species, not the GLOBAL TEMPERATURE RANGE.
• Calculate the molecular weight of each species and display it in the command window.
• Plot the Cp, Enthalpy and Entropy for the local temperature range (low temperature : high temperature) specific for each species.
• Save the plots as images with appropriate names and dump them in separate folders for each species with suitable name. All these should be generated automatically.
• The main program that you will submit will just take the "THERMO.dat" file as an input and generate Cp, H and S plots.

THEORY :

Parsing is the process of analyzing text made of sequence of tokens to determine its grammatical structure with respect to a given formal grammar. The praser then builds a data structure based on the tokens . Parsing is a compiler that is used to break the data into smaller elements coming from lexical analysis phase . A praser takes input in the form of sequence of tokens and produces output in the prase tree.

EXPLANATION :

1. Initially open the given thermodynamic data file and read the data in the file and store the data to the variable file.
2. Then read the header file of the thermodynamic data and calculate the global low, mid and high temperature from the second line of header file.
3. By using fgetl(),strsplit() and str2num() command find out the global low, mid and high temperature from the second line
4. Then calculate the 14 coefficients for all species in the given thermodynamic data file
5. In data file each species has three lines in which 14 coefficients are there
6. By using the for loop calculate the 14 coefficient for all the species.
7. By using strfind() commond find out the number of exponential in first line, by which we can find out the number of coefficient in first line.
8. Once you find out, then convert the coefficient from string from to number from by using str2num() command
9. By doing the above steps for line two and three we can calculate the 14 coefficient for all the species
10. Once you find out the coefficient for a species, then calculate the entropy, enthalpy and specific heat
11. First 7 coefficient are high temperature coefficient and calculate Cp, H and S
12. Assuming universal gas constant R = 8.314J/mol-k and calculating for each species by using the following formulas

    `Cp = R(C1 + C2.T + C3.T^2 +C4.T^3 +C5.T^4`

    `H = RT(C1 + C2. T/2 +C3. T^2/3 + C4. T^3/4 +C5. T^4/5 + (C6)/T)`

    `S = R(C1.ln(T) + C2. T + C3. T^2/2 + C4. T^3/3 + C5.T^4/4 + C7)`

13. Assuming temperatures range from 100 to 5000 take 100 values by linspace() command
14. Once you calculate the entropy, enthalpy and specific heat for all the species ,plot the graph between : T vs H , T vs S ,T vs Cp
15. Store all the three plots for all species in the file . The stored file is shown in the above picture
16. Finally calculating the molecular weight for all the species and displayed in the command window
17. The method for calculating molecular weight is shown in the function program

PROGRAM :

clear all
close all
clc

%reading the data file 
file = fopen('THERMO.dat','r');

%reading data file first line

line1 =  fgetl(file);
%Global temperatures
line2 = fgetl(file);
%skipping comments
line3 = fgetl(file);
line4 = fgetl(file);
line5 = fgetl(file);

%global temperature


A = strsplit(line2,' ' );
global_min_Temp=str2num(A{2});
global_mid_Temp=str2num(A{3});
global_max_Temp=str2num(A{4});



for i=1:53
    species = fgetl(file); %reading the line
    species_name = species(1:6)


%reading next line
Line_1 = fgetl(file);
    L1 = strfind(Line_1,'E');   
    C1 = str2num(Line_1(1:(L1(1)+3)));
    C2 = str2num(Line_1((L1(1)+4):(L1(2)+3)));
    C3 = str2num(Line_1((L1(2)+4):(L1(3)+3)));
    C4 = str2num(Line_1((L1(3)+4):(L1(4)+3)));
    C5 = str2num(Line_1((L1(4)+4):(L1(5)+3)));  
    %extracting low temp coeff
    Line_2 = fgetl(file);
    L2 = strfind(Line_2,'E');
    C6 =str2num(Line_2(1:(L2(1)+3)));
    C7 = str2num(Line_2((L2(1)+4):(L2(2)+3)));
    C8 = str2num(Line_2((L2(2)+4):(L2(3)+3)));
    C9 = str2num(Line_2((L2(3)+4):(L2(4)+3)))
    C10 = str2num(Line_2((L2(4)+4):(L2(5)+3)))
    
     Line_3 = fgetl(file);
    L3 = strfind(Line_3,'E');
    C11 =str2num(Line_3(1:(L3(1)+3)));
    C12 = str2num(Line_3((L3(1)+4):(L3(2)+3)));
    C13 = str2num(Line_3((L3(2)+4):(L3(3)+3)));
   C14 = str2num(Line_3((L3(3)+4):(L3(4)+3)))

 %creating folder for each species
 mkdir(['C:\matlab2\skilllyncplot\parsing\species\',species_name]);
   
   %Temperature range assumed
   Temp = linspace(300,5000,100);
   
   %universal gas constant
   R = 8.314;
   
   %calculations
   for j = 1:length(Temp)
   if Temp(j)<3000
            Cp = (C1+(C2.*Temp)+(C3.*(Temp.^2))+(C4.*(Temp.^3))+(C5.*(Temp.^4)))*R;
            H = (C1+(C2*(Temp/2))+(C3*((Temp.^2)/.3))+(C4*((Temp.^3)/4))+(C5*((Temp.^4)/.5))+(C6./Temp))*R.*Temp;
            S = ((C1.*log(Temp))+(C2.*Temp)+(C3.*((Temp.^2)./2))+(C4.*((Temp.^3)./3))+(C5.*((Temp.^4)./4))+C7)*R;
        else
             Cp = (C8+(C9.*Temp)+(C10.*(Temp.^2))+(C11.*(Temp.^3))+(C12.*(Temp.^4)))*R;
            H = (C8+(C9*(Temp/2))+(C10*((Temp.^2)/.3))+(C11*((Temp.^3)/.4))+(C12*((Temp.^4/.5))+(C13./Temp)))*R.*Temp;
            S = ((C8.*log(Temp))+(C9.*Temp)+(C10.*((Temp.^2)./2))+(C11.*((Temp.^3)./3))+(C12.*((Temp.^4)./4))+C14)*R;
        end
   end
    



%plotting
%Temperatures vs Entropy
            figure(1);
            plot(Temp,S,'linewidth',2,'color','r');
            xlabel('Temperature [K] ');
            ylabel('ENTROPY [J/KG]');
            title('specific heat vs entropy');
 %Temperature vs Enthalpy           
            figure(2);
             plot(Temp,H,'linewidth',2,'color','g');
            xlabel('Temperature [K] ');
            ylabel('ENTHALPY [J/KG)]');
            title('Enthapy vs Temperature');
            
            % Temperature (T) vs specific heat
            
            figure(3);
             plot(Temp,Cp,'linewidth',2,'color','b');
            xlabel('Temperature [K] ');
            ylabel('Specific  heat [J]');
            title('Specific heat vs Temperature');
            
            %changing thefilepath to save the plots
   directory = (['C:\matlab2\skilllyncplot\parsing\species\',species_name]); 
            cd(directory)
                
            %saving fig data
            saveas(1,'Entropy.jpg');
             saveas(2,'Enthalpy.jpg');
              saveas(3,'Specific heat.jpg');
              cd('C:\matlab2\') %switching the directory back to read thermodat file for molecular mass calclation
           
 M = molecularmass(species_name)
             
end

fclose(file);

ERRORS :

Function Program :

%Function for Molecular weight

function M = molecularmass(species_name);
species = upper(species_name);
 Constituents = {'H','C','N','O','A'};
 Atomic_weight = [1.00794 12.0107 14.0067 15.9994 39.948];
 M= 0;
 for i= 1:length(species)
     for j=1:length(Constituents)
         if strcmp(species(i),Constituents(j))
             M= M+ Atomic_weight(j);
             u = j;
         end
     end
     
     %Adding the Weight of the repeating atoms 
     n = str2double(species(i));
     if n>1
         M = M +Atomic_weight(u)*(n-1);         
     end
 end
 disp(['molecular_weight of:',species])
 disp(M)
 cd('C:\matlab2\')
end
 

• To find the molecular weight we are using six elements ['H' 'O' 'C' 'N' 'S' 'AR'] and their atomic weights [1.007,15.999,12.010,14.006,32.159,39.948]
• we are considering initial molecular weight as 0
• Using for loop and string compare command we comparing the species with respect to string array
• When the comparision becomes positive the program calculate the atomic weight to the old atomic weight the current species
• Then using the disp comment we are displaying the molecular weight of the species.

OUTPUT :

plots : O2

plots : N2

PLOTS : CO2

Molecular weight of species :