Challenge

Project 1 - Parsing NASA thermodynamic data

#ANSWER:

AIM:

To write a code in MATLAB to parse the NASA thermodynamic filedata and then calculate the thermodynamic properties of all the species.

objective:

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. ('THERMO.dat') file

The formulae are shown below.

where, R is the universal gas constant, T is the temeprature

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.

4) 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.

THEORY :

THEORY:-
FILE PARSING:-
File parsing is used to read the data from the external text file or any other type of file and to use the data inside the file in the program by extracting the data from the file is known as file parsing.


Specific Heat:-
The heat is required to raise the temperature of the unit mass of a given substance by a given amount. The units of specific heat are usually calorie or joules per gram per Celsius degree.

function code for specific heat:

% ========================= Specific Heat function file =========================== %

function Cp = specificheat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,Temp,R,global_medium_temp)

% comparing with global temperature %
% If less than global temp, min coefficents are taken %
% otherwise max coefficients are taken %

if Temp <= global_medium_temp
    Cp = R*(a8+(a9*Temp)+(a10*(Temp).^2)+(a11*(Temp).^3)+(a12*(Temp).^4));
else
    Cp = R*(a1+(a2*Temp)+(a3*(Temp).^2)+(a4*(Temp).^3)+(a5*(Temp).^4));
end
    


Enthalpy:-
Enthalpy is the measurement of energy in a thermodynamic system. The quantity of enthalpy equals the total content of heat of a system, equivalent to the system’s internal energy plus the product of volume and pressure.

function code for enthalpy:

% =================== Enthalpy function file ======================= %
function H = enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,Temp,R,global_medium_temp)

% comparing with global temperature %
% If less than global temp, min coefficents are taken %
% otherwise max coefficients are taken %

if Temp <= global_medium_temp
    H = R*Temp.*(a8+((a9*Temp)/2)+(a10*((Temp).^2)/3)+(a11*((Temp).^3)/4)+(a12*((Temp).^4)/5)+(a13./Temp));
else
    H = R*Temp.*(a1+((a2.*Temp)/2)+(a3*((Temp).^2)/3)+(a4*((Temp).^3)/4)+(a5*((Temp).^4)/5)+(a6./Temp));
end
 


Entropy:-
A thermodynamic quantity represents the unavailability of a system's thermal energy for conversion into mechanical work often interpreted as the degree of disorder or randomness in the system. The second law of thermodynamics says that entropy always increases with time.

function codr for entropy:

% ===================== Entropy function file ====================== %
function S = entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,Temp,R,global_medium_temp)

% comparing with global temperature %
% If less than global temp, min coefficents are taken %
% otherwise max coefficients are taken %

if Temp <= global_medium_temp
    S = R*(a8.*log(Temp)+(a9*Temp)+(a10*((Temp).^2)/2)+(a11*((Temp).^3)/3)+(a12*((Temp).^4)/4)+a14);
else
    S = R*(a1.*log(Temp)+(a2.*Temp)+(a3*((Temp).^2)/2)+(a4*((Temp).^3)/3)+(a5*((Temp).^4)/4)+a7);
end


Molecular Weight:-
molecular weight, also called molecular mass, is the mass of a molecule of a substance, based on 12 as the atomic weight of carbon-12. It is calculated in practice by summing the atomic weights of the atoms making up the substance’s molecular formula.

function code for molecular weight:

% Molecular weight function file %
function molecularweight = molecular_weight(Species_SP);
% standard atomic names - hydrogen,carbon,oxygen,nitrogen,argon
atomic_name = ['H','C','O','N','A'];

% standard atomic numbers
atomic_weight = [1.008 12.011 15.999 14.007 39.948];

% starting range 
molecularweight = 0;

% adding atomic number to get particular value 
for i = 1:length(Species_SP);
    for j = 1:length(atomic_name);
        if strcmpi(Species_SP(i),atomic_name(j));
            molecularweight = molecularweight + atomic_weight(j);
            position = j;
        end
    end
    n = str2num(Species_SP(i));
    
    if n>1
        molecularweight = molecularweight + atomic_weight(position)*(n-1);
    end
end
%  fprintf('molecular weight of %s = %d',Species_SP,molecularweight);
%  disp(' ');

end

