Challenge

Project 1 - Parsing NASA thermodynamic data

VVedant Meher7 Jun 2022

File Parsing NASA thermodynamic data

NASA come up with polynomials that can be used to evaluate thermodynamic properties such as Cp, H and S. they have also documented the co-efficient that are required to evaluated these polynomials for 1000+ species.

Here we are PARSING 53 species from NASA thermodynamic data file

THEORY

Specific heat(cp) : specific heat is defined by the amount of heat needed to raise the temperature of 1 gram of a substance 1 degree Celsius (). Water has a high Specific heat capacity which we’ll refer to as simply “heat capacity”, meaning it take more energy to increase the temperature of water compared to other substances.

Q = heat energy

M = mass

C = specific heat capacity

Delta T = change in temperature

Entropy(s): Entropy, the measure of system’s thermal energy of unit temperature that is unavailable for doing useful work. Because work is obtained from ordered molecular motion, the amount of entropy is also measure of the molecular disorder, or randomness of a system.

S = entropy

Kb =Boltzmann constant

ln =natural logarithm

Entropy(H): Entropy is a property of a thermodynamic system, define as the system’s internal energy and the product of its pressure and volume. It is a convenient state function standardly used in many measurements in chemical, biological and physical systems at a constant pressure.

H = enthalpy

E = internal energy

P = pressure

V = volume

OBJECTIVE

Write a function that extracts the 114 co-efficient and calculates the enthalpy, entropy, and specific heats for all the specifies in the data file. You will be reading this file NASA thermodynamic data

Calculate the molecular weight of each species
Plot the p, Enthalpy and Entropy for the local temperature range (low temperature: high temperature) species for each species.
Save the plot as image with appropriate name and dump them in separate folders for each species with suitable name. all this should be generated automatically.
Plot for o2, N2 and co2 species

CODE

Main code

clear all
close all
clc

%opening of NASA's thermodynamic data

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

%Obtaining of global, max and medium temperatures

global_temp = fgetl(f1);
temp_range = strsplit(global_temp,' ');
global_min_temp = str2num(temp_range{2});
global_med_temp = str2num(temp_range{3});
global_max_temp = str2num(temp_range{4});

%Range of temperature
T = linspace(global_min_temp,global_max_temp,400);
R = 8.314;

%skipping the next three comments in the file
for i = 3:5
comments = fgetl(f1);
end

%Getting the coefficients of species
for i = 1:53
%Extracting line 1 and the name of the species
line1 = fgetl(f1);
name = strtok(line1);

%Extracting higher and lower temperature coefficients
line2 = fgetl(f1);
line3 = fgetl(f1);
line4 = fgetl(f1);

A = strfind(line2,'E');
B = strfind(line3,'E');
C = strfind(line4,'E');

%Extracting higher temperature coefficients
a1 = str2num(line2(1:A(1)+3));
a2 = str2num(line2(A(1)+4:A(2)+3));
a3 = str2num(line2(A(2)+4:A(3)+3));
a4 = str2num(line2(A(3)+4:A(4)+3));
a5 = str2num(line2(A(4)+4:A(5)+3));
a6 = str2num(line3(1:B(1)+3));
a7 = str2num(line2(B(1)+4:B(2)+3));

%Extracting lower temperature coefficients
b1 = str2num(line3(B(2)+4:B(3)+3));
b2 = str2num(line3(B(3)+4:B(4)+3));
b3 = str2num(line3(B(4)+4:B(5)+3));
b4 = str2num(line4(1:C(1)+3));
b5 = str2num(line4(C(1)+4:C(2)+3));
b6 = str2num(line4(C(2)+4:C(3)+3));
b7 = str2num(line4(C(3)+4:C(4)+3));

%specific heat calculation
cp = specificheat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,T,R,global_med_temp);

%enthalpy calculation
h = enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,T,R,global_med_temp);

%Entropy calculation
s = entropy_calc(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,T,R,global_med_temp);

%calculating molecular weight of species
Molecular_wt(i) = molecular_weight(name)

%calculating individual folder for each species and plotting
plotting(cp,h,s,T,name)

end

To calculate molecular weight

%function of calculate molecular weight
function MW = molecular_weight(name)
atoms = ['H','C','O','N','R'];
fprintf('Molecular Weight of ');
fprintf(name);
fprintf(' : ');
atomic_weight = [1 12 16 14 40];
MW = 0;
for i = 1:length(name)
for j = 1:length(atoms)
if strcmp(name(i),atoms(j))
MW = MW + atomic_weight(j);
position = j;
end
end
n = str2num(name(i));
if n>1
MW = MW + atomic_weight(position)*(n-1);
end
end
fprintf('%d',MW);
disp(' ');
end

