in error-1 image, i have captured one error message where i have recived many times like every time while create the new folder for plots of each species
the errors says that in that folder their NO function programs
for that i have pasted all the function programs each time while creating the new folder
Detail explation for the entire program
Here we are creating programs for that extracts the 14 co-efficients and calculating the enthalpy, entropy, specific heats and Molecular weight. By Parsing NASA thermodynamic data.
We are reading the NASA thermodynamic data using fopen, fgetl and many more syntax’s
Here we are having formulas for Specific heat, enthaply and entropy
Cp = (a1 + a2 T + a3 T2 + a4 T3 +a5 T4) * R
H = (a1 + a2 T/2 + a3 T2/3 + a4 T3/4 + a5 T4/5 + a6/T) * R*T
S = (a1 lnT + a2 T + a3 T2/2 + a4 T3/3 + a5 T4/4 + a7) * R
By using the above formulas we are calculating the Specific heat, enthalpy and entropy
If we read the NASA Thermodynamic data, we having the 53 Species and for every Specie there are 14 temperature coefficients values and there are 3 local temperatures
In every 14 temperature coefficients the first 7 coefficients are high temperature coefficients and second 7 cefficients are low temperature coefficients.
Here R is universal gas constant = 8.314
Function Porgram for the Specific heat
For Specific heat we are having the formula
Cp = (a1 + a2 T + a3 T2 + a4 T3 +a5 T4) * R
function[Cp] = Specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,tempmid)
from above formula we are having the a1,a2,a3,a4,a5, for higher temperature coefficients and a8,a9,a10,a11,a12 are for the low temperature coefficients, Tis for temperature, R is the Universal gas constant and tempmid is the medium temperature for gobal temperature in the NASA thermodynamic data.
if (T<=tempmid)
if the temperature is the less than or equal to tempmid the we ae using the formula with lower temperature coefficients
Cp = (a8 + (a9*T)+(a10*(T).^2)+(a11*(T).^3) + (a12*(T).^4))*R;
Else we are using the higher temperature coefficients
Cp = (a1 + (a2*T)+(a3*(T).^2)+(a4*(T).^3) + (a5*(T).^4))*R;
Function Porgram for the Entropy
For Entropy we are having the formula
S = (a1 lnT + a2 T + a3 T2/2 + a4 T3/3 + a5 T4/4 + a7) * R
function[S] = Entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,tempmid)
from above formula we are having the a1,a2,a3,a4,a5,a7 for higher temperature coefficients and a8,a9,a10,a11,a12,a14 are for the low temperature coefficients, Tis for temperature, R is the Universal gas constant and tempmid is the medium temperature for gobal temperature in the NASA thermodynamic data.
if (T<=tempmid)
if the temperature is the less than or equal to tempmid the we ae using the formula with lower temperature coefficients
S = (a8*log(T) + (a9*T)+((a10*(T).^2)/2)+((a11*(T).^3)/3) + ((a12*(T).^4)/4) + a14)*R;
Else we are using the higher temperature coefficients
S = (a1*log(T) + (a2*T)+((a3*(T).^2)/2)+((a4*(T).^3)/3) + ((a5*(T).^4)/4) + a7)*R;
Function Porgram for the Enthaply
For Enthaply we are having the formula
H = (a1 + a2 T/2 + a3 T2/3 + a4 T3/4 + a5 T4/5 + a6/T) * R*T
function[H] = Enthaply(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,tempmid)
from above formula we are having the a1,a2,a3,a4,a5,a7 for higher temperature coefficients and a8,a9,a10,a11,a12,a14 are for the low temperature coefficients, Tis for temperature, R is the Universal gas constant and tempmid is the medium temperature for gobal temperature in the NASA thermodynamic data.
if (T<=tempmid)
if the temperature is the less than or equal to tempmid the we ae using the formula with lower temperature coefficients
H = (a8+((a9*T)/2)+((a10*(T).^2)/3)+((a11*(T).^3)/4) + ((a12*(T).^4)/5) + (a13./T)).*R.*(T);
Else we are using the higher temperature coefficients
H = (a1+((a2*T)/2)+((a3*(T).^2)/3)+((a4*(T).^3)/4) + ((a5*(T).^4)/5) + (a6./T)).*R.*T;
Function program for the Molecular weight of each Species
function mw = mol_w(Species_SP)
taking the ‘mw’ as molecular weight
as we are having 53 are common with 'H','C','O','N','A', these are taking as the atomic name and these molecular weights are using the molecular weights all the 53 species
by using the for loop we are creating this function
at initial stage we are having molecular weight = 0
for i = 1:length(Species_SP)
for j = 1:length(atomic_name)
if strcmpi(Species_SP(i),atomic_name(j))
mw = mw + atomic_weight(j);
here we are using STRCMPI which means compare with specices and the atomic name if it matches then the atomic weight is taken where the atomic name represents. In the same the molecular weight is calculated for all the species
n = str2num(Species_SP(i));
str2num is used to convert the string to the number
after that if n>1 then
mw = mw+atomic_weight(position)*(n-1)
which helps to find the molecular weight to the species
now we have all function programs
jumping into the main program
Clear all
Close all
Clc
To clear the old data and close all the plots and clc is for the clear the command window
For reading the file we are using the
Fopen(‘file name’,’r’)
Is the syntax for the fopen command
f1 = fopen('THERMO.dat','r');
for moving to next kine we are using the fgetl command
fgetl(assigned file name)
is the syntax for the fgetl command
fgetl(f1);
we are having the gobal temperatures in the THERMO file in the second line, so we are assigning the second line as ‘temp’
temp = fgetl(f1);
splitting the temp using strsplit(temp) at assigning ‘Z’
then writing command str2num(Z{2}) which is in Z line 2 number has to be take which will assign to temp low - templow = str2num(Z{2});
in the same way assigning the 3rd value for tempmid and 4th value for the temphigh
templow = str2num(Z{2});
tempmid = str2num(Z{3});
temphigh = str2num(Z{4});
for the value T we are assigning the 500 value in between the templow and temphigh by using the linspace command
T = linspace(templow,temphigh,500);
The next 3 lines I the data file are not useful to us,we are skiping those lines by writing the I = 1:3 and the going to the next line using the fgetl and assigning to the tline by using the for loop.
for i = 1:3
tline = fgetl(f1);
end
opening the for loop again which here we are assigning for all the 53 species, in this line we are having the local temperatures in every species and we have capture all the local temps for that we are using the
A = strsplit(tline,' ');
Command where in the tline it removes all the spaces for easily capturing the temperature values
For capturing the Species name
Species_SP = (A{1});
Where the ‘A’ is assigned to the tline without spaces and the 1st element the line ‘A’ is the SPECIE which will assigned to the Species_SP
Creating the low,mid and high temps in the local temperatures
tlow = str2num(A{end-3});
where end-3, means the 3rd element or number from the ending of line ‘A’ which was allocated to the tlow
in the same way tmid and thigh was allacated
tlow = str2num(A{end-3});
tmid = str2num(A{end-2});
thigh = str2num(A{end-1});
finding the co-efficients
going to the next line by using fgetl and it will assigned to the tline1
‘2.56942078E+00-8.59741137E-05 4.19484589E-08-1.00177799E-11 1.22833691E-15’
If you see the above line, we are in the same line right now, there are 5 co-efficient in this line and we have assigned the each co-efficient for each value. For that we are using the ‘E’ which is common in every co-efficent
To identify the ‘E’ position in the line we are writing
pos1 = strfind(tline1,'E');
which assigned as pos1
pos1 will be at 12 27 42 57 72
now for the coefficients
a1 = str2num(tline1(1:pos1(1)+3));
explaining the above command which is 1:pos1(1) means pos1(1) is at 12 it will be like 1 no to 12 number and +3 defines the next 3 numbers after ‘E’ then,
str2num is converting the string to number and assigned to a1
a1 will be 2.56942078e+00
and for a2
a2 = str2num(tline1(pos1(1)+4:pos1(2)+3));
explaining the above command which is (tline1(pos1(1)+4:pos1(2)+3))
starting from the pos1(1)+4 means starting from the 12+4 place means 16 position to the pos1(2)+3 means pos1(2) means 27 position is next ‘E’ position and +3
str2num is converting the string to number and assigned to a2
a2 will be 8.59741137E-05
in the same way we have to capture next 3 values from the tline1
and creating the tline2 as like tline1 and capturing all 5 values from tline2
same like tline2, tline3 will be created and 4 values will be capture all the 14 values in the one specie are capture.
In the for loop we will capture all the 53 species in the same way.
Now, adding the all function syntax for Specific heat, Entropy, Enthaply & molecular weight.
Creating the molecular weight with the help of the fopen and fclose
And now creating the folder for every species with the plots in it
Plot1 Cp vs T
Using plot command T,CP, with line width 2 and color red
plot(T,Cp,'linewidth',2,'Color','r')
adding xlabel and y label to the plot
xlabel('Temperature(K)')
ylabel('Specfic heat(KJ)')
adding the title Specfic heat vs Temperature
grid on
saveas command in name figure(1) with specific heat
saveas(figure(1),”Specfic heat.jpg')
same for the plot Entropy vs Temperature and
same for the plot Enthaply vs Temperature
Main Program for thr file Parsing
clear all
close all
clc
%readingthe thermo file
f1 = fopen('THERMO.dat','r');
fgetl(f1);
%temperature from the file
temp = fgetl(f1);
%spliting
Z = strsplit(temp);
%string cells to numbers
templow = str2num(Z{2});
tempmid = str2num(Z{3});
temphigh = str2num(Z{4});
%tempurature range
T = linspace(templow,temphigh,500);
%Universal Gas Constant
R = 8.314;
%skip next 3 lines
for i = 1:3
tline = fgetl(f1);
end
for j = 1:53
tline = fgetl(f1);
A = strsplit(tline,' ');
Species_SP = (A{1});
%local temp's
tlow = str2num(A{end-3});
tmid = str2num(A{end-2});
thigh = str2num(A{end-1});
%finding the coefficents
%line1
tline1 = fgetl(f1);
pos1 = strfind(tline1,'E');
%finding higher coefficents
a1 = str2num(tline1(1:pos1(1)+3));
a2 = str2num(tline1(pos1(1)+4:pos1(2)+3));
a3 = str2num(tline1(pos1(2)+4:pos1(3)+3));
a4 = str2num(tline1(pos1(3)+4:pos1(4)+3));
a5 = str2num(tline1(pos1(4)+4:pos1(5)+3));
%line 2
tline2 = fgetl(f1);
pos2 = strfind(tline2,'E');
a6 = str2num(tline2(1:pos2(1)+3));
a7 = str2num(tline2(pos2(1)+4:pos2(2)+3));
%lower coefficents
a8 = str2num(tline2(pos2(2)+4:pos2(3)+3));
a9 = str2num(tline2(pos2(3)+4:pos2(4)+3));
a10 = str2num(tline2(pos2(4)+4:pos2(5)+3));
%line 3
tline3 = fgetl(f1);
pos3 = strfind(tline3,'E');
a11 = str2num(tline3(1:pos3(1)+3));
a12 = str2num(tline3(pos3(1)+4:pos3(2)+3));
a13 = str2num(tline3(pos3(2)+4:pos3(3)+3));
a14 = str2num(tline3(pos3(3)+4:pos3(4)+3));
%Sepecific heat calculation
Cp = Specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,tempmid);
%entropy Calculation
S = Entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,tempmid);
%enthaply Calculation
H = Enthaply(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,tempmid);
%molecular weight calculation
mw = mol_w(Species_SP)
%writing molecular weight file
mw_file = fopen('molecular weight.txt','w');
fprintf(mw_file,'molecular weight of %sis %d',Species_SP,mw);
fclose(mw_file);
%plot Cp,H,S with temperature
cd('D:\matlab\parsing nasa thermodynamic\species')
mkdir(Species_SP);
cd(Species_SP)
%plot 1 Cp vs T
figure(1)
plot(T,Cp,'linewidth',2,'Color','r')
xlabel('Temperature(K)')
ylabel('Specfic heat(KJ)')
title(sprintf('Specfic heat vs Temperature range of %s',Species_SP))
grid on
saveas(figure(1),'Specfic heat.jpg')
%plot 2 S vs T
figure(2)
plot(T,S,'linewidth',2,'Color','g')
xlabel('Temperature(K)')
ylabel('Entropy(KJ/mol-K)')
title(sprintf('Entropy vs Temperature range of %s',Species_SP))
grid on
saveas(figure(2),'entropy.jpg')
%plot 3 H vs T
figure(3)
plot(T,H,'linewidth',2,'Color','b')
xlabel('Temperature(K)')
ylabel('Enthaply(KJ/mol)')
title(sprintf('Enthaply vs Temperature range of %s',Species_SP))
grid on
saveas(figure(2),'Enthaply.jpg')
end
Function Program For the Molecular weight
%calculation of molecular weight of each species function
function mw = mol_w(Species_SP)
% Creating standard atomic name
atomic_name = ['H','C','O','N','A'];
%Storing atomic weight inatomic number
atomic_weight = [1.0079 12.011 15.999 14.0067 39.948];
%molecular mass is 0
mw = 0;
%creating if pogram
for i = 1:length(Species_SP)
for j = 1:length(atomic_name)
if strcmpi(Species_SP(i),atomic_name(j))
mw = mw + atomic_weight(j);
% Storing the position for furthue use
position = j;
end
end
n = str2num(Species_SP(i));
if n>1
mw = mw+atomic_weight(position)*(n-1);
end
end
fprintf('molecular weight of %s :-',Species_SP);
fprintf('%s',mw);
disp(' ')
end
Function program for Enthaply
%Enthaply function
function[H] = Enthaply(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,tempmid)
if (T<=tempmid)
H = (a8+((a9*T)/2)+((a10*(T).^2)/3)+((a11*(T).^3)/4) + ((a12*(T).^4)/5) + (a13./T)).*R.*(T);
else
H = (a1+((a2*T)/2)+((a3*(T).^2)/3)+((a4*(T).^3)/4) + ((a5*(T).^4)/5) + (a6./T)).*R.*T;
end
end
Function Program for the Entropy
%Entropy function
function[S] = Entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,tempmid)
if (T<=tempmid)
S = (a8*log(T) + (a9*T)+((a10*(T).^2)/2)+((a11*(T).^3)/3) + ((a12*(T).^4)/4) + a14)*R;
else
S = (a1*log(T) + (a2*T)+((a3*(T).^2)/2)+((a4*(T).^3)/3) + ((a5*(T).^4)/4) + a7)*R;
end
end
Function Program for the Specific heat
%Specific heat function
function[Cp] = Specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,tempmid)
if (T<=tempmid)
Cp = (a8 + (a9*T)+(a10*(T).^2)+(a11*(T).^3) + (a12*(T).^4))*R;
else
Cp = (a1 + (a2*T)+(a3*(T).^2)+(a4*(T).^3) + (a5*(T).^4))*R;
end
end