EXPLANATION OF FUNCTION FILE FOR MOLECULAR WEIGHT:

Above Program is for calculating the molecular weight of the elements and their species.

 First, we created the function and name the file as molecularweight, and also gave the same name to the function.

 Then we get the standard elements from the data file then we defined it as a variable called Atomic_name in the array.

 Then we find out the atomic weight of those elements and then and assign it a variable called Atomic_weight.

 Then we wanted to find out the molecular weight of all the species so that we don't need to put it every time so we created such a code so that we can automate the process of calculating the molecular weight.

 For that we used the for-loop before for loop we assigned the value of molecular weight which is molecular weight as MW=0 this is our initial condition. Then for the for loop, we assigned the value of species to the variable i which will go running from 1 to length of Species and the total number of species is 53 so it will run 53 times.

 Then we assigned the variable j to the Atomic_name where our standard atomic elements are presented. And they are 5 in numbers so it is going to run for 5 times in the loop and going to give the required value.

 For starter, we used the if command and compare our both j and i string by using strcmp it is a string comparing command with the help of that command we compared species(i) and Atomic_name(j) so this was our condition.

 Then we used the molecular weight formula which is molecularweight+Atomic_weight(j) so if we run the system now then it will hope to compare with the species 1 and also with the Atomic_name if its value matches then it will going to select the number and will going to find the other number but if its value does not matches then for the species(1) it will going to run for 3 times and will going to match the element each time.

 But this was just for getting the individual elements and combining them but we have numbers in the species and there were no numbers in the Automic_name array so we defined a variable for j as p which is for the position.

 Then we converted the species from string to numbers and assigned it the variable N and then again use the if command and condition was if N>1 then it will going to use the formula MW = MW+Atomic_weight(positin)*(N-1); This formula is going to print the variable for n number of times if, for example, we have species called CH3 so we get the value of c and also H for only one time but we want it again 2 times so for that our formula comes into action and what it will do it will just add the molecular weight of H and get the atomic weight of it and going to repeat it for 3 times and going to get the value as 3 and then it will get subtract one time and will show the CH3. And finally, it will go to an end.

 This program is going to get all the species names and will calculate all species molecular weight.

MAIN CODE:

clear all
close all
clc
% reading given thermo.dat file %
f1 = fopen('THERMO.dat','r');
gt = fgetl(f1);
% reading header information %
global_temp = fgetl(f1);
% splitting the temperature data %
t = strsplit(global_temp);
% converting strings to numbers %
global_low_temp = str2num(t{2});
global_medium_temp = str2num(t{3});
global_high_temp = str2num(t{4});
% temperature range is created b/w low and high temperatures %
Temp = linspace(global_low_temp,global_high_temp,1000);
% universal gas constant
R = 8.314;
% skipping first three lines in the given file which contain comments %
for i = 1:3
 gt = fgetl(f1);
end
% reading the 53 species %
for i = 1:53
 gt = fgetl(f1);
 A = strsplit(gt,' ');
 Species_SP = (A{1});
 
 line_1 = fgetl(f1);
 a = strfind(line_1,'E');
 
 % extracting first 7 high temperature coefficents %
 a1 = line_1(1:a(1)+3);
 a1 = str2num(a1);
 a2 = line_1(a(1)+4:a(2)+3);
 a2 = str2num(a2);
 a3 = line_1(a(2)+4:a(3)+3);
 a3 = str2num(a3);
 a4 = line_1(a(3)+4:a(4)+3);
 a4 = str2num(a4);
 a5 = line_1(a(4)+4:a(5)+3);
 a5 = str2num(a5);
 
 line_2 = fgetl(f1);
 b = strfind(line_1,'E');
 a6 = line_2(1:b(1)+3);
 a6 = str2num(a6);
 a7 = line_2(a(1)+4:b(2)+3);
 a7 = str2num(a7);
 
 % extracting second 7 low temperature coefficents %
 a8 = line_2(b(2)+4:b(3)+3);
 a8 = str2num(a8);
 a9 = line_2(b(3)+4:b(4)+3);
 a9 = str2num(a9);
 a10 = line_2(b(4)+4:b(5)+3);
 a10 = str2num(a10);
 
 line_3 = fgetl(f1);
 c = strfind(line_1,'E');
 
 a11 = line_3(1:c(1)+3);
 a11 = str2num(a11);
 a12 = line_3(c(1)+4:c(2)+3);
 a12 = str2num(a12);
 a13 = line_3(c(2)+4:c(3)+3);
 a13 = str2num(a13);
 a14 = line_3(c(3)+4:c(4)+3);
 a14 = str2num(a14);
 
 % Specific heat %
 Cp = specificheat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,Temp,R,global_medium_temp);
 
 % Enthalpy %
 H = enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,Temp,R,global_medium_temp); 

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

 % Molecular weight
 molecularweight = molecular_weight(Species_SP);
 
