File Parsing of NASA thermodynamic data

File parsing is an important aspect of data handling in MATLAB. In this project the data of the file \'therm.dat\'. The objective is to extract the coefficients of the gas to apply in the equations as follows to obtain the Cp, H and S value over the range of temperature given in the file.

The main code is as follows

clc
close all
clear all

% Opening the thermodynamic file
therm = fopen(\'thermo30.dat\',\'r\');

% Setting up initial conditional values
x = fgetl(therm);
n = fgets(therm);
n = strsplit(n,\' \');

% Setting the temperature list
for i = 1:3
    k(i) = n(i+1);
    t(i) = str2num(cell2mat(k(i)));
end
  min_temp = t(1);
  mid_temp = t(2);
  max_temp = t(3);
% Removing the comments from the file
for i = 1:3
   y = fgetl(therm);
end
clear y;

% Finding the coefficients and the temperatures of the gas
[ C , T , MW ] = get_coeff(therm); 
    

The code basically does the initial conditioning of the file such as removing comments, setting temperature, setting the line to start searching of gas and its coefficients.

To find coefficients the gas is set and empirical formula of the gas is set as the value of the string. The search for gas and extraction of coefficient is done by the function get_coeff defined as follows.

This function also uses the variable gas to

function [ C,T,MW ] = get_coeff( therm )
% Selecting the gas from the file
% Collecting coefficient from the data file
    for m = 1:53
        g = fgetl(therm);
        gs = strsplit(g,\' \');
        gas = cell2mat(gs(1));
        % Setting up the conditions for the selection of data 
        len = length(gs);
        % Setting up the temperature
        T(1) = str2double(cell2mat(gs(len-3)));
        T(2) = str2double(cell2mat(gs(len-1)));
        T(3) = str2double(cell2mat(gs(len-2)));
        % Extraction of coefficients
        Co = zeros(3,5);
        for i = 1:3
            a = fgetl(therm);
            k = strfind(a,\'E\');
            for j = 1:length(k)
                l = 15*(j-1)+1;
                u = 15*j;
                Co(i,j) = str2num(a(l:u));
            end
        end
        C = Co\';
        % Finding the molecular weight
        MW = set_gas_MW(gas);
        fprintf([\'Molecular/Atomic weight of \',gas]);
        fprintf(\' = %0.0d amu \\n\',MW);

        % Finding the properties of the gases
        [ Cp, H, S] = gas_property(C,T);

        % Plotting the Properties
        plot_gas_property(gas,Cp,H,S,T);
        
        cd ..\\..
    end
end

After gas is set Molecular/Atomic weight is calculated from the function set_gas_MW which is defined as follows.

function [ MW ] = set_gas_MW( gas )
% The Function finds the Molecular weight of the gas
% It breaks the string into single characters and generates molecular weight
% based on the database provided by IUPAC.
W = zeros(1,length(gas));
for  i = 1:length(gas)
    k = gas(i);
    % If k is not a character rather a number then the str2num(k) will not
    % be zero and hence the code will take the character before and assign
    % the atomic weight of that character and multiply k-1 as it had been
   if (str2num(k) ~= 0)
      l = gas(i-1);
      n = str2num(k);
      W(i) = Atom_Weight(l)*(n-1);
  else
      W(i) = Atom_Weight(k);
  end
end 
MW = sum(W);
end

MW gives the output as molecular weight and this output is printed in command window by the main program. This function uses Atom_weight function to define atomic weight.

The following function gives Molecular weight out of the string.

function [ W ] = Atom_Weight( s )
% To be used in set_gas function 
%   Sets the Molecular weight based on the string
while (s == \'H\')
    W = 1;
    break
end
while (s == \'C\')
    W = 12;
    break
end
while (s == \'N\')
    W = 14;
    break
end
while (s == \'O\')
    W = 16;
    break
end
while (s == \'A\')
    W = 40;
    break
end
while ( s == \'R\' || s == \'S\' || s == \'(\' || s == \')\')
    W = 0;
    break
end
end

Now the gas properties i.e. Cp, H and S are calculated by the following function gas_property to generate Cp, H ans S values. The C1 and C2 are the two sets of coefficients for two different temperature condtions i.e. low and high as described in the therm.dat file.

function [ Cp, H, S ] = gas_property( coeff_mat,temp_mat )
% Inputs coefficient matrix from the therm.dat file, outputs Cp, H and S 
%   Finds the property of the gas from the Cp/R, H/RT, S/R equations

% Simplifying inputs
C1 = coeff_mat(1:7);
C2 = coeff_mat(8:14);
T2 = linspace(temp_mat(1),temp_mat(2),1000);
T1 = linspace(temp_mat(2),temp_mat(3),1000);
% Universal gas constant R in J/Kg-K
R = 287;

% The properties of the gases are defined by the equation in following
% functions
Cp = get_Cp(C1,C2,T1,T2,R);
H = get_H(C1,C2,T1,T2,R);
S = get_S(C1,C2,T1,T2,R);
end

The calculation is done in the get_Cp, get_H and get_S according to equation as defined respectively in the following functions

function [ Cp ] = get_Cp( C1,C2,T1,T2,R )
% Calculates the Cp for high and low temperatures respectively
% Creating the Cp polynomials
Cpr1 = polyval(C1(5:-1:1),T2);
Cpr2 = polyval(C2(5:-1:1),T1);
Cp = [Cpr1*R;Cpr2*R];    
end



_______________________________________________________________________________________



function [ H ] = get_H( C1,C2,T1,T2,R )
% Find high and low temperature Enthalpies
% Enthalpy has different equation
for i = 1:6
    if (i < 6)
        Ch1(i) = C1(i)/i;
        Ch2(i) = C2(i)/i;
    else
        Ch1(i) = C1(i);
        Ch2(i) = C2(i);
    end
end
Hr1 = Ch1(1)*T2 + Ch1(2)*T2.^2 + Ch1(3)*T2.^3 + Ch1(4)*T2.^4 + Ch1(5)*T2.^5 + Ch1(6);
Hr2 = Ch2(1)*T1 + Ch2(2)*T1.^2 + Ch2(3)*T1.^3 + Ch2(4)*T1.^4 + Ch2(5)*T1.^5 + Ch2(6);
H = [Hr1*R ; Hr2*R];

end



_______________________________________________________________________________________



function [ S ] = get_S( C1,C2,T1,T2,R )
% Finds the high and low temperature Entropy  
% Now for Entropy we have following equations
for i = 1:7
    if  (1

Finally plotting is done as by taking the values from above set of functions and using it in plot_gas_property. This function also saves the plots in new folder named \'Result\'.

function [ W ] = plot_gas_property( gas,Cp,H,S,T )
% Takes in the values and plots the graph of Cp, H and S vs T
% Simple function to plot the gas properties in Matlab from therm.dat file  

% Creating Temperature array
t1 = linspace(T(1),T(2),1000);
t2 = linspace(T(2),T(3),1000);

% For dumping images in new directory we have to create new directory
mkdir(\'Result\',[\'\',gas]);
cd(\'Result\')
cd([\'\',gas])

% Generating the required plots
figure (1)
hold on
plot(t1,Cp(1,:),\'color\',\'b\');
plot(t2,Cp(2,:),\'color\',\'r\');
title([\'The plot of C_{p} vs Temperature of \',gas]);
xlabel(\'Temperature\');
ylabel(\'Specific heat at const. Pressure\');
legend(\'Low Temperature\',\'High Temperature\');
savefig([\'Cp Plot of \',gas]);
figure (2)
hold on
plot(t1,H(1,:),\'color\',\'b\');
plot(t2,H(2,:),\'color\',\'r\');
title([\'H vs Temp of \',gas]);
xlabel(\'Temperature\');
ylabel(\'Enthalpy\');
legend(\'Low Temperature\',\'High Temperature\');
savefig([\'H Plot of \',gas]);
figure (3)
hold on
plot(t1,S(1,:),\'color\',\'b\');
plot(t2,S(2,:),\'color\',\'r\');
title([\'S vs Temp of \',gas]);
xlabel(\'Temperature\');
ylabel(\'Entropy\');
legend(\'Low Temperature\',\'High Temperature\');
savefig([\'S Plot of \',gas]);

close all
end

The close all in the end of this program is important as the plot of the previous reults will superimpose over the recentplot aking an undesirable plot.

Finally the project is complete and all the required results were acheived and the plots were saved in the folders named after the gas in the main folder of \'Result\'.