Parsing NASA thermodynamic data using MATLAB

NOTE- ALL NECESSARY SCREENSHOTS AND PLOTS ARE UPLOADED TO GOOGLE DRIVE.

The objective for this program is to write a function a function that extracts the 14 co-efficients and calculates the enthalpy, entropy and specific heats for all the species in the data file.

NASA came up with polynomials that can be used to evaluate thermodynamic properties such as Cp, H and S. They have also documented the co-efficients that are required to evaluate these polynomials for 1000+ species.

Steps:

1. Started with fopen command to open th provided file i.e \'THERMO.DAT.txt\' file.

2.Initially,fgetl command was used to read the first lines of the file. Later using fscan global low,mid and high temp values are read.

3.To find exponent \'E\', strfind command is used to find value of \'E\' in line to separate coefficients.

4. New directory is made using mkdir command to create folder with species name.

5. To calculate specific heat , enthalpy and entropy respective functions are defined.

6. Later, three graphs are lot viz. T vs Cp, T vs entropy and T vs enthalpy

7. Atlast using cd directory is changed and the molecular mass is calculated using function which is user defined.

The program is as follows;-

MAIN CODE

clear all
close all
clc

% universal gas constant
R = 8.314 ;          

%Reading header
f1= fopen(\'THERMO.dat.txt\',\'r\');
Header = fgetl(f1);

%global temperatures
A = fscanf(f1,\'%f \\n\',[1 3]);
global_low_temp = A(1);
global_mid_temp = A(2);
global_high_temp = A(3);

% skipping 3 lines
line1= fgetl(f1);
line2= fgetl(f1);
line3= fgetl(f1);

for i= 1:53
name_line = fgetl(f1);
%extracting the name of species
species = name_line(1:6);     
% Reading coefficients
tline1 = fgetl(f1);
% Extracting coefficients
y1 = strfind(tline1,\'E\');  
a1 = tline1(1:y1(1)+3);
a1 = str2double(a1);
a2 = tline1((y1(1)+4: (y1(2)+3)));
a2 = str2double(a2);
a3 = tline1(y1(2)+4: (y1(3)+3));
a3= str2double(a3);
a4= tline1 (y1(3)+4: (y1(4)+3));
a4= str2double(a4);
a5= tline1(y1(4)+4: (y1(5)+3));
a5= str2double(a5);
tline2= fgetl(f1);
y2 = strfind(tline2,\'E\');
a6 = tline2(1:(y2(1)+3));
a6 = str2double(a6);
a7 = tline2((y2(1)+4): (y2(2)+3));
a7= str2double(a7);
a8= tline2((y2(2)+4):(y2(3)+3));
a8=str2double(a8);
a9= tline2((y2(3)+4):(y2(4)+3));
a9=str2double(a9);
a10= tline2((y2(4)+4): y2(5)+3);
a10= str2double(a10);
tline3 = fgetl(f1);
y3 = strfind(tline3,\'E\');
a11 = tline3(1: (y3(1)+3));
a11 = str2double(a11);
a12 = tline3((y3(1)+4):(y3(2)+3));
a12=str2double(a12);
a13=tline3((y3(2)+4):(y3(3)+3));
a13=str2double(a13);
a14=tline3((y3(3)+4):(y3(4)+3));
a14=str2double(a14);

% Creating folders with species name
subfolder_name = species;
mkdir([\'C:/Users/Abhishek/Documents/project1\',subfolder_name]);

T = linspace(1,5000,50);
%Function for cp
Cp=specific_heat(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T,R,global_low_temp,global_mid_temp,global_high_temp);

%Function for enthalpy
H=enthalpy(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T,R,global_low_temp,global_mid_temp,global_high_temp);

%Function for entropy
S=entropy(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T,R,global_low_temp,global_mid_temp,global_high_temp);

%graphs

figure(1) %T vs specific heat
plot(T,Cp)
title(subfolder_name)
xlabel(\'Temp\')
ylabel(\'Specific Heat\')

figure(2) %for T vs enthalpy
plot(T,H)
title(species)
xlabel(\'Temp\')
ylabel(\'Enthalpy\')

figure(3) %for T vs entropy
plot(T,S)
title(species)
xlabel(\'Temp\')
ylabel(\'Entropy\')


% Changing the directory
cd(\'C:/Users/Abhishek/Documents\')

% Function for molecular_mass
mm=molecular_mass(species)

% Changing directory to save the graphs in respective species folder
filepath = ([\'C:/Users/Abhishek/Documents/project1\',subfolder_name]);
cd (filepath)
saveas (1,\'T vs specific heat.jpg\')
saveas (2,\'T vs enthalpy.jpg\')
saveas (3,\'T vs entropy.jpg\')

% Saving the molecular mass in .txt file
f2 = fopen(\'Molecular mass.txt\',\'w\');
fprintf(f2, \'%s \\n %s %f \\n\',species,\'Molecular mass\',mm);

fclose(f2);

end

fclose(f1);

Function for enthalpy;-

function H=enthalpy(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T,R,global_low_temp,global_mid_temp,global_high_temp)
if T<=1000
    H=(a1.*T+a2.*(T.^2)./2+a3.*(T.^3)./3+a4.*(T.^4)./4+a5.*(T.^5)./5+a6)*R;
else
    H=(a8.*T+a9.*(T.^2)./2+a10.*(T.^3)./3+a11.*(T.^4)./4+a12.*(T.^5)./5+a13)*R;
end

 

Function for molecular mass;-

function[mm]=molecular_mass(species)

% A program for calculating molecular mass of each species
atoms = {\'H\',\'C\',\'N\',\'O\',\'AR\'};          
% atomic weights
aw = [1,12,14,16,40]; 
% initial value for molecular mass
mm=0;                          

for i=1:4
for j= 1:4
if (strcmp(species(i),atoms(j))==1);

mm= mm+aw(j);
end
end
end

% Adding the weights of the repeating atoms
for g = 1:length(species)
n=str2double(species(g));
if n>1
p=(g-1);
mm=mm+aw(p)*(n-1);
end
end
% Atomic weight of argon
for l=1:length(species)
if strcmp(species(l),\'R\')
mm=mm+38.8;
end
end

Function for specific heat;-

function Cp=specific_heat(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T,R,global_low_temp,global_mid_temp,global_high_temp)
if T<=1000
    Cp=(a1+a2.*(T)+a3.*(T.^2)+a4.*(T.^3)+a5.*(T.^4))*R;
else
    Cp=(a8+a9.*(T)+a10.*(T.^2)+a11.*(T.^3)+a12.*(T.^4))*R;
end

Function for entropy;-

function S=entropy(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T,R,global_low_temp,global_mid_temp,global_high_temp)
if T<=1000
    S=(a1.*log(T)+a2.*T+a3.*(T.^2)./2+a4.*(T.^3)./3+a5.*(T.^4)./4+a7)*R;
else
    S=(a8.*log(T)+a9.*T+a10.*(T.^2)./2+a11.*(T.^3)./3+a12.*(T.^4)./4+a14)*R;
end