Parsing NASA Thermodynamic Data
clear all
close all
clc
%Universal Gas Constant Value
R = 8.314;
%Reading the data file
f1 = fopen(\'THERMO(1).dat\',\'r\');
fgetl(f1);
%Getting the global min, mid and max temperatures
L2 = fgetl(f1);
t = strsplit(L2);
global_min = str2double(t{2});
global_mid = str2double(t{3});
global_max = str2double(t{4});
%Defining the Temperature range
temp = linspace(global_min,global_max,500);
%Ignoring the header Information
for i = 1:3
fgetl(f1);
end
for i = 1:53
l = fgetl(f1);
element_name = strtok(l);
l1 = fgetl(f1);
l2 = fgetl(f1);
l3 = fgetl(f1);
%Finding the position of E
p1 = strfind(l1,\'E\');
p2 = strfind(l2,\'E\');
p3 = strfind(l3,\'E\');
%Finding High Temp Coefficient
hc1 = str2double(l1(1:p1(1)+3));
hc2 = str2double(l1(p1(1)+3:p1(2)+3));
hc3 = str2double(l1(p1(2)+4:p1(3)+3));
hc4 = str2double(l1(p1(3)+4:p1(4)+3));
hc5 = str2double(l1(p1(4)+4:p1(5)+3));
hc6 = str2double(l2(1:p2(1)+3));
hc7 = str2double(l2(p2(1)+4:p2(2)+3));
%Finding Low Temp Coefficient
lc1 = str2double(l2(p2(2)+4:p2(3)+3));
lc2 = str2double(l2(p2(3)+4:p2(4)+3));
lc3 = str2double(l2(p2(4)+4:p2(5)+3));
lc4 = str2double(l3(1:p3(1)+3));
lc5 = str2double(l3(p3(1)+4:p3(2)+3));
lc6 = str2double(l3(p3(2)+4:p3(3)+3));
lc7 = str2double(l3(p3(3)+4:p3(4)+3));
%Functions for Specific Heats, Enthalpy and Entropy
c_p = specific_heat(hc1,hc2,hc3,hc4,hc5,lc1,lc2,lc3,lc4,lc5,R,temp,global_mid);
H = enthalpy(hc1,hc2,hc3,hc4,hc5,hc6,lc1,lc2,lc3,lc4,lc5,lc6,R,temp,global_mid);
S = entropy(hc1,hc2,hc3,hc4,hc5,hc7,lc1,lc2,lc3,lc4,lc5,lc7,R,temp,global_mid);
plotting(c_p,H,S,temp,element_name);
%Finding the Molecular Weight of the species
mw(i) = molecular_weight(element_name);
X = sprintf(\'%s %d\',element_name,mw(i));
disp(X);
end
fclose(f1)
Function for Specific Heat:
function c_p = specific_heat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,R,temp,global_mid)
if temp<=global_mid
c_p = R*(a1 + a2.*temp + a3.*(temp.^2) + a4.*(temp.^3) + a5.*(temp.^4));
else
c_p = R*(b1 + b2.*temp + b3.*(temp.^2) + b4.*(temp.^3) + b5.*(temp.^4));
end
end
Function for Enthalpy:
function H = enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,R,temp,global_mid)
if temp<= global_mid
H = R.*temp.*(a1+a2.*(temp./2)+ a3.*((temp.^2)./3)+ a4.*((temp.^3)./4) + a5.*((temp.^4)./5)+a6./temp);
else
H = R.*temp.*(b1 + b2.*(temp./2) + b3.*((temp.^2)./3)+b4.*((temp.^3)./4)+b5.*((temp.^4)./5)+b6./temp);
end
end
Function for Entropy:
function S = entropy(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,R,temp,global_mid)
if temp<=global_mid
S = R*(a1.*log(temp)+a2.*temp+ a3.*(temp.^2)./2 + a4.*(temp.^3)./3 + a5.*(temp.^4)./4 + a7)
else
S = R*(b1.*log(temp) + b2.*temp + b3.*(temp.^2)./2 + b4.*(temp.^3)./3 +b5.*(temp.^4)./4 + b7);
end
end
Function for Molecular Weight:
function [weight] = molecular_weight(element_name)
s = {\'H\' \'C\' \'N\' \'O\' \'AR\'};
atomic_wt = [1 12 14 16 40];
weight = 0;
for i = 1:length(element_name)
for j = 1:length(s)
if strcmp(element_name(i),s(j))
weight = weight + atomic_wt(j);
pos = j;
end
end
n= str2double(element_name(i));
if n>1
weight = weight + atomic_wt(pos)*(n-1);
end
end
end
Screenshot of the Folders:
Plot for O2:
Plot for N2:
Plot for CO2:
Step-Wise explanation of the code:
- First of all we start with clear all, close all and clc to remove all the previous calculations on the command window and the workspace.
- We take R (Universal Gas Constant) as a given input.
- To read the given text file, ‘fopen’ command is used and is stored in the variable ‘f1’, which acts as a pointer to the file, which can be later used to extract the details from the file.
- Now f1 is used to read all the lines, which can be done by the function ‘fgetl’. Here the first line was redundant hence was let go, then the second line which contains the global temperature range was stored in the 3 different variables, ‘global_min’, ’global_mid’, ’global_max’. Then we use ‘linspace command’, where we define 500 different values in between min and max values.
- Then, 3 more header lines are let go.
- Then, comes the main lines which contains the temperature co-efficients for 53 species. To read each of the coefficient, and to perform calculations on them, the best way was to use the ‘for loop’ command. Inside the ‘for loop’, the first line which contains the element name was stored in a variable then, a function called ‘strtok’ was used which reads the first word of any line, which is basically the element name, and is stored in a variable called ‘element’.
- Then the next three lines contains the different co-efficients, which were first stored in 3 different variables. Now, the main task was to find the position of ‘E’ in all the temperature, and then using this position we needed to extract the coefficient. The position of ‘E’ was found by ‘strfind’ function and stored in another variable, then by using the variable the 14 coefficients were extracted and stored in 14 different variables, one for high temp and other for low temp. Since, everything is a string, the ‘str2double command’ is used to convert the co-efficients into a double character.
- Now 3 separate functions are written to calculate the specific heat, enthalpy and entropy. The co-efficients extracted were passed on as input parameters to these functions and these thermodynamic property were calculated for each species.
- Then these 3 thermodynamic properties along with the temperature range were passed onto the plotting function to plot the necessary graphs.
- Then another function has been created to calculate the molecular weight of each species and is displayed on the command window by the function ‘disp’.
- The file is then closed by the ‘fclose’ command.
Function for specific Heat:
The input parameters for this functions are first 5 low and high coefficients, temp, R, and global mid function. The global mid is used to define the condition whether the temp is less than or greater than it, and to choose the appropriate low and high co-efficient. The condition was checked by if-else command.
Enthapy Function:
The input parameter for this function were the first 6 low and high coefficents, temp, R and global mid function. The global mid is used to define the condition whether the temp is less than or greater than it, and to choose the appropriate low and high co-efficient. The condition was checked by if-else command.
Entropy Function:
The input parameters for his function were, the first 5 and 7th low and high coefficients, temp, R, and global mid function. The global mid is used to define the condition whether the temp is less than or greater than it, and to choose the appropriate low and high co-efficient. The condition was checked by if-else command.
Molecular Weight Function:
The function only takes the element name as the input parameter. First a string is created with the 5 basic element names, along with that an atomic weight function which contains their atomic weights. Then we used ‘strcmp’ command to compare the name of the element with the string array. If there is a repetition of an element, ‘str2double’ is use to find and (n-1) is multiplied to atomic weight of the element to calculate the molecular weight of the respective species.
Plotting function:
This function require the thermodynamic properties, temp, and element name as the input parameters.
- Firstly, ‘cd(pwd)’ command is used to be in the current directory.
- Then ‘mkdir’ is used to make a new folder inside the current directory with the name of the element. This folder is created to save all the plots for the particular species.
- Now 3 different plots are created that are ‘Enthalpy vs temp’, ‘Specific Heat vs Temp’, Entropy vs Temp’.
- Then these plots were then saved by using the ‘saveas command’.
- The cd command is then used to close the folder and to change the back of the original folder.
- As this function is running inside the loop so, 53 folders are created for each particular species.