Project 1 - Parsing NASA thermodynamic data

Aim:

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

Objective:

• To write a function that extracts the 14 coefficients and calculates the enthalpy, entropy, and specific heats for all the species in the NASA thermodynamic data file "Thermo.dat".

The NASA polynomials have the form as shown below:

Where R is the universal gas constant, T is the temperature

• To calculate the molecular weight of each species and display it in the command window.
• To plot the Cp, Enthalpy, and Entropy for the local temperature range (low temperature: high temperature) specific for each species.
• To save the plots as images with appropriate names and dump them in separate folders for each species with suitable name automatically.

File parsing:

File parsing refers to how information can be read from a particular text file or any type of file and how that information can be used in a computer program.

Step 1 - Define the function

Step 2 - Create an Input Parser object

Step 3 - Add inputs to the scheme.

Step 4 - Set properties to adjust parsing.

Step 5 - Parse the inputs.

Step 6 - Use the inputs in function.

Step 7 - Calling function.

NASA thermodynamic file:

NASA thermodynamics file which is THERMO.dat file in which there are temperature range and 53 species are given in these data and each species contains 14 coefficients and temperature ranges and we have to parse these 14 coefficients to find Specific heat, Entropy and Enthalpy using formula as:-

Where,

R = Universal Gas constant

T = Temperature

Cp = Specific heat (KJ)

S = Entropy (KJ/mol-K)

H = Enthalpy (KJ/mol)

a(i) = Coefficients where, i = 1 to 14

Specific heat:

The specific heat is the amount of heat per unit mass required to raise the temperature by one degree Celsius. The relationship between heat and temperature change is usually expressed in the form shown below where c is the specific heat.

Entropy:

Entropy, the measure of a system's thermal energy per unit temperature that is unavailable for doing useful work. Because work is obtained from ordered molecular motion, the amount of entropy is also a measure of the molecular disorder, or randomness, of a system.

Enthalpy:

Enthalpy, a property of a thermodynamic system, is equal to the system's internal energy plus the product of its pressure and volume. In a system contained so as to prevent mass transfer, for processes at constant pressure, the heat absorbed or released equals the change in enthalpy.

MATLAB code:

Main code:

clear all
close all
clc

% input
R = 8.314; % gas constant

% reading data file
% open file and reading THERMO.data file
f1 = fopen('THERMO.dat','r')
get = fgetl(f1)

% now reading temperature from this file
temp = fgetl(f1)

% spliting 
s = strsplit(temp);

% converting strings cells into numbers
templow = str2num(s{2})
tempmid = str2num(s{3})
temphigh = str2num(s{4})

% Temperature range
T = linspace(templow,temphigh,1000);

% skipping comment lines
for i = 1:3
    get = fgetl(f1);
end

