% main code
clear all
close all
clc
% Reading Nasa thermodynamic data file
f1=fopen('THERMO.dat','r');
gt= fgetl(f1);
% reading temperatures
global_temp=fgetl(f1);
%spliting the string
g=strsplit(global_temp);
%converting string to number
t=str2num(global_temp);
global_low_temp=t(1);
global_medium_temp=t(2);
global_high_temp=t(3);
%Gas constant R(J/(mol-k)) and temperature range
temp=linspace(global_low_temp,global_high_temp,1000)
R=8.314
%skipping unwanted comment lines
for i=1:3
gt= fgetl(f1)
end
for i=1:53
gt=fgetl(f1);
A=strsplit(gt,'');
species_sp=(A{1});
line_1=fgetl(f1);
a=strfind(line_1,'E');
%extrating first 7 high temperature coefficients
a1=line_1(1:a(1)+3);
a1=str2num(a1);
a2=line_1(a(1)+4:a(2)+3);
a2=str2num(a2);
a3=line_1(a(2)+4:a(3)+3);
a3=str2num(a3);
a4=line_1(a(3)+4:a(4)+3);
a4=str2num(a4);
a5=line_1(a(4)+4:a(5)+3);
a5=str2num(a5);
line_2=fgetl(f1);
b=strfind(line_2,'E');
a6=line_2(1:b(1)+3);
a6=str2num(a6);
a7=line_2(b(1)+4:b(2)+3);
a7=str2num(a7);
%extractuing low temperature coefficients
a8=line_2(b(2)+4:b(3)+3);
a8=str2num(a8);
a9=line_2(b(3)+4:b(4)+3);
a9=str2num(a9);
a10=line_2(b(4)+4:b(5)+3);
a10=str2num(a10);
line_3=fgetl(f1);
c=strfind(line_3,'E');
a11=line_3(1:c(1)+3);
a11=str2num(a11);
a12=line_3(c(1)+4:c(2)+3);
a12=str2num(a12);
a13=line_3(c(2)+4:c(3)+3);
a13=str2num(a13);
a14=line_3(c(3)+4:c(4)+3);
a14=str2num(a14);
%calcuating specific heat
cp= specificheat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,temp,global_medium_temp,R);
%clcuating enthalpy
s = entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,R,temp,global_medium_temp);
%calcuating entropy
H = Enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,global_medium_temp,R,temp);
%calcuating molecular Weight of species
molecularweight(i) = molecular_weight(species_sp);
%writing and reading moleculrweight file
molecular_list = fopen('molecular weight.txt','a+');
fprintf(molecular_list,'molecular weight of %s is %d n',species_sp,molecularweight );
fclose(molecular_list);
%plotting the graphs for h,cp,temp,species_sp
%creating a folder to save data
cur_dir=pwd;
mkdir(species_sp);
cd (species_sp);
% plotting specific heat and temperature range
figure(1)
plot(temp,cp,'linewidth',3,'color','b');
xlabel('temperature [K]');
ylabel('specific heat [KJ]');
title(sprintf('specific heat changing with temperature range of %s',species_sp));
grid on
saveas(figure(1),'specific heat.jpg');
%plotting graph between enthalpy and temperature range
figure(2)
plot(temp,H,'linewidth',3,'color','r');
xlabel('temperature [K]');
ylabel('Enthalpy [KJ/mol]');
title(sprintf('Enthalpy changing with temperature range of %s',species_sp));
grid on
saveas(figure(2),'Enthalpy.jpg');
%plotting graph between entropy and temperature range
figure(3)
plot(temp,s,'linewidth',3,'color','g');
xlabel('temperature [K]');
ylabel('Entropy [KJ/mol_K]');
title(sprintf('Entropy changing with temperature range of %s',species_sp));
grid on
saveas(figure(3),'Entropy.jpg');
cd(cur_dir)
%molecular weight
N = molecular_weight(species_sp)
f2=fopen('molecular mass','w');
fprintf('%s n %s n %f n','species','molecular mass',N)
fclose(f2)
end
%function for calcuating specific heat
function cp= specificheat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,temp,global_medium_temp,R)
if temp <= global_medium_temp % low temp coefficients are taken
cp=(a8+(a9*temp)+(a10*(temp).^2)+(a11*(temp).^3)+(a12*(temp).^4))*R;
else % high temp coefficients are taken
cp= (a1+(a2*temp)+(a3*(temp).^2)+(a4*(temp).^3)+(a5*(temp).^4))*R;
end
%function for calcuating entropy
function s = entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,R,temp,global_medium_temp)
if temp <= global_medium_temp %low temperature coefficients are taken
s=(a8.*log(temp)+(a9.*temp) + (a10*((temp).^2)/2) + (a11*((temp).^3)/3) + (a12 *((temp).^4)/4) + a14).*R
else %high temperature coefficients are taken
s= (a1.*log(temp)+(a2.*temp)+ (a3*((temp).^2)/2)+ (a4*((temp).^3)/3) + (a5*((temp).^4)/4)+a7).*R
end
%function for calcuating enthalpy
function H = Enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,global_medium_temp,R,temp)
if temp <= global_medium_temp %Minimum temperature coefficients are taken
H=(a8+(a9*temp)/2 + (a10*((temp).^2)/3) + (a11*((temp).^3)/4) + (a12*((temp).^4)/5) + (a13./temp)).*temp*R;
else % maximum temperature coefficients are taken
H=(a1+(a2*temp)/2 + (a3*((temp).^2)/3) + (a4*((temp).^3)/4) + (a5*((temp).^4)/5) + (a6./temp)).*temp*R;
end
% Function for calcuating molecular weight
function molecularweight = molecular_weight(species_sp)
atomic_name=['C','O','H','A','N'];
fprintf('Molecular weight of ');
fprintf(species_sp);
fprintf(' : ');
%Atomic weights
atomic_weight=[12.011 15.999 1.008 39.984 14.007];
% Starting molecular weight with 0
molecularweight = 0;
% Adding atomic number
for i = 1:length(species_sp)
for j = 1:length(atomic_name)
if strcmpi(species_sp(i),atomic_name(j))
molecularweight = molecularweight + atomic_weight(j);
position=j;
end
end
n = str2num(species_sp(i));
if n>1
molecularweight = molecularweight + atomic_weight(position)*(n-1);
end
end
fprintf('molecular weight %s :','species_sp');
fprintf('%d',molecularweight);
disp('');
end
% explanation:
% By using this code the coefficients and other data are read from
% thermo.da file.
%we used fopen and agetl commands to obtain information about the files
%which are open and read the line from file.
%To get each and every values of the species we used for loop which help
%in doing iterations.
%By using strsplit we arranged the co efficients in order and we have used
%strfind to find the location of variable E which is used to extract the
%coefficients from the file.
% functions have been defined to calcuate specific heat,enthalpy,entropy
% and molecular weight.
%By using mkrdir command i have created a folder.
%plot of temp vs(enthaply,entropy,specific heat ) has been plotted and
%images are stored in jpeg format.