Project 1 - Parsing NASA thermodynamic data

Function Code for Specific heat.

function c_p = specificheat_fn(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,t,R,med_temp)
%we compare the temperature t with medtemp and use if statement to proceed
if t<=med_temp
    %second half ie the low temp coeffs are taken
    c_p=R*(a8+(a9*t)+(a10*(t).^2)+(a11*(t).^4)+(a12*(t).^5));
else
    c_p=R*(a1+(a2*t)+(a3*(t).^2)+(a4*(t).^4)+(a5*(t).^5));
    %first half ie the high temp coeffs are taken
end

Function Code for Enthalpy.

function H=enthalpy_fn(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,t,R,med_temp)
%we compare the temperature t with medtemp and use if statement to proceed
if t <= med_temp
    %second half ie the low temp coeffs are taken
    H=R*t.*(a8+((a9.*t)/2)+((a10.*(t).^2)/3)+((a11.*(t).^3)/4)+((a12.*(t).^4)/5)+(a13./t));
else
   %first half ie the high temp coeffs are taken  
    H=R*t.*(a1+((a2.*t)/2)+((a3.*(t).^2)/3)+((a4.*(t).^3)/4)+((a5.*(t).^4)/5)+(a6./t));
end

Function code for Entropy.

function S = entropy_fn(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,t,R,med_temp)
%we compare the temperature t with medtemp and use if statement to proceed
if t <= med_temp
    %second half ie the low temp coeffs are taken
    S = R.*(a8*log(t)+(a9.*t)+(a10.*((t).^2)/2)+(a11.*((t).^2)/2)+(a12.*((t).^2)/2)+a14);
else
    %first half ie the high temp coeffs are taken
    S = R.*(a1*log(t)+(a2.*t)+(a3.*((t).^2)/2)+(a4.*((t).^2)/2)+(a5.*((t).^2)/2)+a7);
end

Main code for Parsing NASA Thermodynamic Data

clear all
close all 
clc

%Reading the thermo dat file and creating file marker for future use
f1=fopen('THERMO.dat','r');
UNW=fgetl(f1);
%reading header for 't' value
global_temp=fgetl(f1);

%string splitting and assigning temp values to variables for comparison
%later in the code
y=strsplit(global_temp,' ');
%Classifying temp ranges as low med and high to use to calc cp,H,S
low_temp= str2num(y{2});
med_temp= str2num(y{3});
high_temp=str2num(y{4});

%inputs
t=linspace(low_temp,high_temp,1000);
R = 8.314 ;                               %Universal Gas constant

%reading and skippinf 3 unnecessary lines
for i=1:3
    UNW=fgetl(f1);
end
 for i = 1:53                  %53 refers to the no of species for which the below code will repeat
    
    UNW=fgetl(f1);
    X=strsplit(UNW,' ');
    species_name=(X{1});
    species=char(species_name);
    %Reading coefficients and assigning them
    line1=fgetl(f1);
    a=strfind(line1,'E');
    
    
    %Extracting the coefficients from the body of the species
    a1= line1(1:a(1)+3);
    a1=str2num(a1);
    a2=line1(a(1)+4:a(2)+3);
    a2=str2num(a2);
    a3=line1(a(2)+4:a(3)+3);
    a3=str2num(a3);
    a4=line1(a(3)+4:a(4)+3);
    a4=str2num(a4);
    a5=line1(a(4)+4:a(5)+3);
    a5=str2num(a5);
    %Moving on to reading the next line
    line2=fgetl(f1);
    b=strfind(line2,'E');
    
    a6=line2(1:b(1)+3);
    a6=str2num(a6);
    a7=line2(b(1)+4:b(2)+3);
    a7=str2num(a7);
    a8=line2(b(2)+4:b(3)+3);
    a8=str2num(a8);
    a9=line2(b(3)+4:b(4)+3);
    a9=str2num(a9);
    a10=line2(b(4)+4:b(5)+3);
    a10=str2num(a10);
    
    %Moving on to reading line 3
    line3=fgetl(f1);
    c=strfind(line3,'E');
    a11=line3(1:c(1)+3);
    a11=str2num(a11);
    a12=line3(c(1)+4:c(2)+3);
    a12=str2num(a12);
    a13=line3(c(2)+4:c(3)+3);
    a13=str2num(a13);
    a14=line3(c(3)+4:c(4)+3);
    a14=str2num(a14);
    
    %CP,H,S,Molecular weight calculations by calling in the function
    c_p = specificheat_fn(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,t,R,med_temp);
    H = enthalpy_fn(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,t,R,med_temp);
    S = entropy_fn(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,t,R,med_temp);
    Mw = molecularweight_fn(species_name);
    w=sprintf('%s',species);
    mkdir(w)
    cd(w)
    %plotting
    %plot 1 (Specific heat vs temp)
    figure (1);
    plot(t,c_p,'linewidth',3,'color','r');
    title ('Specific heat vs Temperature');
    xlabel('temperature');
    ylabel('Specific heat');
    grid on
    saveas(1,'c_p.jpg')
    
    %plot 2 (Enthalpy vs temperature)
    figure (2);
    plot (t,H,'linewidth',3,'color','g');
    title ('Enthalpy vs Temperature');
    xlabel('temperature');
    ylabel('Enthalpy');
    grid on
    saveas(2,'H.jpg')
    
    
    %plot 3 (Entropy vs temperature)
    figure (3)
    plot(t,S,'linewidth',3,'color','b');
    title ('Entropy vs Temperature');
    xlabel('temperature');
    ylabel('Entropy');
    grid on
    saveas(3,'S.jpg')
    cd ..
 end
fclose(f1);

Function code to calculate molecular weight of species

function Mw = molecularweight_fn(species_name)
Mw=0;
elements= ['C','H','O','N','AR'];
Aw= [12 1 16 14 40];

for i=length(species_name)
    for j= length(elements)
     if   strcmpi(species_name(i),elements(j))
        Mw=Mw+Aw(j);
        p = j;
     end
     
    end
    n=str2num(species_name(i));


if  n>1
    Mw=Mw+(Aw(p)*(n-1));
end 
end

sprintf('molecular weight %s',species_name);
sprintf('%d',Mw);
end

The plots of O2,CO2 and N2 as wellas the folders containing plots screenshot is uploaded with the report.