% reading the species
for j = 1:53
    tline = fgetl(f1);

    A = strsplit(tline,'');
    species = (A{1});

    % converting string number
    %lowtemp = str2num(A{end-3});
   % midtemp = str2num(A{end-2});
    %hightemp = str2num(A{end-1});

    tline1 = fgetl(f1);
    a = strfind(tline1,'E');

    % finding first 7 higher and second 7 lower co-efficients
    % hegher co-efficients
    a1 = tline1(1:a(1)+3);
    a1 = str2num(a1);
    a2 = tline1(a(1)+4:a(2)+3);
    a2 = str2num(a2);
    a3 = tline1(a(2)+4:a(3)+3);
    a3 = str2num(a3);
    a4 = tline1(a(3)+4:a(4)+3);
    a4 = str2num(a4);
    a5 = tline1(a(4)+4:a(5)+3);
    a5 = str2num(a5);

    tline2 = fgetl(f1);
    b = strfind(tline2,'E');
    a6 = tline2(1:b(1)+3);
    a6 = str2num(a6);
    a7 = tline2(b(1)+4:b(2)+3);
    a7 = str2num(a7);

    % finding lower co-efficients
    a8 = tline2(b(2)+4:b(3)+3);
    a8 = str2num(a8); 
    a9 = tline2(b(3)+4:b(4)+3);
    a9 = str2num(a9);    
    a10 = tline2(b(4)+4:b(5)+3);
    a10 = str2num(a10);

    tline3 = fgetl(f1);
    c = strfind(tline3,'E');
    a11 = tline3(1:c(1)+3);
    all = str2num(a11);
    a12 = tline3(c(1)+4:c(2)+3);
    a12 = str2num(a12);
    a13 = tline3(c(2)+4:c(3)+3);
    a13 = str2num(a13);    
    a14 = tline3(c(3)+4:c(4)+3);
    a14 = str2num(a14);

    % calculate specific heat
    Cp = specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,tempmid);

    % calculate Entropy
    S = Entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,tempmid);

    % calculate Enthalpy
    H = Enthalpy (a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,tempmid);    

    % calculate molecular weight of the species
    W = Molecular_weight(species);

    % writng molecular weight file
    MW_file = fopen('Molecular Weigth.txt','w');
    fprintf(MW_file,'Molcular weight of %s is %d', species,W);
    fclose(MW_file); % closing the file

    % plotting Cp, H, and S with Temperature
    % saving in a folder
    cd('E:\matlab\project 1')
    mkdir(species);
    cd(species);

    % plot 1 Cp vs T
    figure(1)
    plot(T,Cp,'linewidth',2,'color','b');
    xlabel('Temperature(K)');
    ylabel('Specific Heat(KJ)');
    title(sprintf('specific heat vs Temperature range '));
    grid on
    saveas(figure(1),('specific heat.jpg'));

    % plot 2 S vs T
    figure(2)
    plot(T,S,'linewidth',2,'color','r');
    xlabel('Temperature(K)');
    ylabel('Entropy(KJ/mol-K)');
    title(sprintf(' Entropyvs Temperature range '));
    grid on
    saveas(figure(2),('Entropy.jpg'));

    % plot 3 H vs T
    figure(3)
    plot(T,H,'linewidth',2,'color','g');
    xlabel('Temperature(K)');
    ylabel('Enthalpy(KJ/mol-K)');
    title(sprintf(' Enthalpy vs Temperature range '));
    grid on
    saveas(figure(2),('Enthalpy.jpg'));

    cd('E:\matlab\project 1')

    % Molecular weight

    W = Molecular_weight(species)

end

Steps:

• Firstly, we will Input the parameter R and we will open a file dat file and reading this file using fopen command.
• Then, we will read first THERMO using fgetl command and then we read temperature file and we will split the temperatures using strsplit command to convert it into cell array and then we will convert this string to number using str2num command to find low, high and mid temperatures.
• Then we will assign temperature range using linspace Then we run loop to skip command lines the data file using fgetlcommand inside the loop.
• Now, we will run the loop to find 53 element in the data file by assigning fgetl command and then spliting into the cell array to find the species line and then find the species name and then find the number of the temperatures in the species line.
• Then we read the new line. In order to read the coefficient in the data file we will first find the 'E' string using strfind command to find before and after numbers of the coefficients then we will find remaining 14 coefficients and convert them to numbers from strings.
• Out of which we will find first 7 Coefficients which are for Higher temperature and then Second 7 Coefficients which are for Lower temperature.
• Now we will specify specific heat (Cp), Entropy (S), Enthalpy (H) functions to find these parameters using these 14 coefficients for Lower and higher Temperature. then we will specify molecular weight ᵤ function to find the Molecular weight of the different species.
• Now we will write the Molecular weight.txt file and fprintf all the molecular weights of the species. and then close the txt file using fclose
• Now we change the directory using cd command and then creating directory folder where we have to save our files using mkdir command and then we will change the directory to this created folder.
• Now, we will plot specific heat, entropy and enthalpy with the temperature for different species. We we will use scanf command to in the title command to create different species title.
• Then we will use saveas command to save the plots of different species in the directory folder. Now we will change the directory to which our main code file runs in the MATLAB.

