BlogFor BusinessTalk to a counsellor
CoursesIIT ProgramsPlacementsLog in
BlogFor Business
Talk to a counsellor
Challenge

Project 1 - Parsing NASA thermodynamic data

YYogesh Jagtap10 May 2022
MATLAB

Aim: To write a code in matlab to parse the Nasa thermodynamic data file and calculate the thermodynamic properties of various gas species 

Objective:

  •  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.
  • 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.

Procedure:

File Parsing :

  •  Parsing,Syntax analysis is the process of analyzing the string of symbols,either in natural languages,computer languages or data structures,conforming to the rules of formal grammar.The term parsing comes from latin pars(orationis),meaning part (of speech)  
  • Within computational linguistic the term is used to refer to the formal analysis by a computer of a sentence or other string of words into its constituents, resulting in a parse tree showing their syntactic relation to each other, which may also contain semantic and other information.Some parsing algorithms may generate a parse forest or list of parse trees for a syntactically ambiguous input.
  • Parsing a file means reading the data stream of some sort and building an in memory model of semantic content of the data .This aims at facilitating and performing a transformation on the data.

 Specific heat (cp):   Specific heat means,the quantity of heat required to raise the temperature of one gram of a substance by one Celsius degree.

  the given polynomial for cp

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

where, cp=specific heat  

           R =Universal gas constant

           T= Temperature 

  Enthalpy (H): Enthalpy is the sum of internal energy and the product of pressure and volume.

the given polynomial for H 

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

  Where, H= Enthalpy 

             R =Universal gas constant

             T= Temperature 

  Entropy(S): Entropy, the measure of a system's thermal energy per unit temperature that is unavailable for doing useful work. Because work is obtained from ordered molecular motion, the amount of entropy is also a measure of the molecular disorder, or randomness, of a system.

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

      Where , S= Entropy

                  R =Universal gas constant

                  T= Temperature 

Main Matlab Code:

clear all
close all
clc
f1=fopen('THERMO (1).dat','r');
% reading the first line from data 
first_line=fgetl(f1); 
% Univesal gas constant 
R=8.314;

Comments1=fgetl(f1);
Comments2=fgetl(f1);    
Comments3=fgetl(f1);
comments4=fgetl(f1);
for i= 1:53

    line1= fgetl(f1);
    element_name=strtok(line1);
%obtaining Local max,min.,medium temperature
% obtaining temperature value from splitting data
temp_range=strsplit(line1,' ');
Local_min_temp = str2num(temp_range{end-3});
Local_mid_temp = str2num(temp_range{end-1});
Local_max_temp = str2num(temp_range{end-2});
% create a vector of 500 evenly spaced temperature points

T= linspace(Local_min_temp,Local_max_temp,500);
    line2= fgetl(f1);
    A=  findstr(line2,'E');

    line3= fgetl(f1);
    B=  findstr(line3,'E');
      
    line4= fgetl(f1);   
    C=  findstr(line4,'E');

    a1=str2num(line2(1:A(1)+3));
    a2=str2num(line2(A(1)+4:A(2)+3));
    a3=str2num(line2(A(2)+4:A(3)+3));
    a4=str2num(line2(A(3)+4:A(4)+3));
    a5=str2num(line2(A(4)+4:A(5)+3));

    a6=str2num(line3(1:B(1)+3));
    a7=str2num(line3(B(1)+4:B(2)+3));

    a8=str2num(line3(B(2)+4:B(3)+3));
    a9=str2num(line3(B(3)+4:B(4)+3));
    a10=str2num(line3(B(4)+4:B(5)+3));
    a11=str2num(line4(1:C(1)+3));
    a12=str2num(line4(C(1)+4:C(2)+3));
    a13=str2num(line4(C(2)+4:C(3)+3));
    a14=str2num(line4(C(3)+4:C(4)+3));

    % Specific heat 
    cp=specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,Local_mid_temp);

    % Enthalpy
    H=enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,Local_mid_temp);

    %Entropy
    S=entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,Local_mid_temp);

%calculating molecular weight of the species 
Molecular_wt(i)=molecular_weight(element_name);
warning('off')
% saving the plotting data
cd('C:\Users\yjagt\OneDrive\Documents\MATLAB\File parsing')
mkdir(element_name)
cd(element_name)


figure(1)
plot(T,cp,'linewidth',2,'color','b')
title(sprintf('Specific heat vs temperature for %s',element_name))
xlabel('Temperature')
ylabel('Specific heat')
grid on
saveas(figure(1),'Specific heat.png');

figure(2)
plot(T,H,'linewidth',2,'color','r')
title(sprintf('Enthalpy vs temperature for %s',element_name))
xlabel('Temperature')
ylabel('Enthalpy')
grid on
saveas(figure(2),'Enthalpy.png');

