Challenge

Project 1 - Parsing NASA thermodynamic data in MATLAB

AAnurag Mohanta21 Sept 2020

Objective

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.

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

CODE

close all
clear all
clc

f1=fopen('THERMO.dat','r');       %opening the file for reading
header=fgetl(f1);                 %Reading Header

globaltemp=fgetl(f1);             %reading second line
temp=strsplit(globaltemp,'  ');   %spliting string from the line
t=str2double(temp);               %converting string into number
global_low_temp=t(2);             %Assinging idices of 't' to a variable
global_mid_temp=t(3);
global_high_temp=t(4);

 cd('E:skill-lynchmatlab_projectsNASAfileparsing')
 fileID=fopen('molecular weight.txt','a');
 fwrite(fileID,sprintf('%s %6s','species_name','mol_wt'));
 fwrite(fileID,sprintf('n'));


%Reading the comment lines with for loop
for i=3:5
    fgetl(f1);
end

%SPECIES INFORMATION
for i=1:53
    species=fgetl(f1);
    species_name=species(1:6);
    
    line_1=fgetl(f1);
    a1=strfind(line_1,'E');
    c1=str2num(line_1(1:(a1(1)+3)));
    c2=str2num(line_1(1:(a1(1)+4):(a1(2)+3)));
    c3=str2num(line_1(1:(a1(2)+4):(a1(3)+3)));
    c4=str2num(line_1(1:(a1(3)+4):(a1(4)+3)));
    c5=str2num(line_1(1:(a1(4)+4):(a1(5)+3)));

    
    
    line_2=fgetl(f1);
    a2=strfind(line_2,'E');
    c6=str2num(line_2(1:(a2(1)+3)));
    c7=str2num(line_2(1:(a2(1)+4):(a2(2)+3)));
    c8=str2num(line_2(1:(a2(2)+4):(a2(3)+3)));
    c9=str2num(line_2(1:(a2(3)+4):(a2(4)+3)));
    c10=str2num(line_2(1:(a2(4)+4):(a2(5)+3)));
    
    
    line_3=fgetl(f1);
    a3=strfind(line_3,'E');
    c11=str2num(line_3(1:(a3(1)+3)));
    c12=str2num(line_3(1:(a3(1)+4):(a3(2)+3)));
    c13=str2num(line_3(1:(a3(2)+4):(a3(3)+3)));
    c14=str2num(line_3(1:(a3(3)+4):(a3(4)+3)));
  
    
    %creating folder for each species
    mkdir('E:skill-lynchmatlab_projectsNASAfileparsing',species_name);
    
    %Temperature ranges assumed for calculation purpose
    T=linspace(300,5000,100);
    
    %universial gas constant
    R=8.314;
    
    %calculation of properties based on lower and upper temperature
    for j=1:length(T)
        if T(j)<1000
            cp=(c1+c2.*T+c3.*(T.^2)+c4.*(T.^3)+c5.*(T.^4)).*R;
            H=(c1+(c2.*T)./2+(c3.*(T.^2))./3+(c4.*(T.^3))./4+(c5.*(T.^4))./5+c6./T).*R.*T;
            S=(c1.*log(T)+c2.*T+(c3.*(T.^2))./2+(c4.*(T.^3))./3+c5.*(T.^4)+c7).*R;
        else
            cp=(c8+c9.*T+c10.*(T.^2)+c11.*(T.^3)+c12.*(T.^4)).*R;
            H=(c8+(c9.*T)./2+(c10.*(T.^2))./3+(c11.*(T.^3))./4+(c12.*(T.^4))./5+c13./T).*R.*T;
            S=(c8.*log(T)+c9.*T+(c10.*(T.^2))./2+(c11.*(T.^3))./3+c12.*(T.^4)+c14).*R;
        end
    end
    
    %plotting
    %temperature Vs Entropy
    figure(1)
    plot(T,S,'linewidth',4,'color','b')
    title({'Entropy Vs Temperature','species name= ',species_name})
    xlabel('Temperature')
    ylabel('Entropy')
    
    
    
    
    
    %temperature Vs Enthalpy
    figure(2)
    plot(T,H,'linewidth',4,'color','k')
    title({'Enthalpy Vs Temperature','species name= ',species_name})
    xlabel('Temperature')
    ylabel('Enthalpy')
    
    
    
    
    %temperature Vs Specific Heat
    figure(3)
    plot(T,cp,'linewidth',4,'color','g')
    title({' Temperature Vs Specific Heat','species name= ',species_name})
    xlabel('Temperature')
    ylabel('Specific Heat')
    
    
    directory=(['E:skill-lynchmatlab_projectsNASAfileparsing',species_name]);
    cd(directory)
    
    saveas(figure(1),'Temperature Vs Entropy.png');
    saveas(figure(2),'Temperature Vs Enthaly.png');
    saveas(figure(3),'Temperature Vs Specific_Heat.png');
    
    
    cd('E:skill-lynchmatlab_projectsNASAfileparsing')
    mol_wt=molecular_mass(species_name)
    fileID=fopen('molecular weight.txt','a');
    fwrite(fileID,sprintf('%s %6dn',species_name,mol_wt));
    fclose(fileID);

    
    
end
fclose(f1);

FUNCTION

%molecular weight calculation
function mol_wt = molecular_mass(species_name)
species=upper(species_name)
species_elements={'O','H','C','N','S','A'};%Declaring Species elements
atomic_mass= [16,1,12,14,32,40];%Declaring the atomic mass value
mol_wt=0;%declaring molecular weight is '0' in initial stage

