Project 1 - Parsing NASA thermodynamic data
AIM
Write a function that extracts the 14 co-efficients and calculates the enthalpy, entropy and specific heats for all the species in the data file.
OBJECTIVE
Calculate the molecular weight of each species and display it in the command window.
Plot the Cp, Enthalpy and Entropy for the local temperature range (low temperature : high temperature) specific for each species.
Save the plots as images with appropriate names and dump them in separate folders for each species with suitable name
CODE EXPLANATION
1. main code
clear all
close all
clc
f1=fopen('THERMO.dat','r') // command to open and read the file 'THERMO.dat' which contains the details of species from which the details are to be parsed to the program
prompt = 'enter the specie name '; // command to display
x = input(prompt); // takes the input variable given by the user and returns the entered text, without evaluating the input as an expression. to get the specie name from the user
fgetl(f1); // returns the next line of the specified file,
fgetl(f1);
fgetl(f1);
fgetl(f1);
fgetl(f1);
for i= 1:53 // loop to move the pointer through everyspecie name present in the file, inorder to find the user mentioned specie and extract the required datas related with it
tline=fgetl(f1); // line contains the datas releated with local temperature ranges
temp_data= strsplit(tline," "); //splits the string into pieces wrt to the " " (spaces) bw the characters and stores it into variable temp_data
a=temp_data{1}; // variable a stores the first character of data stored in temp_data
if strcmp(a,x);// compares if the stored in variables a and x are same,
if they are same will go on further with the rest of the codes else will skip next 3lines and moves the pointer to 4th line from curent location and again resumes from start untill it locates the given specie in the file
the logic behind using for , if loops and the above statements
xx = length(temp_data);// variabale xx stores the length or the no of strings in the variable temp_data
low_temp= temp_data{xx-3}; // low_temp willl carry the value present in the location (xx-3)ie length-3, as xx will point towards the last location so xx-3 will locate 3 pos before the last .
low_temp= str2num(low_temp);// the value stored will be in the string format ,to convert the format from string to number and store it back to variable low_temp
high_temp= temp_data{xx-1};
high_temp= str2num(high_temp);
mid_temp= temp_data{xx-2};
mid_temp= str2num(mid_temp);
tline1=fgetl(f1); // gets the pointer to the next line so as to read the datas present in the line - details of coeficients and store them into variable tline1
a=findstr(tline1,'E');// every coeficients have word 'E' as common s as to seperate them using this criteria, will locate where all the word 'E' is present throughout the line and stores the location into the variable 'a'
a1= tline1(1:a(1)+3);// a1 will store the characters from location from 1st to a(1) (which represents the location where the letter 'E' first appears) along with '+' next 3 characters from location of a(1) in the variable tline1
a1= str2num(a1); // the string characters taken are to be converted to numbers for further equating
a2=tline1(a(1)+4:a(2)+3);
a2=str2num(a2);
a3=tline1(a(2)+4:a(3)+3);
a3=str2num(a3);
a4=tline1(a(3)+4:a(4)+3);
a4=str2num(a4);
a5=tline1(a(4)+4:a(5)+3);
a5=str2num(a5);
tline2 = fgetl(f1);// gets to the next line where the next set of coeficent datas are stored
a6=findstr(tline2,'E');
a6= tline2(1:a(1)+3);
a6= str2num(a6);
a7=tline2(a(1)+4:a(2)+3);
a7=str2num(a7);
a8=tline2(a(2)+4:a(3)+3);
a8=str2num(a8);
a9=tline2(a(3)+4:a(4)+3);
a9=str2num(a9);
a10=tline2(a(4)+4:a(5)+3);
a10=str2num(a10);
tline3 = fgetl(f1);
a11=findstr(tline3,'E');// gets to the next line where the next set of datas are stored
a11= tline3(1:a(1)+3);
a11= str2num(a11);
a12=tline3(a(1)+4:a(2)+3);
a12=str2num(a12);
a13=tline3(a(2)+4:a(3)+3);
a13=str2num(a13);
a14=tline3(a(3)+4:a(4)+3);
a14=str2num(a14);// the 14 set of coeficents are extracted as string and conv into numbers from the file
lt = linspace(low_temp,high_temp,100); // creates an array with values starting from low to high temp into 100 parts, so as to find cp,h and s with respect to temperature values
cp = daset(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,low_temp,high_temp,mid_temp);// function calling statement to call the function that calculates and returns the specific heat into the array cp for each temp value from the array 'lt'
H= daset1(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,low_temp,high_temp,mid_temp);//function calling statement to call the function that calculates and returns the enthalpy calculated
S=daset2(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,low_temp,high_temp,mid_temp);//function calling statement to call the function that calculates and returns the entropy caluclated
mol_wt=fre(x); // function calling statement to call the function that calculates and returns and displays the value of molecular weight of the given specie
figure(1) // commands to plot the graph
plot( lt,cp,'linewidth',3,'color','b');
xlabel('temperature')
ylabel('specific_heat')
figure(2)// commands to plot the graph
plot( lt,H,'linewidth',3,'color','b');
xlabel('temperature')
ylabel('enthalpy')
figure(3)// commands to plot the graph
plot( lt,S,'linewidth',3,'color','b');
xlabel('temperature')
ylabel('entropy')
disp('moleuclar weight of the given specie = ')// command to disp the molecular weight calculated
disp(mol_wt)
mkdir (['C:\Users\Student\Documents\MATLAB\PLOTS\',x]) // command that creates a newfolder in the specified address with the name specified in the variable 'x'
cd (['C:\Users\Student\Documents\MATLAB\PLOTS\',x])// after the folder has been created the
saveas(1,'SP vs TMP .png') // command to save the plots generated as an image
saveas(2,'H vs TMP .png')
saveas(3,'S vs TMP .png')
break// command to break the conditional statement(if i= 1:53) so that it wont loop again after the specie has been found in the file
else // else condition for the conditional statement (if i= 1:53)
fgetl(f1);
fgetl(f1);
fgetl(f1);
end
end
2. function to calculate specific heat
function[cp]=daset(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,low_temp,high_temp,mid_temp)
R= 8.314;
lt = linspace(low_temp,high_temp,100);
for i=1:length(lt)
if lt(i)> mid_temp// conditional statement to search if the temp value taken from the temp array 'lt' is greater than mid_temp then first 7 coeficents ie a1 to a7 are used in the equation else the coeficents a8 to a14 are used
cp(i) = ( a1+ a2*lt(i) + a3*(lt(i)^2) + a4*(lt(i)^3) + a5*(lt(i)^4) )*R;// equation to calculate the specific heat using 1st 7 coeficents
else
cp(i) = ( a8+ a9*lt(i) + a10*(lt(i)^2) + a11*(lt(i)^3) + a12*(lt(i)^4) )*R;//equation to calculate the specific heat using 2nd set of coeficents
end
end
end
3. function to calculate enthalpy
function[H]=daset1(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,low_temp,high_temp,mid_temp)
R= 8.314;
lt = linspace(low_temp,high_temp,100);
for i=1:length(lt)
if lt(i)> mid_temp// conditional statement to search if the temp value taken from the temp array 'lt' is greater than mid_temp then first 7 coeficents ie a1 to a7 are used in the equation else the coeficents a8 to a14 are used
H (i) = (a1+ a2*lt(i)/2 + a3*(lt(i)^2)/3 + a4*(lt(i)^3)/4 + a5*(lt(i)^4)/5 + (a6/lt(i))) * (R*lt(i));// equation to calculate enthalpy using 1st 7 coeficents
else
H (i) = (a8+ a9*lt(i)/2 + a10*(lt(i)^2)/3 + a11*(lt(i)^3)/4 + a12*(lt(i)^4)/5 + (a13/lt(i))) * (R*lt(i));// equation to calculate enthalpy using 2nd set of coeficents
end
end
end
4. function to calculate entropy
function[S]=daset2(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,low_temp,high_temp,mid_temp)
R= 8.314;
lt = linspace(low_temp,high_temp,100);
for i=1:length(lt)
if lt(i)> mid_temp
S (i) = (a1*log(lt(i)) + a2*lt(i) + a3*(lt(i)^2)/2 + a4*(lt(i)^3)/3 + a5+(lt(i)^4)/4 + a7 )* R;
else
S (i) = (a8*log(lt(i)) + a9*lt(i) + a10*(lt(i)^2)/2 + a11*(lt(i)^3)/3 + a12+(lt(i)^4)/4 + a14 )* R;
end
end
end
5. function to calculate molecular weight
function[mol_wt] = fre(x)
name=(x); // gets the species name that is stored in the variable x
specie=["O", "H","C","N","A"];// array containing letters that form the name of the species
mass=[16,1,12,14,40];// array containing the numbers or individual molecular mass of each element and are arranged so as that it has the same corresponding location as the letters, so if any letter is taken its location or position if located in this array would give
represent the molecular weight of the letter or element
mol_wt=0; // initially the value of variable molecular weight is defined as zero
for i= 1:length(name); // stores the lenght of string stored (name of the specie defined by the user)
a= name(i); // variable 'a' stores letter in the position defined by name(i)
q = str2num(a);// command to convert a string to number
c=length(q);// c stores the value of lenght of the element in the variable q
if c == 0;// if the element present in q is a string, after the str2num command the value of c will be null as an alphabet string cannot be converted into a number format, this conditional statement will ensure that the value of 'a' is character
for j= 1:length(specie) ; // looping statement so as to check 'a' throughout with all elements present in the array specie
v =strcmp(a,specie(j));// command to check wether the strings a and the string indicated in the position'specie(j)' are same,
if v == 1;// if the above statement is true the logical operator will return value as 1 to the variable 'v'
mol_wt= mol_wt+mass(j); // equation, mol_wt+ mass(j)- will return the value of molecular weight corresponding of the element
end
end
else
p = name(i-1);/ if its not a letter, it'll be a number then the no is related to the letter just before the position of the number present , so the value of the letter present before the location of the number is to be found,
to find its value its again compared or searched in the array 'specie' where the letters are stored and the corresponding loction in the array mass would return the value of its molecular mass which is tobe multiplied
with the number
for b= 1:length(specie) ;// conditional statement
u=strcmp(p,specie(b));// comparing command to check wether p and the element in the location specie(b)
if u==1// if true logical statement will return 1 , using the concept confirming that the elements compared are same
rs = mass(b)*q; // the number in the name of specie is multiplied with the molecular mass of the element present before the number
mol_wt=mol_wt-mass(b);// as the molecular mass of the letter before the location of the number in the specie name as it is a letter, would be accepted and its molecular mass would be added in the above statement
mol_wt=mol_wt+rs;// value of rs to be added to variable mol_wt which stores the value of molecular weight of specie mentioned by the user
end
end
end
end
end
% main code
clear all
close all
clc
f1=fopen('THERMO.dat','r')
prompt = 'enter the specie name ';
x = input(prompt);
fgetl(f1);
fgetl(f1);
fgetl(f1);
fgetl(f1);
fgetl(f1);
for i= 1:53
tline=fgetl(f1);
temp_data= strsplit(tline," ");
a=temp_data{1};
if strcmp(a,x);
xx = length(temp_data);
low_temp= temp_data{xx-3};
low_temp= str2num(low_temp);
high_temp= temp_data{xx-1};
high_temp= str2num(high_temp);
mid_temp= temp_data{xx-2};
mid_temp= str2num(mid_temp);
tline1=fgetl(f1);
a=findstr(tline1,'E');
a1= tline1(1:a(1)+3);
a1= str2num(a1);
a2=tline1(a(1)+4:a(2)+3);
a2=str2num(a2);
a3=tline1(a(2)+4:a(3)+3);
a3=str2num(a3);
a4=tline1(a(3)+4:a(4)+3);
a4=str2num(a4);
a5=tline1(a(4)+4:a(5)+3);
a5=str2num(a5);
tline2 = fgetl(f1);
a6=findstr(tline2,'E');
a6= tline2(1:a(1)+3);
a6= str2num(a6);
a7=tline2(a(1)+4:a(2)+3);
a7=str2num(a7);
a8=tline2(a(2)+4:a(3)+3);
a8=str2num(a8);
a9=tline2(a(3)+4:a(4)+3);
a9=str2num(a9);
a10=tline2(a(4)+4:a(5)+3);
a10=str2num(a10);
tline3 = fgetl(f1);
a11=findstr(tline3,'E');
a11= tline3(1:a(1)+3);
a11= str2num(a11);
a12=tline3(a(1)+4:a(2)+3);
a12=str2num(a12);
a13=tline3(a(2)+4:a(3)+3);
a13=str2num(a13);
a14=tline3(a(3)+4:a(4)+3);
a14=str2num(a14);
lt = linspace(low_temp,high_temp,100);
cp = daset(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,low_temp,high_temp,mid_temp);
H= daset1(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,low_temp,high_temp,mid_temp);
S=daset2(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,low_temp,high_temp,mid_temp);
mol_wt=fre(x);
figure(1)
plot( lt,cp,'linewidth',3,'color','b');
xlabel('temperature')
ylabel('specific_heat')
figure(2)
plot( lt,H,'linewidth',3,'color','b');
xlabel('temperature')
ylabel('enthalpy')
figure(3)
plot( lt,S,'linewidth',3,'color','b');
xlabel('temperature')
ylabel('entropy')
disp('moleuclar weight of the given specie = ')
disp(mol_wt)
mkdir (['C:\Users\Student\Documents\MATLAB\PLOTS\',x])
cd (['C:\Users\Student\Documents\MATLAB\PLOTS\',x])
saveas(1,'SP vs TMP .png')
saveas(2,'H vs TMP .png')
saveas(3,'S vs TMP .png') in
break
else
fgetl(f1);
fgetl(f1);
fgetl(f1);
end
end
% function to calculate specific heat
function[cp]=daset(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,low_temp,high_temp,mid_temp)
R= 8.314;
lt = linspace(low_temp,high_temp,100);
for i=1:length(lt)
if lt(i)> mid_temp
cp(i) = ( a1+ a2*lt(i) + a3*(lt(i)^2) + a4*(lt(i)^3) + a5*(lt(i)^4) )*R;
else
cp(i) = ( a8+ a9*lt(i) + a10*(lt(i)^2) + a11*(lt(i)^3) + a12*(lt(i)^4) )*R;
end
end
end
%function to calculate enthalpy of given specie
function[H]=daset1(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,low_temp,high_temp,mid_temp)
R= 8.314;
lt = linspace(low_temp,high_temp,100);
for i=1:length(lt)
if lt(i)> mid_temp
H (i) = (a1+ a2*lt(i)/2 + a3*(lt(i)^2)/3 + a4*(lt(i)^3)/4 + a5*(lt(i)^4)/5 + (a6/lt(i))) * (R*lt(i));
else
H (i) = (a8+ a9*lt(i)/2 + a10*(lt(i)^2)/3 + a11*(lt(i)^3)/4 + a12*(lt(i)^4)/5 + (a13/lt(i))) * (R*lt(i));
end
end
end
%function to calculate molecular weight of given specie
function[mol_wt] = fre(x)
name=(x);
specie=["O", "H","C","N","A"];
mass=[16,1,12,14,40];
mol_wt=0;
for i= 1:length(name);
a= name(i);
q = str2num(a);
c=length(q);
if c == 0;
for j= 1:length(specie) ;
v =strcmp(a,specie(j));
if v == 1;
mol_wt= mol_wt+mass(j);
end
end
else
p = name(i-1);
for b= 1:length(specie) ;
u=strcmp(p,specie(b));
if u==1
rs = mass(b)*q;
mol_wt=mol_wt-mass(b);
mol_wt=mol_wt+rs;
end
end
end
end
end
%function to calculate entropy of given specie
function[S]=daset2(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,low_temp,high_temp,mid_temp)
R= 8.314;
lt = linspace(low_temp,high_temp,100);
for i=1:length(lt)
if lt(i)> mid_temp
S (i) = (a1*log(lt(i)) + a2*lt(i) + a3*(lt(i)^2)/2 + a4*(lt(i)^3)/3 + a5+(lt(i)^4)/4 + a7 )* R;
else
S (i) = (a8*log(lt(i)) + a9*lt(i) + a10*(lt(i)^2)/2 + a11*(lt(i)^3)/3 + a12+(lt(i)^4)/4 + a14 )* R;
end
end
end