Project 1 - Parsing NASA thermodynamic data
AIM:
To write a code in matlab to parse the NASA thermodynamic data file and the calculate thermodynamics properties , molecular weight and plot nature of the thermodynamic properties of various species.
FUNCTION – SPECIFIC HEAT
*This function uses the high temperature co efficient , low temperature coefficients , universal gas constant, and mintemp,midtempo,maxtemp values as the input parameters.
*Using ‘if’,’else ’ statement to check if the temp values are less tan or equal to the midtemp values and select the appropriate equation with high and low temperature co efficients is used to calculate the specific heat value denoted by c_p and return that value.
functionc_p=specific_heat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,R,temp,tmid)
if temp<= tmid
c_p=R*(b1+(b2.*temp)+(b3*(temp.^2))+(b4*((temp).^3))+(b5*((temp).^4)));
else
c_p=R*(a1+(a2*temp)+a3*((temp).^2)+a4*((temp).^3)+a5*((temp).^4));
end
end
FUNCTION-ENTHALPY
Then using the ‘if’ ‘else ’ conditional statement to check if the temp values aree less than or equal to the midtemp values and the appropriate equation with the proper maximum and minimum temperature co efficients is selected and the enthalpy is calculated and its value represented by H is returned.
function H = enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,R,temp,tmid)
if temp <= tmid
H = R*temp.*(b1+b2*(temp./2)+b3.*((temp.^2)./3)+b4.*((temp.^3)./4)+b5.*((temp.^4)./5)+b6./temp);
else
H = R.*temp.*(a1+a2.*(temp./2)+a3.*((temp.^2)./3)+a4.*((temp.^3)./4)+a5.*((temp.^4)./5)+a6./temp);
end
end
FUNCTION-MOLECULAR WEIGHT
1.A string array is created which contains the name of the atoms.
2.A variable atomic_wt is created which contains the atomoic weight of respective atoms in the string array.
3.then used ‘strcmp’ command to compare the name f the element with the string array.
function [w] = molecular_weight(element_name)
s = {'H' 'C' 'N' 'O' 'AR'};
%atomic weight of each elements
atomic_wt = [1 12 14 16 40];
w = 0;
for i =1:length(element_name)
for j = 1:length(s)
if strcmp(element_name(i),s(j))
w=w+atomic_wt(j);
% storing the position of matched element
position = j;
end
end
n = str2double(element_name(i));
if n>1
w=w+atomic_wt(position)*(n-1);
end
end
fprintf('Molecular weight of %s = %f',element_name,w)
disp(' ')
end
EXPLANATION-MAIN PROGRAM
1.The value of the universal gas constant R is given as 8.314 .
2.The given data is formated as therm.dat file
3.We use the 'fopen' command .the file is used to open the file in the read mode specified by'r' keyword 'fopen ' command.
4.The file is read line by line using 'fgetl' command.
5.The temperature array is defined using 'linspace ' command.
6.The 'for' loop is created loop through the species 53 data given in the file and calculate the required values.
7. The name of the element stored in a variable using 'strtok' command.
8.The next three line of the file is read which contain the high and low temp data of the particular specis.
9.The low temperature values are-a1,a2,a3,a4,a5,a6,a6
10.The high temperature are -b1,b2,b3,b4,b5,b6,b7
11.The next line is the reference to a function called specific heat which calculated and returns the value of specific heat using the high and low temp values as inputs.
12.The next line is the reference to a function called enthalpy which calculated and returns the value of enthalpy using the high and low temp values as inputs.
13.The next line is the reference to a function called entropy which calculated and returns the value of entropy using the high and low temp values as inputs.
14.The next line is the reference to a function calculated moleculae weight of the different element.
close all
clear
clc
%Inputs
R=8.314;
%Reading the data file
f1 = fopen('THERMO.dat','r');
%Reading the first line of the file
L1 = fgetl(f1);
%Getting the global minimum,maximum and 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 comments in the file
for i=1:3
fgetl(f1);
end
for i =1:53
%Reading the first line which includes the formation for species
l=fgetl(f1);
%Reading the nam eof the species
% element_name=strtok(1);
A=strsplit(l,' ');
element_name=(A{1});
tlow=str2num(A{end-3});
thigh=str2num(A{end-2});
tmid=str2num(A{end-1});
%Reading the three lines which contains 7 high temp coefficients & 7 low temp coefficients
l1 = fgetl(f1);
l2 = fgetl(f1);
l3 = fgetl(f1);
%Finding the position of E so that we can extract the coefficients
a=strfind(l1,'E');
b=strfind(l1,'E');
c=strfind(l3,'E');
% Extracting high temp coefficient
a1=str2double(l1(1:a(1)+3));
a2=str2double(l1(a(1)+4:a(2)+3));
a3=str2double(l1(a(2)+4:a(3)+3));
a4=str2double(l1(a(3)+4:a(4)+3));
a5=str2double(l1(a(4)+4:a(5)+3));
a6=str2double(l2(1:b(1)+3));
a7=str2double(l2(b(1)+4:b(2)+3));
%Extracting low temperature coefficient
b1=str2double(l2(b(2)+4:b(3)+3));
b2=str2double(l2(b(3)+4:b(4)+3));
b3=str2double(l2(b(4)+4:b(5)+3));
b4=str2double(l3(1:c(1)+3));
b5=str2double(l3(c(1)+4:c(2)+3));
b6=str2double(l3(c(2)+4:c(3)+3));
b7=str2double(l3(c(3)+4:c(4)+3));
%Calculating specific heat
c_p = specific_heat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,R,temp,tmid);
%Calculating enthalpy
H = enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,R,temp,tmid);
%Calculating entropy
S =entropy(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,R,temp,tmid);
%Plotting graphs
mkdir(element_name)
warning('off')
cd(element_name)
figure(1)
plot(temp,c_p)
xlabel('temperature')
ylabel('Specific heat')
title('Specific heat vs temperature of ',sprintf('%S ',element_name))
saveas(figure(1),'c_p vs temp.jpg')
figure(2)
plot(temp,H)
xlabel('temperature')
ylabel('Enthalpy')
title('Enthalpy vs temperature of ',sprintf('%S ',element_name))
saveas(figure(2),'H vs temp.jpg')
figure(3)
plot(temp,S)
xlabel('temperature')
ylabel('Entropy')
title('Entropy vs temperature of ',sprintf('%S ',element_name))
saveas(figure(3),'S vs temp.jpg')
cd ..
%Calculating molecular weights
M = molecular_weight(element_name);
end
O2:
N2:
CO2:
