Project 1 - Parsing NASA thermodynamic data

Main Program

clear all
close all
clc

% Universal gas constant (J/mol-K)

R = 8.314 ;
number_of_species = 53;

% Opening file from directory

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

fgetl(f1) % Skipping the first line.

tline = fgetl(f1); % Storing the second line in a variable tline

A = strsplit(tline,' ');

global_low_temp = str2num(A{2});

global_mid_temp = str2num(A{3});

global_high_temp = str2num(A{4});

global_temp_array = linspace(global_low_temp, global_high_temp,50);

% Skiping next three lines
fgetl(f1)
fgetl(f1)
fgetl(f1)

% Collecting information from each line for the next 53 times.

for i = 1:number_of_species
species_details = fgetl(f1);
species_name = strtok(species_details); % Name of each species
temp_find = findstr(species_details,'.');
local_temp_1 = str2double(species_details((temp_find(1)-3):(temp_find(1)+3)));
local_temp_2 = str2double(species_details((temp_find(2)-4):(temp_find(2)+3)));
local_temp_3 = str2double(species_details((temp_find(3)-4):(temp_find(3)+3)));

% Creating an array of local temperature values

local_temp_array = linspace(local_temp_1,local_temp_3,50);

% Collecting coefficients of every species

cline1 = fgetl(f1);
A1 = findstr(cline1,'E');
coeff_l1_str1 = cline1(1:A1(1)+3);
coeff_l1_str2 = cline1(A1(1)+4:A1(2)+3);
coeff_l1_str3 = cline1(A1(2)+4:A1(3)+3);
coeff_l1_str4 = cline1(A1(3)+4:A1(4)+3);
coeff_l1_str5 = cline1(A1(4)+4:A1(5)+3);
x(1) = str2num(coeff_l1_str1);
x(2) = str2num(coeff_l1_str2);
x(3) = str2num(coeff_l1_str3);
x(4) = str2num(coeff_l1_str4);
x(5) = str2num(coeff_l1_str5);


cline2 = fgetl(f1);
A2 = findstr(cline2,'E');
coeff_l2_str1 = cline2(1:A2(1)+3);
coeff_l2_str2 = cline2(A2(1)+4:A2(2)+3);
coeff_l2_str3 = cline2(A2(2)+4:A2(3)+3);
coeff_l2_str4 = cline2(A2(3)+4:A2(4)+3);
coeff_l2_str5 = cline2(A2(4)+4:A2(5)+3);
x(6) = str2num(coeff_l2_str1);
x(7) = str2num(coeff_l2_str2);
x(8) = str2num(coeff_l2_str3);
x(9) = str2num(coeff_l2_str4);
x(10) = str2num(coeff_l2_str5);

cline3 = fgetl(f1);
A3 = findstr(cline3,'E');
coeff_l3_str1 = cline3(1:A3(1)+3);
coeff_l3_str2 = cline3(A3(1)+4:A3(2)+3);
coeff_l3_str3 = cline3(A3(2)+4:A3(3)+3);
coeff_l3_str4 = cline3(A3(3)+4:A3(4)+3);
x(11) = str2num(coeff_l3_str1);
x(12) = str2num(coeff_l3_str2);
x(13) = str2num(coeff_l3_str3);
x(14) = str2num(coeff_l3_str4);

% Calculating molecular weight of each element
mol_wg(i) = cal_mol_wg(species_name);

% Calculating specific heat
specific_heat = cal_cp(x,local_temp_2,local_temp_array,R);

% Calculating enthalpy
enthalpy = cal_H(x,local_temp_2,local_temp_array,R);

% Calculating entropy
entropy = cal_s(x,local_temp_2,local_temp_array,R);

% Creating Folders to save the plots.

FD = pwd;
mkdir(species_name)
cd(species_name)


% Plotting

% Plotting Temperature(K) vs Cp(J/Kg.K)

figure(1)
plot(local_temp_array,specific_heat,'lineWidth',3)
xlabel('Temperature(K)')
ylabel('Specific Heat (J/Kg.K)')
title(sprintf('Temperature vs Specific Heat of %s',species_name))

% Saving the plot as jpg file.

saveas(figure(1),'T vs Cp.jpg')

% Plotting Plotting Temperature(K) vs H(J)

