Project 1 - Parsing NASA thermodynamic data

AIM:

To write a code in MATLAB to parse NASA thermodynamic data file and then calculate thermodynamic properties of various gas species.

OBJECTIVE:

• Write a function that extracts the 14 coefficients and calculates the enthalpy,entropy and specific heats for all the species in the data file.
• Calculating molecular weight of each species.
• Plotting the specific heat,enthalpy and entropy for the local temperature range.
• Saving the plots as images with appropriate names and dumping them in seperate folders for each species.

THEORY:

File Parsing :

File parsing is a string or symbol either in natural language, computer language or data structure confirming to the rules of the former grammer.

MATHEMATICAL EQUATION:

For specific heat (Cp),

CpR=a1+a2⋅T+a3⋅T2+a4⋅T3+a5⋅t4

For enthalpy (H),

HR=a1+a2⋅T2+a3⋅T23+a4⋅T34+a5⋅T45+a6T

For entropy (S),

SR=a1⋅ln(T)+a2.T+a3⋅T22+a4⋅T33+a5⋅T44+a7

where,

T=local temprature

a1 - a7 = high temperature coefficient

a8 - a14 = low temprature coefficient

R = Universal gas constant = 8.314

14 coefficients are obtained from the thermodynamic data file of which the first 7 are used in low temperature conditions and the latter 7 are used for high temperature conditions. The temperature for each species is mentioned in the data file with a maximum and minimum limit as well as a median limit. For temperatures above the median limit, the high temperature coefficients are used and for temperatures below the median limit, the low temperature coefficients are used.

MATLAB CODE :

MAIN CODE :

clear all
close all
clc

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

%Reading header information
fgetl(f1);

%Global temperature
t_line = fgetl(f1); 

a=strsplit(t_line,' ');

global_low_temp=str2double(a{2});
global_mid_temp=str2double(a{3});
global_high_temp=str2double(a{4});

%Neglecting the 3 comments lines
fgetl(f1);
fgetl(f1);
fgetl(f1);


for i=1:53
%Reading the species block

t_line2 = fgetl(f1);
s=strsplit(t_line2,' ');
species_name = s{1};

A = strfind(t_line2,'G');
temp_values = t_line2(A(1)+4:length(t_line2));
B = strsplit(temp_values,' ');

local_low_temp=str2double(B{1});
local_high_temp=str2double(B{2});
local_mid_temp=str2double(B{3});

%Line1
line1 = fgetl(f1);
b = strfind(line1,'E');
%1st coefficient
b1 = line1(1:b(1)+3);
b1 = str2double(b1);
%2nd coefficient
b2 = line1(b(1)+4:b(2)+3);
b2 = str2double(b2);
%3rd coefficient
b3 = line1(b(2)+4:b(3)+3);
b3 = str2double(b3);
%4th coefficient
b4 = line1(b(3)+4:b(4)+3);
b4 = str2double(b4);
%5th coefficient
b5 = line1(b(4)+4:b(5)+3);
b5 = str2double(b5);

%Line2
line2 = fgetl(f1);
%6th coefficient
b6 = line2(1:b(1)+3);
b6 = str2double(b6);
%7th coefficient
b7 = line2(b(1)+4:b(2)+3);
b7 = str2double(b7);
%8th coefficient
b8 = line2(b(2)+4:b(3)+3);
b8 = str2double(b8);
%9th coefficient
b9 = line2(b(3)+4:b(4)+3);
b9 = str2double(b9);
%10th coefficient
b10 = line2(b(4)+4:b(5)+3);
b10 = str2double(b10);

%Line3
line3 = fgetl(f1);
%11th coefficient
b11 = line3(1:b(1)+3);
b11 = str2double(b11);
%12th coefficient
b12 = line3(b(1)+4:b(2)+3);
b12 = str2double(b12);
%13th coefficient
b13 = line3(b(2)+4:b(3)+3);
b13 = str2double(b13);
%14th coefficient
b14 = line3(b(3)+4:b(4)+3);
b14 = str2double(b14);

%Universal gas constant
R = 8.314;

t = linspace(local_low_temp,local_high_temp,100);

