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.
- To Calculate the molecular weight of each species and to display it in the command window.
- To Plot the Cp, Enthalpy and Entropy for the local temperature range specific for each species.
- To Save the plots as images with appropriate names and dump them in separate folders for each species with suitable name.
PROCEDURE:
CODE:
close all
clear all
clc
%intial values
R =8.314;
%reading the file
data = fopen(\'THERMO.dat\',\'r\')
%header information
line1 = fgetl(data);
line2 = fgetl(data);
line3 = fgetl(data);
line4 = fgetl(data);
line5 = fgetl(data);
A = strsplit(line2,\' \');
Global_low_temperature = str2num(A{2})
Global_mid_temperature = str2num(A{3})
Global_high_temperature = str2num(A{4})
%species information
for i=1:53
Species= fgetl(data);
species_name = Species(1:6)
first_line = fgetl(data);
E1 = strfind(first_line,\'E\');
coefficient1 = first_line(1:(E1(1)+3));
coefficient1 = str2num(coefficient1);
coefficient2 = first_line((E1(1)+4):(E1(2)+3));
coefficient2 = str2num(coefficient2);
coefficient3 = first_line((E1(2)+4):(E1(3)+3));
coefficient3 = str2num(coefficient3);
coefficient4 = first_line((E1(3)+4):(E1(4)+3));
coefficient4 = str2num(coefficient4);
coefficient5 = first_line((E1(4)+4):(E1(5)+3));
coefficient5 = str2num(coefficient5);
second_line = fgetl(data);
E2 = strfind(second_line,\'E\');
coefficient6 = second_line(1:(E2(1)+3));
coefficient6 = str2num(coefficient6);
coefficient7 = second_line((E2(1)+4):(E2(2)+3));
coefficient7 = str2num(coefficient7);
coefficient8 = second_line((E2(2)+4):(E2(3)+3));
coefficient8 = str2num(coefficient8);
coefficient9 = second_line((E2(3)+4):(E2(4)+3));
coefficient9 = str2num(coefficient9);
coefficient10 = second_line((E2(4)+4):(E2(5)+3));
coefficient10 = str2num(coefficient10);
third_line = fgetl(data);
E3 = strfind(third_line,\'E\');
coefficient11 = third_line(1:(E3(1)+3));
coefficient11 = str2num(coefficient11);
coefficient12 = third_line((E3(1)+4):(E3(2)+3));
coefficient12 = str2num(coefficient12);
coefficient13 = third_line((E3(2)+4):(E3(3)+3));
coefficient13 = str2num(coefficient13);
coefficient14 = third_line((E3(3)+4):(E3(4)+3));
coefficient14 = str2num(coefficient14);
%creating folder for each species
mkdir([\'C:\\Users\\Home\\Desktop\\matlab\\challenges\\nasa file parsing\\species\\\',species_name]);
Temp = linspace(1,5000,100);
%calculation of properties based on the lower and upper temperature
for k=1:length(Temp)
if Temp(k)<3000
Cp = R*(coefficient1+coefficient2.*Temp+coefficient3.*Temp.^2+coefficient4.*Temp.^3+coefficient5.*Temp.^4);
H = R.*Temp.*(coefficient1+(coefficient2.*Temp)./2+(coefficient3.*Temp.^2)/3+(coefficient4.*Temp.^3)/4+(coefficient5.*Temp.^4/5)+coefficient6./Temp);
S = R.*(coefficient1*log(Temp)+coefficient2.*Temp+(coefficient3.*Temp.^2)/2+(coefficient4.*Temp.^3)/3+(coefficient5.*Temp.^4)/4+coefficient7);
else
Cp = R.*(coefficient8+coefficient9.*Temp+coefficient10.*Temp.^2+coefficient11.*Temp.^3+coefficient12.*Temp.^4);
H = R.*Temp.*(coefficient8+(coefficient9.*Temp)/2+(coefficient10.*Temp.^2)/3+(coefficient11.*Temp.^3)/4+(coefficient12.*Temp.^4/5)+coefficient13./Temp);
S = R.*(coefficient8.*log(Temp)+coefficient9.*Temp+(coefficient10.*Temp.^2)/2+(coefficient11.*Temp.^3)/3+(coefficient12.*Temp.^4)/4+coefficient14);
end
end
%plotting
%temperature vs entropy
figure(1);
plot(Temp,S,\'linewidth\',3,\'color\',\'r\');
title(species_name);
xlabel(\'Temperature\');
ylabel(\'Entropy\');
%temperature vs specific heat
figure(2);
plot(Temp,Cp,\'linewidth\',3,\'color\',\'g\');
title(species_name);
xlabel(\'Temperature\');
ylabel(\'Specific Heat\');
%Temperature vs Enthalpy
figure(3);
plot(Temp,H,\'linewidth\',3,\'color\',\'b\');
title(species_name);
xlabel(\'Temperature\');
ylabel(\'Enthalpy\');
%changing the directory to save the plots
address = ([\'C:\\Users\\Home\\Desktop\\matlab\\challenges\\nasa file parsing\\species\\\',species_name]);
cd(address);
%saving the plots to its respective species
saveas(1,\'Entropy.jpg\');
saveas(2,\'specific_heat.jpg\');
saveas(3,\'Enthalpy.jpg\');
%changing the directory for calculating the molecular_weight
cd(\'C:\\Users\\Home\\Desktop\\matlab\\challenges\\nasa file parsing\\\');
Molecular_weight = molecular_mass(species_name)
end
fclose(data)
function[mass] = molecular_mass(species_name)
species = upper(species_name);
atoms = \'HCNO\';
Atomic_Weight = [1.00784 12.0096 14.00643 15.99903];
mass = 0;
for i = 1:length(species)
for j = 1:4
if strcmp(species(i), atoms(j))
mass = mass +Atomic_Weight(j);
end
end
end
for q = 1:length(species)
n = str2num(species(q));
if n>1
p = (q-1);
mass = mass + Atomic_Weight(p)*(n-1);
end
end
for l = 1:length(species)
if strcmp(species(l),\'R\')
mass = mass + 38.8;
end
end
endStep 1: commands \"close all\" ; \"clear all\" ; \"clc\" these commands close all the previous files,clear all the previous contents and clears the content in command window respectivly.
Step 2 : commands \"R =8.314;\" this command assigns the value for the variable \'R\'.
Step 3 : Commands \"data = fopen(\'THERMO.dat\',\'r\')\" this command opens the file \'thermo.dat\' and reads value from it and stores it in the variable \'data\'.
Step 4 : Commands \"line1 = fgetl(data);\"
\"line2 = fgetl(data);\"
\"line3 = fgetl(data);\"
\"line4 = fgetl(data);\"
\"line5 = fgetl(data);\"
\"A = strsplit(line2,\' \');\"
\"Global_low_temperature = str2num(A{2})\"
\"Global_mid_temperature = str2num(A{3})\"
\"Global_high_temperature = str2num(A{4})\"
In the above commands ,the command \'fgetl(data)\' reads a line in the file and stores the value in the variable \'line1\'. The command \'A = strsplit(line2,\' \');\' splits the \'line2\' into segments whenever a \' \' comes. The values of the \'A\' are converted into numerical entity and stored in variables \'Global_low_temperature\' ; \'Global_mid_temperature\' ; \'Global_high_temperature\'.
Step 5 : commands \"for i=1:53
\"Species= fgetl(data);\"
\"species_name = Species(1:6)\"
The above commands makes a loop for 53 lines and reads next lines each time it runs .The command \"species_name = Species(1:6)\" stores the first six character in variable \'species_name\' as species name.
Step 6 : commands \"first_line = fgetl(data);\"
\"E1 = strfind(first_line,\'E\');\"
\"coefficient1 = first_line(1:(E1(1)+3));\"
\"coefficient1 = str2num(coefficient1);\"
\"coefficient2 = first_line((E1(1)+4):(E1(2)+3));\"
\"coefficient2 = str2num(coefficient2);\"
\"coefficient3 = first_line((E1(2)+4):(E1(3)+3));\"
\"coefficient3 = str2num(coefficient3);\"
\"coefficient4 = first_line((E1(3)+4):(E1(4)+3));\"
\"coefficient4 = str2num(coefficient4);\"
\"coefficient5 = first_line((E1(4)+4):(E1(5)+3));\"
\"coefficient5 = str2num(coefficient5);\"
The set of commands reads first line of the species and stores in the variable \'first_line\' and then searches and find for the \'E\' in the line and stores the position of \'E\' in the variable \'E1\' and then calculates the coefficients of the species.
Step 7 : commands \"second_line = fgetl(data);\"
\"E2 = strfind(second_line,\'E\');\"
\"coefficient6 = second_line(1:(E2(1)+3));\"
\"coefficient6 = str2num(coefficient6);\"
\"coefficient7 = second_line((E2(1)+4):(E2(2)+3));\"
\"coefficient7 = str2num(coefficient7);\"
\"coefficient8 = second_line((E2(2)+4):(E2(3)+3));\"
\"coefficient8 = str2num(coefficient8);\"
\"coefficient9 = second_line((E2(3)+4):(E2(4)+3));\"
\"coefficient9 = str2num(coefficient9);\"
\"coefficient10 = second_line((E2(4)+4):(E2(5)+3));\"
\"coefficient10 = str2num(coefficient10);\"
The above set of commands is the same as Step 6.Expect that it reads the second line of the species and calculates the coefficients.
Step 8 : commands \"third_line = fgetl(data);\"
\"E3 = strfind(third_line,\'E\');\"
\"coefficient11 = third_line(1:(E3(1)+3));\"
\"coefficient11 = str2num(coefficient11);\"
\"coefficient12 = third_line((E3(1)+4):(E3(2)+3));\"
\"coefficient12 = str2num(coefficient12);\"
\"coefficient13 = third_line((E3(2)+4):(E3(3)+3));\"
\"coefficient13 = str2num(coefficient13);\"
\"coefficient14 = third_line((E3(3)+4):(E3(4)+3));\"
\"coefficient14 = str2num(coefficient14);\"
The above set of commands is the same as Step 6.Expect that it reads the third line of the species and calculates the coefficients.
Step 9 : commands \"mkdir([\'C:\\Users\\Home\\Desktop\\matlab\\challenges\\nasa file parsing\\species\\\',species_name]);\" .The command creates folder in the path mentioned with \'species_name\' as the folder name.
Step 10 : Command\"Temp = linspace(1,5000,100);\"
\"for k=1:length(Temp)\"
\"if Temp(k)<3000\"
\"Cp = R*(coefficient1+coefficient2.*Temp+coefficient3.*Temp.^2+coefficient4.*Temp.^3+coefficient5.*Temp.^4);\"
\"H = R.*Temp.*(coefficient1+(coefficient2.*Temp)./2+(coefficient3.*Temp.^2)/3+(coefficient4.*Temp.^3)/4+(coefficient5.*Temp.^4/5)+coefficient6./Temp);\"
\"S = R.*(coefficient1*log(Temp)+coefficient2.*Temp+(coefficient3.*Temp.^2)/2+(coefficient4.*Temp.^3)/3+(coefficient5.*Temp.^4)/4+coefficient7);\"
\"else\"
\"Cp = R.*(coefficient8+coefficient9.*Temp+coefficient10.*Temp.^2+coefficient11.*Temp.^3+coefficient12.*Temp.^4);\"
\"H = R.*Temp.*(coefficient8+(coefficient9.*Temp)/2+(coefficient10.*Temp.^2)/3+(coefficient11.*Temp.^3)/4+(coefficient12.*Temp.^4/5)+coefficient13./Temp);\"
\"S = R.*(coefficient8.*log(Temp)+coefficient9.*Temp+(coefficient10.*Temp.^2)/2+(coefficient11.*Temp.^3)/3+(coefficient12.*Temp.^4)/4+coefficient14);\"
\"end\"
\"end\"
The above commands first assings value for temperature between \'1-5000\' with equal spacing.Then a loop is created for calculating the properties of species for each temperature value assinged and we assume that \'2500\' to be the transient temperature(i.e below it is lower temperature and above it as higher temperature). If the temperature is a low temperature then, first 7 coefficients of the species is used to calculate its properties and if the temperature is a high temperature then, last 7 coefficients is used to calculate the properties.The properties are calculated based on the formulas:
Step 11 : Commands \" figure(1);\"
\"plot(Temp,S,\'linewidth\',3,\'color\',\'r\');\"
\"title(species_name);\"
\"xlabel(\'Temperature\');\"
\"ylabel(\'Entropy\');\"
\"figure(2);
\"plot(Temp,Cp,\'linewidth\',3,\'color\',\'g\');\"
\"title(species_name);\"
\"xlabel(\'Temperature\');\"
\"ylabel(\'Specific Heat\');\"
\"figure(3);\"
\"plot(Temp,H,\'linewidth\',3,\'color\',\'b\');\"
\"title(species_name);\"
\"xlabel(\'Temperature\');\"
\"ylabel(\'Enthalpy\');\"
The commands create three plots namely \'Temperature vs Entropy\' ; \'Temperature vs specific heat\' and \'Temperature vs Enthalpy\' for the species with the title as species name.
Step 12 : commands \"address = ([\'C:\\Users\\Home\\Desktop\\matlab\\challenges\\nasa file parsing\\species\\\',species_name]);
cd(address);\"
\"saveas(1,\'Entropy.jpg\');
\"saveas(2,\'specific_heat.jpg\');
\"saveas(3,\'Enthalpy.jpg\');
The commands changes the directory path to the species folder created .so that the plots of the species are stored in their respective folders.And the \'saveas\' command saves the plot.
Step 13 : Commands \"cd(\'C:\\Users\\Home\\Desktop\\matlab\\challenges\\nasa file parsing\\\');\"
\"Molecular_weight = molecular_mass(species_name)\"
We again change the directory where main program is stored as to execute the function to calculate the molecular weight of the sepcies.Then the function \'molecular_mass(species_name)\' is called as to calculate the molecular_weight of the species.
Step 14 : commands \"species = upper(species_name);\"
\"atoms = \'HCNO\';\"
\"Atomic_Weight = [1.00784 12.0096 14.00643 15.99903];\"
\"mass = 0;\"
The above commands first converts the character in the \'species_name\' to uppercase then atoms and their atomic weights are initialised and the mass is initialised to zero.
Step 15 : commands \"for i = 1:length(species)\"
\"for j = 1:4\"
\"if strcmp(species(i), atoms(j))\"
\"mass = mass +Atomic_Weight(j);\"
\"end\"
\"end\"
\"end\"
\"for q = 1:length(species)\"
\"n = str2num(species(q));\"
\"if n>1\"
\"p = (q-1);\"
\"mass = mass + Atomic_Weight(p)*(n-1);\"
\"end\"
\"end\"
\"for l = 1:length(species)\"
\"if strcmp(species(l),\'R\')\"
\"mass = mass + 38.8;\"
\"end\"
\"end\"
\"end\"
In the above commands we use loop to compare the character that represents the atoms and add their atomic weights to calculate the molecular weight of the species and in the end we used a special loop for argon as it involves character \'R\'.
OUTPUT:
RESULT :
Thus NASA thermodynamic data is parsed and the plots for the properties of each species are created and the molecular weights of the each species is calculated by adding the atomic weight of the atoms present in the species.
Attachements:
https://drive.google.com/open?id=1ALLSTf_LAgGJiS4x9eXr3fLcPnbtDny5
https://drive.google.com/open?id=1WkWXDkMw5JCJo4kHwgopGKw_pRX39hxC
https://drive.google.com/open?id=1c171cyuEFMYJ_wQOv_l1COC2hQ4YmeEI
https://drive.google.com/open?id=1bYJdPrMc9f2KD_skF8eHEWsiNbTsGBfO
https://drive.google.com/open?id=1Rok0HstLJlH9aS_yrn8kGzgGmzQKcDod
https://drive.google.com/open?id=1RjSpp5kY7yVoQR1LOXl_hLOdWqEcuA_e
https://drive.google.com/open?id=1SNXsh9FaJimvj_wh82CN1ROaaVVQoU8z
https://drive.google.com/open?id=1So4v8WuBMqjsPgn_pPuqExaokzymJvld
https://drive.google.com/open?id=1_duqShwi00U1_CZVKp53OsEH7xJqJlPd
https://drive.google.com/open?id=1nzGe5HBe4SdCNh7NGgvuPr5KXXKGt37M
https://drive.google.com/open?id=1i2Yd_EqQtNDLtoZE8iFPIlC6IXdGSo5K
https://drive.google.com/open?id=1MsxSnuqNrK_ZC0VA-ui4x7jZkPsJRbOO
https://drive.google.com/open?id=11SpClJ7pZP1SGhin6Kj1_Khky0ru7KZF
https://drive.google.com/open?id=1lclQu2Vh1_-jHnJAOX_vHICs1vLcar0n
Referrence:
https://www.mathworks.com/help/matlab/ref/mkdir.html