figure(3)
plot(T,S,'linewidth',2,'color','m')
title(sprintf('Entropy vs temperature for %s',element_name))
xlabel('Temperature')
ylabel('Entropy')
grid on
saveas(figure(3),'Entropy.png');

cd('C:\Users\yjagt\OneDrive\Documents\MATLAB\File parsing')
end

Explanation of above matlab code is 

Step 1) First use fopen command to open NASA  thermodynamic data .

Step 2) Use fgetl command to read the first line of the data ie ,THERMO, 

Step 3) Use the value of Universal gas constant (R),and by using the fgetl command access the comment line used in NASA  thermodynamic data .

Step 4) Use for loop to calculate the cp.H,S for all the 53 species and we use strtok command to get the selected part of the string (Species name)

Step 5)Obtaining Local max,minimun,medium temperature,obtaining temperature value from splitting data for that purpose we use strsplit command .convert the string into the number use str2num command  and we create variable for local minimun,medium,max temperature.

Step 6) Create a vector of 500 evenly spaced temperature points,by using linspace command .

Step 7)In line no 2 of the nasa thermodynamic data there is coefficient values ,In order to access this values we use fgetl command and assign variable line2 , findstr command for finding the position of  'E' in line 2 we use the variable A to store this values. similarly for line3,line4 we use variable B,C to store the value of position of 'E'.

Step 8) To extract the a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14 coefficient we use str2num command.

Step 9)To calculate the Specific heat,Enthalpy,Entropy,molecular weight of each species , we call the function code which is shown in above code.

Step 10) To save the plotting data we use mkdir and cd command .

Step 11)To create the plot of specific heat,enthalpy,entropy verses temperature we use plot command ,we use sprintf command in title to create the different species title.

Function code for specific heat calculation :

function cp = specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,Local_mid_temp)
% comparing with Local temperature
if T <= Local_mid_temp
   cp=R.*(a8+a9.*T+a10.*T.^2+a11.*T.^3+a12.*T.^4);
else 
   cp=R.*(a1+a2.*T+a3.*T.^2+a4.*T.^3+a5.*T.^4);
    
end

Function code for Enthalpy Calculation:     

function H=enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,Local_mid_temp);
% comparing with Local temperature
if T<= Local_mid_temp
    H=R.*(a8.*T+(a9.*T.^2)./2+(a10.*T.^3)./3+(a11.*T.^4)./4+(a12.*T.^5)./5+a13); 
else
    H=R.*(a1.*T+(a2.*T.^2)./2+(a3.*T.^3)./3+(a4.*T.^4)./4+(a5.*T.^5)./5+a6);
   
end

Function code for Entropy Calculation:

function S=entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,Local_mid_temp);
% comparing with Local temperature
if T<= Local_mid_temp
 S=R*(a8*log(T)+a9.*T+a10.*(T.^2)./2 +a11.*(T.^3)./3+a12.*(T.^4)./4+a14); 
else
  S=R*(a1*log(T)+a2.*T+a3.*(T.^2)./2 +a4.*(T.^3)./3+a5.*(T.^4)./4+a7);
  
end

Explanation of above three function code is :

  •   We use the if else command to calculate the value of specific heat,enthalpy,entropy as per the formula mentioned in the above code .

 Function Code For Molecular weight calculation:

function MW= molecular_weight(element_name)
atoms=['H','C','N','O','A'];
fprintf('Molecular weight of');
fprintf(element_name);
fprintf(':');
atomic_weight=[1 12 14 16 39];
MW=0;
for i=1:length(element_name)
    for j=1:length(atoms)
        if strcmp(element_name(i),atoms(j))
            MW=MW+atomic_weight(j);
            position =j;
        end
    end
    n=str2num(element_name(i));
    if n>1
        MW=MW+atomic_weight(position)*(n-1);
    end
end
fprintf('%d',MW);
disp(' ');
end

Explanation of above code :

  • For molecular weight calculation we create the function which takes atomic weight of hydrogen,carbon,nitrogen,oxygen,argon as inputs because all the elements are consist of all these atoms.
  • Now we use for loop and strcmp command is used to compare elements with species if they both match then molecular weight will be calculated.
  • Now ,take the string of species name and position to calculate their weight.here n is the number of molecules of the elements present in the species.  

Folder in which all the plots are saved:

 Plot of O2,N2,CO2:

Molecular Weight of each species :

References:

https://en.wikipedia.org/wiki/Parsing

https://in.mathworks.com/help/matlab/ref/fgetl.html

https://in.mathworks.com/help/matlab/ref/str2num.html

https://in.mathworks.com/help/matlab/ref/findstr.html

Your portfolio, next

Ship work like this
in your own course.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.