% To calculate Specific heat
c_p = specific_heat(b1,b2,b3,b4,b5,b8,b9,b10,b11,b12,local_mid_temp,t,R);

%To calculate enthalpy
H = enthalpy(b1,b2,b3,b4,b5,b6,b8,b9,b10,b11,b12,b13,local_mid_temp,t,R);

%To calculate entropy
S = entropy(b1,b2,b3,b4,b5,b7,b8,b9,b10,b11,b12,b14,local_mid_temp,t,R);

%To calculate molecular weight
molecular_weight = molecular(species_name);
sprintf('Molecular Weight of %s = %d',species_name,molecular_weight)

%TO Create directory
file_name=sprintf('%s',species_name);
mkdir(file_name)
cd (file_name)

%Plotting

%Temperature vs Specific heat
figure(1)
plot(t,c_p,'linewidth',3)
xlabel('Temperature (K)')
ylabel('Specific heat (J/kg.K)')
title(species_name,'Temperature vs Specific heat')
saveas(gcf,'cp.png')

%Temperature vs Enthalpy
figure(2)
plot(t,H,'linewidth',3)
xlabel('Temperature (K)')
ylabel('Enthalpy (J/kg)')
title(species_name,'Temperature vs Enthalpy')
saveas(gcf,'H.png')

%Temperature vs Entropy
figure(3)
plot(t,S,'linewidth',3)
xlabel('Temperature (K)')
ylabel('Entropy (J/K)')
title(species_name,'Temperature vs Entropy')
saveas(gcf,'S.png')
cd ..
end 

FUNCTION CODE :

• For Molecular weight

function [out]=molecular(species_name)
%Defining the elements
elements=['H','C','N','O','R']; 

%Atomic weights of each element
Atomic_Weight = [1 12 14 16 40]; 
out=0; 

%Looping for each element inside each species_name and finding the atomic weight
for i=1:length(species_name) 
    for j=1:length(elements) 
        if strcmp(species_name(i),elements(j)) 
            out= out + Atomic_Weight(j);
            p=j; 
        end
    end
    
    %loop to multiply the coefficients of the elements if greater than 1.
    n=str2double(species_name(i)); 
    if n>1 
        out = out + (Atomic_Weight(p))*(n-1);
    end
end

end

• For Specific heat

function [c_p] = specific_heat(b1,b2,b3,b4,b5,b8,b9,b10,b11,b12,local_mid_temp,t,R)
if (t >local_mid_temp)
    c_p = R*(b1+b2*t+b3*(t.^2)+b4*(t.^3)+b5*(t.^4));
else
    c_p = R*(b8+b9*t+b10*(t.^2)+b11*(t.^3)+b12.*(t.^4));

end

• For Enthalpy

function [H] = enthalpy(b1,b2,b3,b4,b5,b6,b8,b9,b10,b11,b12,b13,local_mid_temp,t,R)
if (t >local_mid_temp)
    H = R*(b1+b2*(t/2)+b3*((t.^2)/3)+b4*((t.^3)/4)+b5*((t.^4)/5)+(b6*t));
else
    H = R*(b8+b9*(t/2)+b10*((t.^2)/3)+b11*((t.^3)/4)+b12*((t.^4)/5)+(b13*t));

end

• For Entropy

function [S] = entropy(b1,b2,b3,b4,b5,b7,b8,b9,b10,b11,b12,b14,local_mid_temp,t,R)
if (t >local_mid_temp)

 %Entropy
 S = R*((b1*log(t)) + (b2*t) + (b3*(t.^2)/2) + (b4*(t.^3)/3) + (b5*(t.^4)/4) + b7);

 else

 S = R*((b8*log(t)) + (b9*t) + (b10*(t.^2)/2) + (b11*(t.^3)/3) + (b12*(t.^4)/4) + b14);

end

end

RESULTS :

• Molecular weight for all species

• Screenshot of the window where all the folders are saved

• PLOTS :

For O2

For N2

For CO2

CONCLUSION :

• The characteristic curves for all the species are obtained by using MATLAB.
• Cp and S are inversely proportional to temperature they follow a parabolic path and as H is directly proportional to temperature it shows linear plot.