Parsing NASA Thermodynamic Data
1. 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.
2. GOVERNING EQUATIONS IF ANY:
1. `Cp=R⋅(a_1+a_2T+a_3T^2+a_4T^3+a_5T^4)`
2.`H=R⋅T⋅(a_1+a_2T^2+a_3(T^2)/3+a_4(T^3)/4+a_5(T^4)/5+a_6/T)`
3. `S = R ⋅ ( a_1 ⋅ ln ( T ) + a_2 T + a_3 (T^ 2)/ 2 + a_4( T^ 3)/ 3 + a_5 (T^ 4)/ 4 + a_7 )`
3. OBJECTIVES OF THE PROJECT:
1. To calculate the molecular weight of each species and display it in the command window.
2. To plot the Cp, Enthalpy and Entropy for the local temperature range (low temperature : high temperature) specific for each species.
3. 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.
4. FLOW OF THE PROJECT:
1. Opening the data file and check for key things like unnecessary lines and line of start i.e from where actual data can be extracted.
2. Extracting global and local minimum, medium and maximum temperatures.
3. Looping to extract all the 14 coefficients line wise and assigning then to a variable.
4. Using governing equations, finding the values of Specific heat, Enthalpy and Entropy.
5. Plotting the graphs for `O_2,N_2 and CO_2`.
6. Saving all plots in individual folder.
7. Calculating molecular weight and displaying it in command window.
5.MATLAB CODE:
1.Main Code
clc
clear all
close all
%Opening and accessing the THERMO.dat file
nasa_data = fopen('THERMO.dat','r');
file_title = fgetl(nasa_data);
global_temp = str2num(fgetl(nasa_data));
T_low = global_temp(1);
T_medium = global_temp(2);
T_high = global_temp(3);
for i = 1:3
fgetl(nasa_data);
end
%end_line = 'END';
R = 8.314;
output = main_function(nasa_data,R)
Explanation:
- We opened the file named THERMO.dat giving only read permission.
- We have skipped lines consisting of the comments in the THERMO.dat file so that we can extract global temperatures i.e. Low, Medium and High from the second line.We assigned it to variables using string to number conversion.
- Created a variable named as = 'END' so as to determine the end of the dat file.
6. FUNCTIONS
1. main_function(nasa_data,R)
function[final_plots,molecular_weights] = main_function(nasa_data,R)
%Implementation of for loop for extracting coefficients of the compounds
for count = 1: 53
o = fgetl(nasa_data);
species1 = findstr(o,' ');
o_local_temp = findstr(o,'G');
o1 = str2num(o(48:75));
o_low_temp = o1(1);
o_high_temp = o1(2);
o_med_temp = o1(3);
lower_temp = linspace(o_low_temp,o_med_temp,100);
higher_temp = linspace(o_med_temp,o_high_temp,100);
o_val1 = fgetl(nasa_data);
%for first line coefficients, we extract them by finding the suitable
%character in the string, i.e. E
a = findstr(o_val1,'E');
bspecies=[o(1:species1(1)-1)]
a1 = str2num(o_val1((1:a(1)+3)));
a2 = str2num(o_val1(a(1)+4:a(2)+3));
a3 = str2num(o_val1(a(2)+4:a(3)+3));
a4 = str2num(o_val1(a(3)+4:a(4)+3));
a5 = str2num(o_val1(a(4)+4:a(5)+3));
%Extracting 2nd line coefficients
o_val2 = fgetl(nasa_data);
b = findstr(o_val2,'E');
a6 = str2num(o_val2(1:b(1)+3));
a7 = str2num(o_val2(b(1)+4:b(2)+3));
a8 = str2num(o_val2(b(2)+4:b(3)+3));
a9 = str2num(o_val2(b(3)+4:b(4)+3));
a10 = str2num(o_val2(b(4)+4:b(5)+3));
%Extracting third line coefficients
o_val3 = fgetl(nasa_data);
c = findstr(o_val2,'E');
a11 = str2num(o_val3(1:c(1)+3));
a12 = str2num(o_val3(c(1)+4:c(2)+3));
a13 = str2num(o_val3(c(2)+4:c(3)+3));
a14 = str2num(o_val3(c(3)+4:c(4)+3));
final_plots = plots(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,R,higher_temp,lower_temp,bspecies);
molecular_weights = MoleWt(bspecies)
end
Explanation:
- The main function accepts the nasa data from the file
- The function ouptuts the final plots of the specific heat, enthalpy and entropy curve versus the temperature.
- After this, we ran the for loop for 53 times, as the total number of the compounds in the database file. By using the fgetl command, the appropriate liine is selected from the file for extracting the data from that particular line.
- Thus, the values of local temperatures and thermal coefficients are thus extracted from the loop.
- The final line of the code basically calls the other function to calculate the main plots and finding the molecular weights, which will be explained below.
2. plots(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,R,higher_temp,lower_temp,bspecies)
function [Cp_high,Cp_low,H_high,H_low,S_high,S_low] = plots(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,R,higher_temp,lower_temp,bspecies)
%Specific Heat values of species
Cp_high = R*(a1+a2*higher_temp+a3*(higher_temp).^2+a4*(higher_temp).^3+a5*(higher_temp).^4);
Cp_low = R*(a8+a9*lower_temp+a10*(lower_temp).^2+a11*(lower_temp).^3+a12*(lower_temp).^4);
figure(1)
plot(lower_temp,Cp_low,'linewidth',2,'color','r')
hold on
plot(higher_temp,Cp_high,'linewidth',2,'color','b')
hold off
xlabel('Temperature')
ylabel('Specific heat')
%Creating a folder and saving the files
folder = 'D:\Matlab practice\NASA 2';
cd(folder)
cd('species')
mkdir(bspecies)
cd (bspecies)
saveas(gcf,'Specific Heat.png')
cd(folder);
%Enthalpy values of the species
H_high = R*(a1+a2.*higher_temp./2+a3.*(higher_temp).^2./3+a4.*(higher_temp).^3./4+a5.*(higher_temp).^4./5+a6./(higher_temp)).*(higher_temp);
H_low = R*(a8+a9.*lower_temp./2+a10.*(lower_temp).^2./3+a11.*(lower_temp).^3./4+a12.*(lower_temp).^4./5+a13./(lower_temp)).*(lower_temp);
figure(2)
plot(higher_temp,H_high,'linewidth',2,'color','r')
hold on
plot(lower_temp,H_low,'linewidth',2,'color','b')
hold off
xlabel('Temperature')
ylabel('Enthalpy')
cd('species')
cd(bspecies);
saveas(gcf,'Enthalpy.png')
cd(folder);
%Entropy values of the species
S_high = R*(a1.*log(higher_temp)+a2.*higher_temp+a3.*(higher_temp).^2./2+a4.*(higher_temp).^3./3+a5.*(higher_temp).^4./4+a7);
S_low = R*(a8.*log(lower_temp)+a9.*lower_temp+a10.*(lower_temp).^2./2+a11.*(lower_temp).^3./3+a12.*(lower_temp).^4./4+a14);
figure(3)
plot(higher_temp,S_high,'linewidth',2,'color','r')
hold on
plot(lower_temp,S_low,'linewidth',2,'color','b')
hold off
xlabel('Temperature')
ylabel('Entropy')
cd('species')
cd(bspecies);
saveas(gcf,'Entropy.png')
cd(folder);
end
Explanation:
- The local temperatures are divided into 2 temperature ranges by linspace() command. One range is from Low Temperature to Medium Temperature and second range is from Medium Temperature to High Temperature.
- The next portion of the code consists of the formulae to find out the values of the Specific Heat, Enthalpy and the Entropy.
- Next, a figure of the plots, i.e. the graphs of the various values of the Specific Heat, Enthalpy and the Entropy are plotted against the Temperature.
- The plots are needed to be saved into a directory. Hence the path where the plots are to be saved is specified. Then we use the makedir() command to make the directory and save the file. The plots are now saved as a .png file with saveas() command.
- Same process is repeated for plotting all the graphs.
3. molecularWt(bspecies)
function [molecular_wt] = MoleWt(bspecies)
%Initialization of the elements and its atomic weights
species = ['C','H','O','N','AR'];
mol_wts = [12.0096,1.00784,15.99903,14.00643,39.948];
%Calculation of the molecular weights
molecular_wt = 0;
for i = 1:length(bspecies) %2
for j = 1:5
if strcmp(bspecies(i),species(j))
molecular_wt = molecular_wt + mol_wts(j);
a=j;
end
end
n = str2num(bspecies(i));
if n>1
molecular_wt = molecular_wt+mol_wts(a)*(n-1);
end
end
molecular_wt
%Writing the values of molecular weights in a notepad file
f3 = fopen('Molecular Weights.txt','a+')
fprintf(f3,'Molecular weight of %s is %f \n',bspecies,molecular_wt)
fclose(f3)
end
Explanation:
- Since each compound is made of some parent element, we have pre-defined some parent elements such as H, C, O, N, Ar along with their molecular weights.Now, the initial value of molecular weight is assigned as zero.
- We now use the for loop and nested if loops to determine the chemical formula of the compound.
- Within the nested if loop, we use the strcmp() command, which basically compares the two values submitted in the strcmp syntax.
- Then, we initialize the variable n and convert it to number or integer.In teh next if loop, if the vlue of the variable is greater than 1, the actual molecular weight of the compound is the weight of the element plus the number of times the element is present in the compound.
- Thus, molecular weights of all the compounds are found are saved in a separate notepad file, by using the fopen command.
7. RESULTS
- Plot For `N_2`
2.Plot for `O_2`
3 Plot for `CO_2`
The screenshot of the folder consisting graphs of every compound in the THERMO.dat file
The screenshot of notepad file containing molecular weight of all the compounds
8.CONCLUSION:
Thus we have successfully parsed the thermodynamic data from a data file and plotted the graphs of different thermodynamic properties and calculated the molecular weights of different species present in the data file.