Project 1 - Parsing NASA thermodynamic data

Main Program:

clear all
close all
clc

%Input
R=8.314;
% Here R is Gas Constant

%open file and reading THERMO.dat file
f1=fopen('THERMO.dat', 'r')
get=fgetl(f1);

%Now Reading Temperature from this file
temp=fgetl(f1);

%spliting
s=strsplit(temp,' ');

%converting strings cells into numbers
templow=str2num(s{2});
tempmid=str2num(s{3});
temphigh=str2num(s{4});

%Temperature Range
T=linspace(templow,temphigh,1000);

%skipping comment lines
for i=1:3
get=fgetl(f1);
end

%writing molecular weight file
M_file=fopen('Molecular_weight.txt', 'w');


%Reading the Species
for j=1:53
tline=fgetl(f1);
%spliting into cell array
A=strsplit(tline,' ');
%reading the name of the species
species=(A{1});
%Converting string to number
lowtemp=str2num(A{end-3});
midtemp=str2num(A{end-2});
hightemp=str2num(A{end-1});


%getting next line
tline1=fgetl(f1);
%Finding the string till point E
a=strfind(tline1, 'E');

%Finding First 7 higher and second 7 lower coefficients

%Finding Higher coefficients

%Finding the Value on basis of position of E
%And Converting the string obtained to number
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);

%getting next line
tline2=fgetl(f1);
%Finding the string till point E
b=strfind(tline2, 'E');
a6=tline2(1:b(1)+3);
a6=str2num(a6);
a7=tline2(b(1)+4:b(2)+3);
a7=str2num(a7);

%Finding Lower coefficients

a8=tline2(b(2)+4:b(3)+3);
a8=str2num(a8);
a9=tline2(b(3)+4:b(4)+3);
a9=str2num(a9);
a10=tline2(b(4)+4:b(5)+3);
a10=str2num(a10);

%getting next line
tline3=fgetl(f1);
%Finding the string till point E
c=strfind(tline3, 'E');
a11=tline3(1:c(1)+3);
a11=str2num(a11);
a12=tline3(c(1)+4:c(2)+3);
a12=str2num(a12);
a13=tline3(c(2)+4:c(3)+3);
a13=str2num(a13);
a14=tline3(c(3)+4:c(4)+3);
a14=str2num(a14);


%calculation specific heat
Cp=Specific_heat(a1, a2, a3, a4, a5, a8, a9, a10, a11, a12,T,R, tempmid);

%calculation Entropy
S=Entropy(a1, a2, a3, a4, a5, a7, a8, a9, a10, a11, a12, a14,T,R, tempmid);

%claculation Enthalpy
H=Enthalpy(a1, a2, a3, a4, a5, a6, a8, a9, a10, a11, a12, a13,T,R, tempmid);

%Calculation Molecular weight of the Species
M(j)= Molecular_weight(species);

%changing directory
cd('/MATLAB Drive/Therm_data');
%creating directory
mkdir(species); 
%changing directory
cd(species);


%Writing on to file
fprintf(M_file, 'Molecular weight of %s is %d ',species,M(j));
fprintf(M_file,'\n');


%Plot1 T vs Cp
figure(1)
plot(T,Cp, 'linewidth', 2, 'color','r');
xlabel ('Temperature(K)');
ylabel('Specific Heat(KJ)');
title(sprintf('Specific heat vs Temperature range of %s', species));
grid on
%Saving Plot1
saveas (figure(1), 'specific Heat.jpg');

%Plot2 T vs S
figure(2)
plot (T,S, 'linewidth',2,'color','g');
xlabel('Temperature');
ylabel('Entropy(KJ/mol-K)');
title(sprintf('Entropy vs Temperature range of %s',species));
grid on
%Saving Plot2
saveas (figure(2), 'Entropy.jpg');

%Plot3 T vs H
figure(3)
plot (T,H, 'linewidth',2, 'color','b');
xlabel('Temperature');
ylabel('Enthalpy (KJ/mol)');
title(sprintf('Enthalpy vs Temperature range of %s', species));
grid on
%Saving Plot3
saveas(figure(3),'Enthalpy.jpg');

%Changing directory to Original matlab path
cd('/MATLAB Drive');


end

%Closing the file
fclose(M_file);