Plotting

%creat plotting function
function plotting(cp,h,s,T,name)
cd(pwd)
dir_name = name;
mkdir(dir_name);

figure(1)
plot(T,cp,'linewidth',2,'color','b')
title(sprintf('Specific heat Vs Temperature for %s',name))
xlabel('TEMPERATURE')
ylabel('SPECIFIC HEAT')

figure(2)
plot(T,h,'linewidth',2,'color','r')
title(sprintf('Enthalpy Vs Temperature for %s',name))
xlabel('TEMPERATURE')
ylabel('ENTHALPY')

figure(3)
plot(T,s,'linewidth',2,'color','g')
title(sprintf('Entropy Vs Temperature for %s',name))
xlabel('TEMPERATURE')
ylabel('ENTROPY')

cd(dir_name)

saveas(1,'cp_vs_t.png')
saveas(2,'h_vs_t.png')
saveas(3,'s_vs_t.png')

cd ..

end

Calculate entropy

%Calculation of entropy
function s = entropy_calc(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,T,R,global_med_temp)
if T<=global_med_temp
s = R*(a1.*log(T)+a2.*T+a3.*(T.^2)./2+a4.*(T.^3)./3+a5.*(T.^4)./4+a7);
else
s = R*(b1.*log(T)+b2.*T+b3.*(T.^2)./2+b4.*(T.^3)./3+b5.*(T.^4)./4+b7);
end

Calculate enthalpy

%Calculation of enthalpy
function h = enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,T,R,global_med_temp)
if T >global_med_temp
h = R.*(a1.*T+(a2.*T.^2)./2+(a3.*T.^3)./3+(a4.*T.^4)./4+(a5.*T.^5)./5+a6);
else
h = R.*(b1.*T+(b2.*T.^2)./2+(b3.*T.^3)./3+(b4.*T.^4)./4+(b5.*T.^5)./5+b6);

end

Calculate specific heat

%calculation of specific heat
function cp = specificheat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,T,R,global_med_temp)
%comparing with global temp
%If less then global it will take min coefficients
if T > global_med_temp
cp = R.*(a1+a2.*T+a3.*T.^2+a4.*T.^3+a5.*T.^4);
%If temp is more it will take max coefficients
else
cp = R.*(b1+b2.*T+b3.*T.^2+b4.*T.^3+b5.*T.^4);
end

CODE EXPLINATION

Starting with clear all, close all and clc
By using fopen command we are opening the file of THERMO.dat and ‘r’ represents to read the file

fgetl: is used to read the line in the file of which we have given and for every fgetl it will be reading the next line in the file so the fgetl is used to read the line in the file

By using fgetl we will be reading the first line of the file which is THERMO and for next fgetl we will be getting the second line which is the Global_temperature and we are assigning to the G_temp

global_temp = fgetl(f1);

And we are skipping the next three lines by creating the loop because those are the comments

And then in the 4^TH line we have species name in the first position and then the local_temperature and we are assigning this line to the line
By using the strsplit we are splitting the line with the reference of ‘space’

Strsplit: split string or character vector at specified delimiter

temp_range = strsplit(global_temp,' ');

Delimiting characters, specified as a character vector, a 1-by-n cell array of character vectors, or a 1-by-n sting array. Text specified in delimiter dose not appear in the output c.

Specify multiple delimiters in a cell array or a string array. The strsplit function splits str on the elements of delimiter. The order in which delimiters appear in delimiter dose not matter unless multiple delimiters begin a match at the same character in str. In that case strsplit splits on the first matching delimiter in delimiter.

Strsplit is used to split the line and it is assigned to the variable a so the A (1) is the species name and for the local_temp we will be taking from the end
global_min_temp = str2num(temp_range{2});
global_med_temp = str2num(temp_range{3});
global_max_temp = str2num(temp_range{4});

str2num: Convert character array or string to numeric array

And then we have to read the co-efficient which are from the next line and then by using fgetl read the next line and assigned to a variable
To assign each value to a variable we used findstr

Findstr: find strings within other strings

We will be finding the E in the string (line) and by refiring that we will be taking the separate value and assign to a variable as shown
a_line = fgetl(f1);
a = findstr(a_line,'E');
a1 = str2num(line2(1:A(1)+3));

Here a will be getting the position of the E in the line and by using the position and number of numbers to left and right we will be getting the co-efficient and for the second number we have to start from the first number position as shown

a2 = str2num(line2(A(1)+4:A(2)+3));