%  read & write molecular weights 
 molecular_list = fopen('molecular weight.txt','a+');
 fprintf(molecular_list,'molecular weight of %s is %d',Species_SP,molecularweight);

 fclose(molecular_list);
 
 % plotting (Cp,H,S,Temp,Species_SP) and saving as folder %
 cur_dir = pwd;
 mkdir(Species_SP);
 cd(Species_SP);
 
 % plot 1 - Cp vs Temp %
 figure(1)
 plot(Temp,Cp,'linewidth',2,'color','r');
 xlabel('Temperature');
 ylabel('Specific Heat');
 title(sprintf('Variation of Specific Heat with Temperature Range -%s',Species_SP));

 grid on
 saveas(figure(1),'specific heat.jpg');
 
 % plot 2 - E vs Temp %
 figure(2)
 plot(Temp,H,'linewidth',2,'color','b');
 xlabel('Temperature')
 ylabel('Enthalpy')
 title(sprintf('Variation of Enthalpy with Temperature Range -%s',Species_SP));
%s',Species_SP));
 grid on
 saveas(figure(2),'enthalpy.jpg')
 
 % plot 3 - S vs Temp %
 figure(3)
 plot(Temp,S,'linewidth',2,'color','g');
 xlabel('Temperature')
 ylabel('Entropy')
 title(sprintf('Variation of Entropy with Temperature Range - %s',Species_SP));
 grid on
 saveas(figure(3),'entropy.jpg')
 
 cd(cur_dir);
 
 % molecuar weight %
 M = molecular_weight(Species_SP);
 f2 = fopen('molecular mass','w');
 sprintf('Molecular weight of %s is %f',Species_SP,M)

end

CODE EXPLANATION:

1. The code was started with the basic commands clear all, close all and clc.
2. Open the given THERMO.dat file using the 'fopen' command and use 'r' to read the file.
3. 'fgetl' command is used to read the line in the text file.
4. To get the global temperature ranges from the given text file, use the 'strsplit' command to split the temperature into a particular string.
5. Convert the obtained string cell into the numbers using 'str2num' command.
6. A temperature range cell is created between the global low-temperature and global high-temperature using 'linspace' command.
7. Universal gas constant (R) value input given as 8.314.
8. 'For loop' is used to skip the first three lines in the given text file which containing the comments.
9. Again 'for loop' is used to read all 53 species contained in the text file.
10.To find the strings, 'strfind' command is used.
11.All the 14 coefficients are found using the 'strfind' and they are converted into the numbers using 'str2num'.
12.Saved function files of Specific heat, Enthalpy, Entropy and Molecular weight are called to find them.
13.Creating directory using 'mkdir' command.
14.Plots are created between temperature and specific heat, enthalpy and entropy, and they are saved in the folder with their species name in jpg format.
15.Molecular weight is displayed on the command window along with species name using 'fprintf' command. used during the coding

SYNTAX USED:

  • pwd
  • 'fopen'
  • 'r'
  • 'fgetl'
  • 'fgetl'
  • 'str2num'
  • 'linspace'
  • 'strfind'
  • 'for loop'
  • 'strssplit'
  • 'fprintf'

MOLECULAR WEIGHT OUTPUT:

ans =

    'Molecular weight of O is 15.999000'


ans =

    'Molecular weight of O2 is 31.998000'


ans =

    'Molecular weight of H is 1.008000'


ans =

    'Molecular weight of H2 is 2.016000'


ans =

    'Molecular weight of OH is 17.007000'


ans =

    'Molecular weight of H2O is 18.015000'


ans =

    'Molecular weight of HO2 is 33.006000'