Function for calculation of Entropy:

function S=Entropy(a1, a2,a3,a4,a5,a7,a8, a9,a10, a11,a12, a14,T,R, tempmid)
   
   %calculating entropy based on Temperature.
    if (T <= tempmid)
       S=R*(a8.*log(T)+(a9*T)+((a10*(T).^2)/2)+((a11*(T).^3)/3)+((a12*(T).^4)/4)+a14);
    else
       S=R*(a1.*log(T)+(a2*T)+((a3*(T).^2)/2)+((a4*(T).^3)/3)+((a5*(T).^4)/4)+a7);
    
    end
end

Function for calculation of Enthalpy:

function H=Enthalpy(a1, a2, a3, a4, a5, a6, a8, a9, a10, a11, a12, a13,T,R, tempmid)
   
   %calculating enthalpy based on Temperature
   if (T <= tempmid)
       H=R*T.*(a8+((a9*T)/2)+((a10*(T).^2)/3)+((a11*(T).^3)/4)+((a12*(T).^4)/5)+((a13./T)));
   else
       H=R*T.*(a1+((a2*T)/2)+((a3*(T).^2)/3)+((a4*(T).^3)/4)+((a5*(T).^4)/5)+((a6./T)));
   end
end

Function for calculation of Molecular Weight:

function M=Molecular_weight(species)
%standard-Atomic Name
atomic_name=['H', 'C', 'O','N','Ar'];
%Hydrogen, Carbon, Oxygen, Nitrogen, Argon
%standard Atomic weights
atomic_weight=[1.008 12.011 15.999 14.007 39.948];
%starting Range of Molecule
M=0;
%@ for loops are used for comparison of each value og species to
for i=1:length(species)
for j=1:length(atomic_name)
if strcmp(species(i),atomic_name(j))
    M=M+atomic_weight(j);
    position=j;
 end
end

n=str2num(species(i));
if n>1
    M=M+atomic_weight(position)*(n-1);
end
end
fprintf('molecular weight %s :',species);
fprintf('%f',M);
disp(' ');
end

Function for calculation of Specific Heat:

function Cp=Specific_heat(a1,a2, a3, a4, a5, a8, a9, a10, a11, a12, T,R, tempmid)
   
   %calculating specific heat based on Temperature.
    if (T <= tempmid)
       Cp=R*(a8+(a9*T)+(a10*(T).^2)+(a11*(T).^3)+(a12*(T).^4));
   else
       Cp=R*(a1+(a2*T)+(a3*(T).^2)+(a4*(T).^3)+(a5*(T).^4));
   end
end

AIM:

Developing a program to extract 14 coefficient values for the corresponding species provided in a file Thermo.dat so as to calculate various value like Entropy, Enthalpy, Specific Heat and Molecular Weight and applying these values obtained to plot graphs with Temperature for Entropy, Enthalpy and Specific Heat values.

EQUATIONS:

Cp/R=al+a2T+a3 T^2+a4 T^3+a5 T^4

H/RT=al+a2T/2+a3 T^2 /3+a4 T^3 /4+a5 T^4 /5+a6/T

S/R=al InT+a2T+a3 T^2 /2+a4 T^3 /3+a5 T^4 /4+a7

Here R is Gas Constant = 8.314 J⋅K^(−1)⋅mol^(−1)

Cp is the Specific Heat, H is Enthalpy and S is Entropy Values.

STEPS:

1. Obtaining the Thermo.dat file and saving it in same directory as the MATLAB code
2. Development of Code to read the data from Thermo.dat File by applying necessary looping statements.
3. Development of Individual Functions for Calculation of Molecular Weight, Entropy, Enthalpy and Specific heats in accordance to the mathematical equations provided.
4. Changing directories and creating a folder using the species names so as to store their respective plots for Temperature vs Enthalpy, Entropy and Specific Heat.
5. Obtained data is Checked for accuracy

CODE Explanation:

_________________________________________________________________________

clear all

close all

clc

%Input

R=8.314;

% Here R is Gas Constant

%open file and reading THERMO.dat file

f1=fopen('THERMO.dat', 'r');

get=fgetl(f1);

%Now Reading Temperature from this file

temp=fgetl(f1);

%spliting

s=strsplit(temp,' ');