lly to all co-efficient as shown

A = strfind(line2,'E');
B = strfind(line3,'E');
C = strfind(line4,'E');

%Extracting higher temperature coefficients
a1 = str2num(line2(1:A(1)+3));
a2 = str2num(line2(A(1)+4:A(2)+3));
a3 = str2num(line2(A(2)+4:A(3)+3));
a4 = str2num(line2(A(3)+4:A(4)+3));
a5 = str2num(line2(A(4)+4:A(5)+3));
a6 = str2num(line3(1:B(1)+3));
a7 = str2num(line2(B(1)+4:B(2)+3));

%Extracting lower temperature coefficients
b1 = str2num(line3(B(2)+4:B(3)+3));
b2 = str2num(line3(B(3)+4:B(4)+3));
b3 = str2num(line3(B(4)+4:B(5)+3));
b4 = str2num(line4(1:C(1)+3));
b5 = str2num(line4(C(1)+4:C(2)+3));
b6 = str2num(line4(C(2)+4:C(3)+3));
b7 = str2num(line4(C(3)+4:C(4)+3));

We have to write a formula for low and high temperature the FIRST 7 coefficients are HIGH-temperature coefficients and the SECOND 7 coefficients are LOW-temperature coefficients.

%Specific heat, Enthalpy, Entropy @ High temperature

%specific heat calculation
cp = specificheat(a1,a2,a3,a4,a5,b1,b2,b3,b4,b5,T,R,global_med_temp);

%enthalpy calculation
h = enthalpy(a1,a2,a3,a4,a5,a6,b1,b2,b3,b4,b5,b6,T,R,global_med_temp);

%Entropy calculation
s = entropy_calc(a1,a2,a3,a4,a5,a7,b1,b2,b3,b4,b5,b7,T,R,global_med_temp);

And then we have to create a folder to save the plots in their respective folders to do that we will be using the mkdir

mkdir: to create a folder

cd: change current folder

%creat plotting function
function plotting(cp,h,s,T,name)
cd(pwd)
dir_name = name;
mkdir(dir_name);

Plotting the enthalpy, entropy and specific heat for a low and high temperature coefficients

figure(1)
plot(T,cp,'linewidth',2,'color','b')
title(sprintf('Specific heat Vs Temperature for %s',name))
xlabel('TEMPERATURE')
ylabel('SPECIFIC HEAT')

figure(2)
plot(T,h,'linewidth',2,'color','r')
title(sprintf('Enthalpy Vs Temperature for %s',name))
xlabel('TEMPERATURE')
ylabel('ENTHALPY')

figure(3)
plot(T,s,'linewidth',2,'color','g')
title(sprintf('Entropy Vs Temperature for %s',name))
xlabel('TEMPERATURE')
ylabel('ENTROPY')

cd(dir_name)

saveas(1,'cp_vs_t.png')
saveas(2,'h_vs_t.png')
saveas(3,'s_vs_t.png')

cd ..

After every plot saveas will save the plot as fig in the respective folder

saveas: Save figure to specific file format gcf is the handle to save the plot image

And we have to give the .. command to come out from the particular folder

cd..: change current folder

A path relative to the current folder. cd assumes that the path is relative to the current folder even if the ./ is not present.

Calculating the molecular weight

function MW = molecular_weight(name)
atoms = ['H','C','O','N','R'];
fprintf('Molecular Weight of ');
fprintf(name);
fprintf(' : ');
atomic_weight = [1 12 16 14 40];
MW = 0;

for i = 1:length(name)
for j = 1:length(atoms)
if strcmp(name(i),atoms(j))
MW = MW + atomic_weight(j);
position = j;
end
end
n = str2num(name(i));
if n>1
MW = MW + atomic_weight(position)*(n-1);
end
end
fprintf('%d',MW);
disp(' ');
end

Create two strings of atom symbol and atomic weight as shown
atoms = ['H','C','O','N','R'];
atomic_weight = [1 12 16 14 40];
Creating a for loop to the species convert the species to number by using str2num command with 2 results
In species of it is number status will be 1 it is a letter(element) the status will be 0
By applying if condition status == 0
ind = find(elements == species(k)) we are finding the elements in species in the element array by using find command
Then the molecular_we = atomic_weight(ind) else
molecular_wt(k-1) so that the previous position to the number in the species will be changed
and then display the molecular_wt

RESULTS:

workspace

molecular weight of each species display it in the command window.

save the plot as images with appropriate names and dump them in separate folders for each species with suitable name is done AUTOMATICALLY
plot for o2,n2 and co2 species

PLOTS OF AR