ans =

    'Molecular weight of H2O2 is 34.014000'


ans =

    'Molecular weight of C is 12.011000'


ans =

    'Molecular weight of CH is 13.019000'


ans =

    'Molecular weight of CH2 is 14.027000'


ans =

    'Molecular weight of CH2(S) is 14.027000'


ans =

    'Molecular weight of CH3 is 15.035000'


ans =

    'Molecular weight of CH4 is 16.043000'


ans =

    'Molecular weight of CO is 28.010000'


ans =

    'Molecular weight of CO2 is 44.009000'


ans =

    'Molecular weight of HCO is 29.018000'


ans =

    'Molecular weight of CH2O is 30.026000'


ans =

    'Molecular weight of CH2OH is 31.034000'


ans =

    'Molecular weight of CH3O is 31.034000'


ans =

    'Molecular weight of CH3OH is 32.042000'


ans =

    'Molecular weight of C2H is 25.030000'


ans =

    'Molecular weight of C2H2 is 26.038000'


ans =

    'Molecular weight of C2H3 is 27.046000'


ans =

    'Molecular weight of C2H4 is 28.054000'


ans =

    'Molecular weight of C2H5 is 29.062000'


ans =

    'Molecular weight of C2H6 is 30.070000'


ans =

    'Molecular weight of CH2CO is 42.037000'


ans =

    'Molecular weight of HCCO is 41.029000'


ans =

    'Molecular weight of HCCOH is 42.037000'


ans =

    'Molecular weight of H2CN is 28.034000'


ans =

    'Molecular weight of HCN is 27.026000'


ans =

    'Molecular weight of HNO is 31.014000'


ans =

    'Molecular weight of N is 14.007000'


ans =

    'Molecular weight of NNH is 29.022000'


ans =

    'Molecular weight of N2O is 44.013000'


ans =

    'Molecular weight of NH is 15.015000'


ans =

    'Molecular weight of NH2 is 16.023000'


ans =

    'Molecular weight of NH3 is 17.031000'


ans =

    'Molecular weight of NO is 30.006000'


ans =

    'Molecular weight of NO2 is 46.005000'


ans =

    'Molecular weight of HCNO is 43.025000'


ans =

    'Molecular weight of HOCN is 43.025000'


ans =

    'Molecular weight of HNCO is 43.025000'


ans =

    'Molecular weight of NCO is 42.017000'


ans =

    'Molecular weight of CN is 26.018000'


ans =

    'Molecular weight of HCNN is 41.033000'


ans =

    'Molecular weight of N2 is 28.014000'


ans =

    'Molecular weight of AR is 39.948000'


ans =

    'Molecular weight of C3H8 is 44.097000'


ans =

    'Molecular weight of C3H7 is 43.089000'


ans =

    'Molecular weight of CH3CHO is 44.053000'


ans =

    'Molecular weight of CH2CHO is 43.045000'

PLOTS FOR O2 ,N2 , CO2:

PLOT FOR O2:

a) temperature vs enthalpy

b) temperature vs entropy:

c) temperature vs specificheat:

PLOT FOR N2:

A) temperature vs enthalpy:

b) temperature vs entropy:

c) temperature vs specificheat:

PLOT FOR CO2:

A) TEMPERATURE VS ENTHALPY

B) TEMPERATURE VS ENTROPY:

C) TEMPERATURE VS SPECIFICHEAT:

We have plotted all the plots of all the 53 species but We have uploaded some of them for reference. In the plot, we can see that the values are changing concerning temperature.

By applying the formula in the code we get all the values of Specific heat, enthalpy, and entropy concerning temperature. As for both the high and low temperatures.

R is the universal gas constant which is 8.314

SCREENSHOT OF THE WINDOW WHERE ALL THE FOLDERS ARE LOCATED :

CONCLUSION:

With the help of MATLAB, we calculated the specific heat, enthalpy, Entropy, and molecular weight of all species and made the code in
such a way that it will be able to calculate for every species and plotted the data on the curve.
And learn how to make the file and insert the data in the file with Matlab code.

GOOGLE DRIVE LINK:

https://drive.google.com/drive/folders/1ssCCoY7c1bTRNzDbGklY3dL2eAONef5T?usp=sharing

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.