%converting strings cells into numbers

templow=str2num(s{2});

tempmid=str2num(s{3});

temphigh=str2num(s{4});

%Temperature Range

T=linspace(templow,temphigh,1000);

%skipping comment lines

for i=1:3

get=fgetl(f1);

end

%writing molecular weight file

M_file=fopen('Molecular_weight.txt', 'w');

%Reading the Species

for j=1:53

tline=fgetl(f1);

%spliting into cell array

A=strsplit(tline,' ');

%reading the name of the species

species=(A{1});

%Converting string to number

lowtemp=str2num(A{end-3});

midtemp=str2num(A{end-2});

hightemp=str2num(A{end-1});

%getting next line

tline1=fgetl(f1);

%Finding the string till point E

a=strfind(tline1, 'E');

%Finding First 7 higher and second 7 lower coefficients

%Finding Higher coefficients

%Finding the Value on basis of position of E

%And Converting the string obtained to number

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);

%getting next line

tline2=fgetl(f1);

%Finding the string till point E

b=strfind(tline2, 'E');

a6=tline2(1:b(1)+3);

a6=str2num(a6);

a7=tline2(b(1)+4:b(2)+3);

a7=str2num(a7);

%Finding Lower coefficients

a8=tline2(b(2)+4:b(3)+3);

a8=str2num(a8);

a9=tline2(b(3)+4:b(4)+3);

a9=str2num(a9);

a10=tline2(b(4)+4:b(5)+3);

a10=str2num(a10);

%getting next line

tline3=fgetl(f1);

%Finding the string till point E

c=strfind(tline3, 'E');

a11=tline3(1:c(1)+3);

a11=str2num(a11);

a12=tline3(c(1)+4:c(2)+3);

a12=str2num(a12);

a13=tline3(c(2)+4:c(3)+3);

a13=str2num(a13);

a14=tline3(c(3)+4:c(4)+3);

a14=str2num(a14);

%calculation specific heat

Cp=Specific_heat(a1, a2, a3, a4, a5, a8, a9, a10, a11, a12,T,R, tempmid);

%calculation Entropy

S=Entropy(a1, a2, a3, a4, a5, a7, a8, a9, a10, a11, a12, a14,T,R, tempmid);

%claculation Enthalpy

H=Enthalpy(a1, a2, a3, a4, a5, a6, a8, a9, a10, a11, a12, a13,T,R, tempmid);

%Calculation Molecular weight of the Species

M(j)= Molecular_weight(species);

%changing directory

cd('/MATLAB Drive/Therm_data');

%creating directory

mkdir(species);

%changing directory

cd(species);

%Writing on to file

fprintf(M_file, 'Molecular weight of %s is %d ',species,M(j));

fprintf(M_file,'\n');

%Plot1 T vs Cp

figure(1)

plot(T,Cp, 'linewidth', 2, 'color','r');

xlabel ('Temperature(K)');

ylabel('Specific Heat(KJ)');

title(sprintf('Specific heat vs Temperature range of %s', species));

grid on

%Saving Plot1

saveas (figure(1), 'specific Heat.jpg');

%Plot2 T vs S

figure(2)

plot (T,S, 'linewidth',2,'color','g');

xlabel('Temperature');

ylabel('Entropy(KJ/mol-K)');

title(sprintf('Entropy vs Temperature range of %s',species));

grid on

%Saving Plot2

saveas (figure(2), 'Entropy.jpg');

%Plot3 T vs H

figure(3)

plot (T,H, 'linewidth',2, 'color','b');

xlabel('Temperature');

ylabel('Enthalpy (KJ/mol)');

title(sprintf('Enthalpy vs Temperature range of %s', species));

grid on

%Saving Plot3

saveas(figure(3),'Enthalpy.jpg');

%Changing directory to Original matlab path

cd('/MATLAB Drive');

end

%Closing the file

fclose(M_file);

1. Default Commands for obtained cleared output windows like clc etc are provided.
2. Constant Value R is provided
3. f1=fopen('THERMO.dat', 'r');

get=fgetl(f1);

These command are used to store the file Thermo.dat into f1 and provide it with reading access and the fgetl command is used to obtain the first line of f1.

4. %spliting

s=strsplit(temp,' ');