Function to find specific heat:

function [Cp] = specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,T,R,tempmid);

    % calculating specific heat
    % calculating specific heat for higher and lower temperature repectively
if (T <= tempmid)
    Cp = R * (a8 + (a9*T) + (a10*(T).^2) + (a11*(T).^3) + (a12 * (T).^4));

else
    Cp = R*(a1 + (a2*T) + (a3*(T).^2) + (a4*(T).^3) + (a5*(T).^4));
end
end

Steps:

• Now we will find the Specific heat by creating Specific heat function and assigning its relted formula coefficient from 14 coefficients, Gas constant, Temperature and mid temperature.
• Now we will create an condition using if loop and make the temperature less than and equal to mid temperature to specify lower temperature coefficient in the formula.
• Then we will use else condition to specify coefficient of Higher temperature coefficient in the formula.
• Then end the loop and function to calculate Specific Heat.

Function to find Entropy:

function S = Entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,T,R,tempmid)

    % calculating entropy
    % calculating entropy for higher and lower temperature repectively.

    if (T <= tempmid)
        S = R*(a8.*log(T) + (a9*T) + ((a10*(T).^2)/2) + ((a11*(T).^3)/3) + ((a12*(T).^4)/4) + a14);
    else
        S = R*(a1.*log(T) + (a2*T) + ((a3*(T).^2)/2) + ((a4*(T).^3)/3) + ((a5*(T).^4)/4) + a7);
    
    end
end

Steps:

• Now we will find the Entropy by creating Entropy function and assigning its relted formula coefficient from 14 coefficients, Gas constant, Temperature and mid temperature.
• Now we will create an condition using if loop and make the temperature less than and equal to mid temperature to specify lower temperature coefficient in the formula.
• Then we will use else condition to specify coefficient of Higher temperature coefficient in the formula.
• Then end the loop and function to calculate Entropy.

Enthalpy:

function H = Enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,T,R,tempmid)

    % calculating enthalpy
    % calculating enthalpy for higher and lower temperature repectively.

    if (T <= tempmid)
        H = R*T.*(a8 + ((a9*T)/2) + ((a10*(T).^2)/3) + ((a11*(T).^3)/4) + ((a12*(T).^4)/5) + ((a13./T)));
    else
        H = R*T.*(a1 + ((a2*T)/2) + ((a3*(T).^2)/3) + ((a4*(T).^3)/4) + ((a5*(T).^4)/5) + ((a6./T)));
    end
end

Steps:

• Now we will find the Enthalpy by creating Enthalpy function and assigning its relted formula coefficient from 14 coefficients, Gas constant, Temperature and mid temperature.
• Now we will create an condition using if loop and make the temperature less than and equal to mid temperature to specify lower temperature coefficient in the formula.
• Then we will use else condition to specify coefficient of Higher temperature coefficient in the formula.
• Then end the loop and function to calculate Enthalpy.

Function to find molecular weight of the different species:

function W = Molecular_weight(species)

    % Defining Atomic name of the species
    atomic_name = ['H','C','O','N','Ar'];
    
    % Defining Atomic masses of main species
    atomic_weight = [1.008 12.011 15.999 14.007 39.948];
    
    W=0;  % Molecular weight starting range
    
    for i=1:length(species)   % Loop runs to the length of the species
        for j=1:length(atomic_name)  % Loop runs to the length of the atomic names
            if strcmp(species(i), atomic_name(j))  %function to compare 2 strings
                W = W + atomic_weight(j);
                position = j;
            end
        end
        n = str2num(species);   % Now we will find if there is more element in the species then incrementing atomic weight of the species
        
        if n>1
            W = W + atomic_weight(position)*(n-1);
        end
    end
    
    fprintf('Molecular Weight of %s :- ',species);
    fprintf('%f',W);
    disp(' ')
end

Plots:

Plots of specie O2:

specific heat:

Entropy:

Enthalpy: