Project 1 - Parsing NASA thermodynamic data
Aim of project :-
1. Calculate thermodynamic properties like Cp(specific heat), H(enthalpy) and S(entropy)
3.Generate plot between - a)Temprature vs Cp.
b) Temprature vs H.
C)Temprature vs S.
attach only for o2, co2 and n2.
2.Calculate moleculer weight as well.
Equations used:-
1)calculation of Cp , H, and S in high temp
Cp = (a1+a2*T+a3*T.^2+a4*T.^3+a5*T.^4)*R
H = (a1+0.5*a2*T+0.333*a3*T.^2+0.25*a4*T.^3+0.2*a5*T.^4+a6/T).*R.*T
S = (a1*log(T)+a2*T+0.5*a3*T.^2+0.333*a4*T.^3+0.25*a5*T.^4+a7)*R
2)calculation of Cp , H, and S in low temp
Cp = (a8+a9*T+a10*T.^2+a11*T.^3+a13*T.^4)*R
H = (a8+0.5*a9*T+0.333*a10*T.^2+0.25*a11*T.^3+0.2*a12*T.^4+a13).*R.*T
S = (a8*log(T)+a9*T+0.5*a10*T.^2+0.333*a11*T.^3+0.25*a13*T.^4+a14)*R
where
R = Characterstic gas constant
T = local temprature
a1 to a14 = cofficients
CODE WTHE EXPLATION
clear all
close all
clc
f1 = fopen('THERMO.dat','r') % open"THERMO" file for reading
for i=1:5 %creating for loop for reading 1st five lines
fgetl(f1) %reading line
end
for j = 1:53 %initiation of for loop for all 53 species
line_1 = fgetl(f1); %reading line
a = strsplit(line_1,' '); %split the string into cell
species_name = a{1} %extract character
len = length(a); %define length
localmid = str2num(a{len-1}) % extracting local mid temp and converted into numbers
localhigh = str2num(a{len-2}) % extracting local high temp and converted into numbers
locallow = str2num(a{len-3}) % extracting local low temp and converted into numbers
coff_line1 = fgetl(f1); % reading line 1 of cofficients
poss_E = strfind(coff_line1,'E') % finding the position of 'E'
% finding out all cofficients in line 1
a1 = str2num(coff_line1(1:poss_E(1)+3));
a2 = str2num(coff_line1(poss_E(1)+4:poss_E(2)+3));
a3 = str2num(coff_line1(poss_E(2)+4:poss_E(3)+3));
a4 = str2num(coff_line1(poss_E(3)+4:poss_E(4)+3));
a5 = str2num(coff_line1(poss_E(4)+4:poss_E(5)+3));
coff_line2 = fgetl(f1); % reading line 2 of cofficients
poss_E2 = strfind(coff_line2,'E'); % finding the position of 'E'
% finding out all cofficients in line 2
a6 = str2num(coff_line2(1:poss_E2(1)+3))
a7 = str2num(coff_line2(poss_E2(1)+4:poss_E2(2)+3));
a8 = str2num(coff_line2(poss_E2(2)+4:poss_E2(3)+3));
a9 = str2num(coff_line2(poss_E2(3)+4:poss_E2(4)+3));
a10 = str2num(coff_line2(poss_E2(4)+4:poss_E2(5)+3));
% finding out all cofficients in line 3
coff_line3 = fgetl(f1);
poss_E3 = strfind(coff_line3,'E');
a11 = str2num(coff_line3(1:poss_E3(1)+3));
a12 = str2num(coff_line3(poss_E3(1)+4:poss_E3(2)+3));
a13 = str2num(coff_line3(poss_E3(2)+4:poss_E3(3)+3));
a14 = str2num(coff_line3(poss_E3(3)+4:poss_E3(4)+3));
T = linspace(locallow,localhigh,10); % defining array of temp
R= 8.134; % defining characterstic gas constant
for i = 1:length(T) % creating for loop for selecting temp from temp array
% initiation of if else loop
if T>localmid % defining the condition
% calculation of Cp , H, and S for if condition
Cp = (a1+a2*T+a3*T.^2+a4*T.^3+a5*T.^4)*R
H = (a1+0.5*a2*T+0.333*a3*T.^2+0.25*a4*T.^3+0.2*a5*T.^4+a6/T).*R.*T
S = (a1*log(T)+a2*T+0.5*a3*T.^2+0.333*a4*T.^3+0.25*a5*T.^4+a7)*R
else
% calculation of Cp , H, and S for else condition
Cp = (a8+a9*T+a10*T.^2+a11*T.^3+a13*T.^4)*R
H = (a8+0.5*a9*T+0.333*a10*T.^2+0.25*a11*T.^3+0.2*a12*T.^4+a13).*R.*T
S = (a8*log(T)+a9*T+0.5*a10*T.^2+0.333*a11*T.^3+0.25*a13*T.^4+a14)*R
end
cd('E:\Nasa') % changing the directory for new path
mkdir(species_name) % creating directory
cd(species_name)
figure(1) % generating seperate figure
plot(T,Cp) %plot between T and Cp
saveas(figure(1),'T vs cp.jpg') %saving figure and specify its name
figure(2)
plot(T,H)
saveas(figure(2),'T vs H.jpg')
figure(3)
plot(T,S)
saveas(figure(1),'T vs s.jpg')
end
end
main command and function used:-
1.fopen - use to make the given file readable to matlab.
2.fgetl - to read single line
3.strsplit - to split given string into cell
4.str2num - to convert string itno numbers
5.strfind - to find out the position of a charater of string
6.cd - to change the directory
7.mkdir - to create the directory
8.saveas - to save a given plot
Code output
1. for O2
T vs H
T vs S
FOR CO2
T vs Cp
T vs H
T vs S
FOR N2
T vs Cp
T vs H
T vs S
CODE FOR CALCULATING MOLECULER WEIGHT
atomic_names = ['H','C','O','N','A'];
atomic_weight = [1,12,16,14,40];
w = 0;
for i = 1:length(species_name)
for j =1:length(atomic_names)
if species_name(i)==atomic_weight(j)
w=w+atomic_weight(j);
position = j;
end
end
n = str2num(species_name(i));
if n>=1
w = w+atomic_weight(position).*(n-1)
end
end
end