%converting strings cells into numbers

templow=str2num(s{2});

tempmid=str2num(s{3});

temphigh=str2num(s{4});

These command are used to split the string rray into cell array and convert

Cell array into numbers.

5. %Temperature Range

T=linspace(templow,temphigh,1000);

%skipping comment lines

for i=1:3

get=fgetl(f1);

end

These commands are used for obtaining a temperature range and also skipping the comment lines in the file using for loop

6. %writing molecular weight file

M_file=fopen('Molecular_weight.txt', 'w');

%Reading the Species

for j=1:53

tline=fgetl(f1);

%spliting into cell array

A=strsplit(tline,' ');

%reading the name of the species

species=(A{1});

%Converting string to number

lowtemp=str2num(A{end-3});

midtemp=str2num(A{end-2});

hightemp=str2num(A{end-1});

%getting next line

tline1=fgetl(f1);

%Finding the string till point E

a=strfind(tline1, 'E');

These commands are used for extracting data using for loop. The command fgetl is used here for skipping the present line. Storing it on to a variable and splitting the according the blank spaces and determining the lower, higher and mid temperature values then the next line is spitted and the next line is stored in as tline1 and the value of the string are sorted out using by finding the presence of ‘E’ in the data and applying the logic accordingly. Molecular weight text file is opened before for loop so as to prevent the rewriting of file each time the loop runs.

7. %Finding First 7 higher and second 7 lower coefficients

%Finding Higher coefficients

%Finding the Value on basis of position of E

%And Converting the string obtained to number

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);

These commands use the logic obtained by determining the position of ‘E’ to Calculate the coefficients similar logic is applied for determining the entire coefficients.

8. %calculation specific heat

Cp=Specific_heat(a1, a2, a3, a4, a5, a8, a9, a10, a11, a12,T,R, tempmid);

%calculation Entropy

S=Entropy(a1, a2, a3, a4, a5, a7, a8, a9, a10, a11, a12, a14,T,R, tempmid);

%claculation Enthalpy

H=Enthalpy(a1, a2, a3, a4, a5, a6, a8, a9, a10, a11, a12, a13,T,R, tempmid);

%Calculation Molecular weight of the Species

M(j)= Molecular_weight(species);

These Commands are used for calling the respective functions so as to determine the Specific Heat, Entropy, Enthalpy and Molecular Weight.

9. %changing directory

cd('/MATLAB Drive/Therm_data');

%creating directory

mkdir(species);

%changing directory

cd(species);

%Writing on to file

fprintf(M_file, 'Molecular weight of %s is %d ',species,M(j));

fprintf(M_file,'\n');

%Plot1 T vs Cp

figure(1)

plot(T,Cp, 'linewidth', 2, 'color','r');

xlabel ('Temperature(K)');

ylabel('Specific Heat(KJ)');

title(sprintf('Specific heat vs Temperature range of %s', species));

grid on

%Saving Plot1

saveas (figure(1), 'specific Heat.jpg');

%Plot2 T vs S

figure(2)

plot (T,S, 'linewidth',2,'color','g');

xlabel('Temperature');

ylabel('Entropy(KJ/mol-K)');

title(sprintf('Entropy vs Temperature range of %s',species));

grid on

%Saving Plot2

saveas (figure(2), 'Entropy.jpg');

%Plot3 T vs H

figure(3)

plot (T,H, 'linewidth',2, 'color','b');

xlabel('Temperature');

ylabel('Enthalpy (KJ/mol)');

title(sprintf('Enthalpy vs Temperature range of %s', species));

grid on

%Saving Plot3

saveas(figure(3),'Enthalpy.jpg');

%Changing directory to Original matlab path

cd('/MATLAB Drive');

end

These Commands are used for writing data of molecular weight into a text file and also to change directories using ’Cd’ command to generate folders in name of species and store their plots of Temperature vs Enthalpy, Entropy and Specific Heat as the for loop runs.

10. Function for calculation of Entropy:

function S=Entropy (a1, a2, a3, a4, a5, a7, a8, a9, a10, a11, a12, a14, T,R, tempmid)

%calculating entropy based on Temperature.

if (T <= tempmid)

S=R* (a8.*log(T)+(a9*T)+((a10* (T).^2)/2)+((a11*(T).^3)/3)+((a12*(T).^4)/4)+a14);

else

S=R* (a1.*log(T)+(a2*T)+((a3*(T).^2)/2)+((a4*(T).^3)/3)+((a5*(T).^4)/4)+a7);

end

end

This Function returns Entropy Value on basis of the equation: “S/R=al InT+a2T+a3 T^2 /2+a4 T^3 /3+a5 T^4 /4+a7”

11. Function for calculation of Enthalpy:

function H=Enthalpy(a1, a2, a3, a4, a5, a6, a8, a9, a10, a11, a12, a13,T,R, tempmid)

%calculating enthalpy based on Temperature

if (T <= tempmid)

H=R*T.*(a8+((a9*T)/2)+((a10*(T).^2)/3)+((a11*(T).^3)/4)+((a12* (T).^4)/5)+((a13./T)));

else

H=R*T.*(a1+((a2*T)/2)+((a3*(T).^2)/3)+((a4*(T).^3)/4)+((a5*(T).^4)/5)+((a6./T)));

end

end

This Function returns Enthalpy Value on basis of the equation: “H/RT=al+a2T/2+a3 T^2 /3+a4 T^3 /4+a5 T^4 /5+a6/T”

12. Function for calculation of Specific Heat:

function Cp=Specific_heat(a1, a2, a3, a4, a5, a8, a9, a10, a11, a12, T,R, tempmid)

%calculating specific heat based on Temperature.

if (T <= tempmid)

Cp=R* (a8+(a9*T)+(a10*(T).^2)+(a11*(T).^3)+(a12*(T).^4));

else

Cp=R* (a1+(a2*T)+(a3*(T).^2)+(a4*(T).^3)+(a5*(T).^4));

end

end

This Function returns Enthalpy Value on basis of the equation: “Cp/R=al+a2T+a3 T^2+a4 T^3+a5 T^4”

13. Function for calculation of Molecular Weight.

function M=Molecular_weight(species)

%standard-Atomic Name

atomic_name=['H', 'C', 'O','N','Ar'];

%Hydrogen, Carbon, Oxygen, Nitrogen, Argon

%standard Atomic weights

atomic_weight=[1.008 12.011 15.999 14.007 39.948];

%starting Range of Molecule

M=0;

%@ for loops are used for comparison of each value og species to

for i=1:length(species)

for j=1:length(atomic_name)

if strcmp(species(i),atomic_name(j))

M=M+atomic_weight(j);

position=j;

end

end

n=str2num(species(i));

if n>1

M=M+atomic_weight(position)*(n-1);

end

end

fprintf('molecular weight %s :',species);

fprintf('%f',M);

disp(' ');

end

This Function Uses 2 for loops to compare the species character by character with the atomic names and used the value from the atomic weight table to along with the str2num converter to sum and multiply the elements according to their chemical formula thus providing the total molecular weight. This result is added to the Molecular_weight text file.

RESULT:

Applied the data read from Thermo.dat file in determining the constants which were applied on to the provided equations in the problem statements to calculate specie specific Entropy, Enthalpy, Specific heat and Molecular weight values and Temperature vs Enthalpy, Temperature vs Entropy and Temperature vs specific heat graphs were plotted and saved onto folders by utilizing the data provided.

MATLAB Drive Link (In case of requirement):

https://drive.matlab.com/sharing/fb3b1e76-d88b-4fe3-a8b8-62fcb0a456cc

Errors:

1. Typo errors.
2. File not getting written text file used to store molecular weight as required

Error Resolutions:

1. Typo errors were rectified by referring to syntaxes.
2. Rectified and modified the logic and other errors in the code.

I have attached the Screenshots for your reference.

Resources:

1. https://in.mathworks.com/help/matlab/ref/saveas.html
2. https://in.mathworks.com/help/matlab/import_export/writing-to-text-data-files-with-low-level-io.html
3. https://in.mathworks.com/help/matlab/ref/cd.html#:~:text=Use%20the%20cd%20command%20to%20display%20the%20new%20current%20folder.&text=Change%20the%20current%20folder%20back,display%20the%20new%20current%20folder.

# Have done it using MATLAB Web and there were certain limitations but have done it to my best also detailed code explanations have been already been provided with as comments with the code.