for i=1:length(species)
    for u=1:length(species_elements)
        if strcmp(species(i),species_elements(u))
            mol_wt=mol_wt+atomic_mass(u);
            y=u;%position is stored
        end
    end
    n=str2double(species(i));
    if n>1
        mol_wt=mol_wt+atomic_mass(y)*(n-1);
    end
end
end

CODE EXPLANATION

close all
clear all
clc

f1=fopen('THERMO.dat','r'); %opening the file for reading
header=fgetl(f1); %Reading Header

globaltemp=fgetl(f1); %reading second line
temp=strsplit(globaltemp,' '); %spliting string from the line
t=str2double(temp); %converting string into number
global_low_temp=t(2); %Assinging idices of 't' to a variable
global_mid_temp=t(3);
global_high_temp=t(4);

cd('E:skill-lynchmatlab_projectsNASAfileparsing') %go to the file directory
fileID=fopen('molecular weight.txt','a'); %Open or create new file for writing. Append data to the end of the file.
fwrite(fileID,sprintf('%s %6s','species_name','mol_wt')); %writing first line as string
fwrite(fileID,sprintf('n')); %crating new line

%Reading the comment lines with for loop
for i=3:5 %loop to read comment
fgetl(f1);
end

%SPECIES INFORMATION
for i=1:53 %loop to read elements
species=fgetl(f1); %reading sixth line which is first element name
species_name=species(1:6); %after that six numbers are cofficients of euation

line_1=fgetl(f1); %fetching data from single line
a1=strfind(line_1,'E'); %finding character "E"
c1=str2num(line_1(1:(a1(1)+3))); %converting string to number
c2=str2num(line_1(1:(a1(1)+4):(a1(2)+3)));
c3=str2num(line_1(1:(a1(2)+4):(a1(3)+3)));
c4=str2num(line_1(1:(a1(3)+4):(a1(4)+3)));
c5=str2num(line_1(1:(a1(4)+4):(a1(5)+3)));

line_2=fgetl(f1); %fetching data from single line
a2=strfind(line_2,'E'); %finding character "E"
c6=str2num(line_2(1:(a2(1)+3))); %converting string to number
c7=str2num(line_2(1:(a2(1)+4):(a2(2)+3)));
c8=str2num(line_2(1:(a2(2)+4):(a2(3)+3)));
c9=str2num(line_2(1:(a2(3)+4):(a2(4)+3)));
c10=str2num(line_2(1:(a2(4)+4):(a2(5)+3)));

line_3=fgetl(f1); %fetching data from single line
a3=strfind(line_3,'E'); %finding character "E"
c11=str2num(line_3(1:(a3(1)+3))); %converting string to number
c12=str2num(line_3(1:(a3(1)+4):(a3(2)+3)));
c13=str2num(line_3(1:(a3(2)+4):(a3(3)+3)));
c14=str2num(line_3(1:(a3(3)+4):(a3(4)+3)));

%creating folder for each species
mkdir('E:skill-lynchmatlab_projectsNASAfileparsing',species_name);

%Temperature ranges assumed for calculation purpose
T=linspace(300,5000,100);

%universial gas constant
R=8.314;

%calculation of properties based on lower and upper temperature
for j=1:length(T)
if T(j)<1000
cp=(c1+c2.*T+c3.*(T.^2)+c4.*(T.^3)+c5.*(T.^4)).*R;
H=(c1+(c2.*T)./2+(c3.*(T.^2))./3+(c4.*(T.^3))./4+(c5.*(T.^4))./5+c6./T).*R.*T;
S=(c1.*log(T)+c2.*T+(c3.*(T.^2))./2+(c4.*(T.^3))./3+c5.*(T.^4)+c7).*R;
else
cp=(c8+c9.*T+c10.*(T.^2)+c11.*(T.^3)+c12.*(T.^4)).*R;
H=(c8+(c9.*T)./2+(c10.*(T.^2))./3+(c11.*(T.^3))./4+(c12.*(T.^4))./5+c13./T).*R.*T;
S=(c8.*log(T)+c9.*T+(c10.*(T.^2))./2+(c11.*(T.^3))./3+c12.*(T.^4)+c14).*R;
end
end

%plotting
%temperature Vs Entropy
figure(1)
plot(T,S,'linewidth',4,'color','b')
title({'Entropy Vs Temperature','species name= ',species_name})
xlabel('Temperature')
ylabel('Entropy')

%temperature Vs Enthalpy
figure(2)
plot(T,H,'linewidth',4,'color','k')
title({'Enthalpy Vs Temperature','species name= ',species_name})
xlabel('Temperature')
ylabel('Enthalpy')

%temperature Vs Specific Heat
figure(3)
plot(T,cp,'linewidth',4,'color','g')
title({' Temperature Vs Specific Heat','species name= ',species_name})
xlabel('Temperature')
ylabel('Specific Heat')

directory=(['E:skill-lynchmatlab_projectsNASAfileparsing',species_name]);%going to this directory
cd(directory)

%saving figures
saveas(figure(1),'Temperature Vs Entropy.png');
saveas(figure(2),'Temperature Vs Enthaly.png');
saveas(figure(3),'Temperature Vs Specific_Heat.png');

cd('E:skill-lynchmatlab_projectsNASAfileparsing')
mol_wt=molecular_mass(species_name) %calling function molecularmass to calucalte molecular weight
fileID=fopen('molecular weight.txt','a'); %opening file molecular weight.txt as read write append format
fwrite(fileID,sprintf('%s %6dn',species_name,mol_wt)); %writing each species name and its molecular weight
fclose(fileID);%close file molecular weight.txt