figure(2)
plot(local_temp_array,enthalpy,'LineWidth',3)
xlabel('Temperature(K)')
ylabel('Enthalpy(J)')
title(sprintf('Temperature vs Enthalpy of %s',species_name))

% Saving the plot as jpg file.

saveas(figure(2),'T vs H.jpg')

 % Plotting Plotting Temperature(K) vs S(J/K)

figure(3)
plot(local_temp_array,entropy,'LineWidth',3)
xlabel('Temperature(K)')
ylabel('Entropy(J/K)')
title(sprintf('Temperature vs Entropy of %s',species_name))

% Saving the plot as jpg file.

saveas(figure(3),'T vs S.jpg')

cd(FD)

end
fclose(f1)


% Function for calculating Molecular weights

function mol_wg = cal_mol_wg(species_name)

mol_wg = 0;

key_species = ['H','C','N','O','S','A'];

at_mass = [1,12,14,16,32,40];

for i=1:length(species_name)
    for j= 1:length(key_species)
        if strcmp(species_name(i),key_species(j))
            mol_wg = mol_wg + at_mass(j);
            position = j;
        end
    end

    n = str2num(species_name(i));

    %checking the validity of n

    if n>1  
        mol_wg = mol_wg+(at_mass(position)*(n-1));
    end
end
% To print the molecular weight values of all species in the command window
fprintf('Molecular weight of species %s is %f \n',species_name,mol_wg)

end



% Function to Calculate specific heat

function [specific_heat] = cal_cp(x,local_temp_2,local_temp_array,R)

% Cp/R = x1 + x2*T + x3*T^2 + x4*T^3 + x5*T^4

for i = 1:length(local_temp_array)

    if local_temp_array(i) <= local_temp_2 
        % For low temperature range use last 7 coefficients, since they are
        % low temperature coefficients

        specific_heat(i) = R*(x(8) + (x(9)*local_temp_array(i)) + (x(10)*local_temp_array(i)^2) + (x(11)*local_temp_array(i)^3) + (x(12)*local_temp_array(i)^4));

    else specific_heat(i) = R*(x(1) + (x(2)*local_temp_array(i)) + (x(3)*local_temp_array(i)^2) + (x(4)*local_temp_array(i)^3) + (x(5)*local_temp_array(i)^4));
        
        %For high temperature range use first 7 coefficients, since they
        % are high temperature coefficients.
    end
end
end


% Function to calculate Enthalpy

function enthalpy = cal_H(x,local_temp_2,local_temp_array,R)

% H/RT = x1 + x2*T/2 + x3*T^2/3 + x4*T^3/4 + x5*T^4/5 + x6/T

for i = 1:length(local_temp_array)
    if local_temp_array(i) <= local_temp_2

% For low temperature range use last 7 coefficients, since they are low
% temperature coefficients.
 
  enthalpy(i) = R*local_temp_array(i)*(x(8) + x(9)*local_temp_array(i)/2 + x(10)*local_temp_array(i)^2/3 + x(11)*local_temp_array(i)^3/4 + x(12)*local_temp_array(i)^4/5 + x(13)/local_temp_array(i));

  %For high temperature range use first 7 coefficients, since they
        % are high temperature coefficients.

    else enthalpy(i) = R*local_temp_array(i)*(x(1) + x(2)*local_temp_array(i)/2 + x(3)*local_temp_array(i)^2/3 + x(4)*local_temp_array(i)^3/4 + x(5)*local_temp_array(i)^4/5 + x(6)/local_temp_array(i));

    end
end
end

% Function to calculate Entropy.

function entropy = cal_s(x,local_temp_2,local_temp_array,R)


for i = 1:length(local_temp_array)

    if local_temp_array(i) <= local_temp_2

        %  For low temperature range use last 7 coefficients, since they are low
% temperature coefficients.

entropy(i) = R*(x(8)*log(local_temp_array(i)) + x(9)*local_temp_array(i) + x(10)*local_temp_array(i)^2/2 + x(11)*local_temp_array(i)^3/3 + x(12)*local_temp_array(i)^4/4 + x(14));

% For high temperature range use first 7 coefficients, since they
        % are high temperature coefficients.
 
    else entropy(i) = R*(x(1)*log(local_temp_array(i)) + x(2)*local_temp_array(i) + x(3)*local_temp_array(i)^2/2 + x(4)*local_temp_array(i)^3/3 + x(5)*local_temp_array(i)^4/4 + x(7));

    end
end
end