Project 1 - Parsing NASA thermodynamic data

Aim: 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.

Objective:

Writing a function that extracts the 14 co-efficients and calculates the enthalpy, entropy and specific heats for all the species in the data file. By reading this file NASA Thermodynamic data file. (Tip: use fopen and fgetl). The formulae are shown below.

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

The FIRST 7 coefficients are HIGH-temperature coefficients and the SECOND 7 coefficients are LOW-temperature coefficients.

By using LOCAL TEMPERATURE RANGES which is unique for each individual species, not the GLOBAL TEMPERATURE RANGE.

2. Calculate the molecular weight of each species and display it in the command window.
3. Plotting the Cp, Enthalpy and Entropy for the local temperature range (low temp high temp) specific for each species.
4. Saving the plots as images with appropriate names and dump them in separate folders for each species with suitable name.

Theroy:

Specific heat:

The specific heat capacity of a substance is the heat capacity of a sample of the substance divided by the mass of the sample. Informally, it is the amount of energy that must be added, in the form of heat, to one unit of mass of the substance in order to cause an increase of one unit in its temperature. The SI unit of specific heat is joule per kelvin and kilogram, J/.

The specific heat constants for constant pressure and constant volume processes are related to the gas constant for a given gas. This rather remarkable result has been derived from thermodynamic relations, which are based on observations of physical systems and processes. Using the kinetic theory of gases, this same result can be derived from considerations of the conservation of energy at a molecular level.

We can define an additional variable called the specific heat ratio, which is given the Greek symbol "gamma", which is equal to cp divided by cv:

gamma = cp / cv

Enthalpy:

Enthalpy, a property of a thermodynamic system, is the sum of the system's internal energy and the product of its pressure and volume. It is a state function used in many measurements in chemical, biological, and physical systems at a constant pressure, which is conveniently provided by the large ambient atmosphere.

If we are dealing with a gas, it is most convenient to use forms of the thermodynamics equations based on the enthalpy of the gas. From the definition of enthalpy:

h = e + p * v

Entropy:

Entropy is a scientific concept as well as a measurable physical property that is most commonly associated with a state of disorder, randomness, or uncertainty. The term and the concept are used in diverse fields, from classical thermodynamics, where it was first recognized, to the microscopic description of nature in statistical physics, and to the principles of information theory.

What is File parsing ?

File parsing means reading a file from the outer folder inside the matlab or extracting the value of a folder provided extracting the value inside the folder. The name Parse refers in matlab reading a function, script or internal code inside the memory. Some of the common commands used are (fopen, fgetl, fgets etc).

We are using this technique to read the Thermodynamic data by using file parsing to get the value of Cp, H, and S.

Explanation:

• Getting the heading of the file.
• The given program will get the header of the file which start with clear all close all and clc
• The command fopen will used for reading the thermodynamic data.
• The fgetl command returns the next line of the specified file, removing the newline characters.
• If the file is nonempty, then fgetlreturns tline as a character vector.
• If the file is empty and contains only the end-of-file marker, then fgetlreturns tline as a numeric value -1.

Code:

%% NASA Thermodynamic Data and calculating the molecular weight
clear all 
close all
clc

diary Results_NASA_Thermo_Data.txt;
% Getting the header info
f1 = fopen('THERMO.data', 'r');

header = fgetl(f1);

Getting the Global Temperature:

• The global temperature info is stored in string so we use the command str2double so that it returns the value in {} as a string.
• Which gives us the cell array. From which we can get the global temperature value.
• For skipping the next three comments we use fgetl command.

Code:

% Getting the temperature data

global_temp = fgetl(f1);

x = strsplit(global_temp);

global_temp_low = str2double(x{2});
global_temp_mid = str2double(x{3});
global_temp_high = str2double(x{4});

fgetl(f1);
fgetl(f1);
fgetl(f1);

Getting the Temperature info of individual species:

• Species name is store in length of 6 for that I used line1(1:6) matlab will get me the value and so on for getting the low_temp value which is on 49 to 55^th location I used the command line1(49:55) matlab will get me the value. As well for higher temp and mid temp.
• Here in this section I have use the terms local_temp which show I’m targeting the temperature of the species.
• And then converting to the string by using again str2double command. With giving temperature name low_tmp, high_temp, mid_temp.

Code:

% getting the line
for i = 1:53
    line1 = fgetl(f1);
    species_name = line1(1:6);
    
    % getting local temp
    local_temp1 = line1(49:55);
    low_temp = str2double(local_temp1);
    
    local_temp2 = line1(58:65);
    mid_temp = str2double(local_temp2);
    
    local_temp3 = line1(68:75);
    high_temp = str2double(local_temp3);

Extracting the co-efficient values:

• In order to extract co-efficient value I got the lines by fgetl command and marked the lines as line2, line3, line4 and I use these line to extract the co-efficient.
• By using the line2, line3, line4 I will try to separate the co-efficient taking a constant E and I will use the command strfind to find E, which lies at same position at every length.
• So I target E as the extractor of my co-efficient plus the extra numbers. And same for all the 14 co-efficient.
• And to converted string I have used the command str2double.

Code:

%% Extracting the main coefficients from the data
    
    % coefficient
    line2 = fgetl(f1);
    line3 = fgetl(f1);
    line4 = fgetl(f1);
    
    a = strfind(line2, 'E');
    b = strfind(line3, 'E');
    c = strfind(line4, 'E');
    
    % For low temp
    a1 = str2double(line2(1:a(1)+3));
    a2 = str2double(line2(a(1)+4:a(2)+3));
    a3 = str2double(line2(a(2)+4:a(3)+3));
    a4 = str2double(line2(a(3)+4:a(4)+3));
    a5 = str2double(line2(a(4)+4:a(5)+3));
    a6 = str2double(line3(1:b(2)+3));
    a7 = str2double(line3(b(1)+4:b(2)+3));
    a8 = str2double(line3(b(2)+4:b(3)+3));
    a9 = str2double(line3(b(3)+4:b(4)+3));
    a10 = str2double(line3(b(4)+4:b(5)+3));
    a11 = str2double(line4(1:c(1)+3));
    a12 = str2double(line4(c(1)+4:c(2)+3));
    a13 = str2double(line4(c(2)+4:c(3)+3));
    a14 = str2double(line4(c(3)+4:c(4)+3));

Specific Heat function:

• Here I have used the same formula as mention above. To calculate the cp, h and s value.
• For loop is for repeating these condition for all three at once cp, h and s.
• If loop is used here because there are condition applied for finding value all 7 co-efficient of higher temperature. And value of 7 co-efficient of lower temperature.

Code:

function [cp, h, s] = coeff(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T)

% Given conditions
R = 8.314;
for j = 1:length(T)
    if T(j) < 300
        
        cp = (a1+a2.*T+a3.*T.^2+a4.*T.^3+a5.*T.^4).*R;
        h = (a1+(a2.*T)./2+(a3.*T.^2)./3+(a4.*T.^3)./4+(a5.*T.^4)./5+(a6./T)).*R.*T;
        s = (a1.*log(T)+a2.*T+(a3.*T.^2)./2+(a4.*T.^3)./3+(a5.*T.^4)./4+a7).*R;
    else
        cp = (a8+a9.*T+a10.*T.^2+a11.*T.^3+a12.*T.^4).*R;
        h = (a8+(a9.*T)./2+(a10.*T.^2)./3+(a11.*T.^3)./4+(a12.*T.^4)./5+(a13./T)).*R.*T;
        s = (a8.*log(T)+a9.*T+(a10.*T.^2)./2+(a11.*T.^3)./3+(a12.*T.^4)./4+a14).*R;

    end
end
end

Molecular Weight function:

• Hence we can see the thermodynamic data contains 5 species that are oxygen, carbon, nitrogen, Hydrogen, Argon.
• These are named as elements and are stored in array.
• Where as their atomic masses are stored in atomic mass array.
• Here I have used For loop so that the process is to be repeated for each and every species.
• Strcmp command is used to compare the elements and corresponding atomic mass.
• And the molecular weight is stored in mol_wt.
• And Hence the molecular weight is calculated.
• Now I have to pass this function to the main code.

Code:

function [molecular_weight_cal] = mole_wt(species_name)

elements = {'H','N','O','C','A'};
atomic_mass = [1 14 16 12 40];

mol_wt = 0;
for i = 1:length(species_name)
    for j = 1:length(elements)
        if strcmp(species_name(i),elements(j))
            mol_wt = mol_wt + atomic_mass(j);
            m = j;
        end
    end

   s = str2double(species_name(i));
    if s > 1
        mol_wt = mol_wt+(atomic_mass(m)*(s-1));
    end
end
molecular_weight_cal = mol_wt
fprintf('molecular_weight_cal of %s:- ',species_name)
fprintf('%d',molecular_weight_cal)
fprintf(' ')
fprintf('amu')
fprintf('\n')
end

Creating a folder:

• Inorder to create a folder I have used the command mkdir which allows matlab to create a new folder to store the species data.
• Strcat command is used to concatenate a string to identify the name of the folder.
• Pwd command is used to identify the nature of the folder. For example .txt .file etc.

Code:

 % Creating a folder 'SPECIES' to save all the species_name
    species_name = strcat('Property ...',species_name);
    mkdir(species_name);
    file = pwd; 

Plotting the graph:

• Plot the graph between cp and temperature, enthalpy vs temperature, entropy vs temperature.
• Saving the graph by using fullfill command which store the entire data generated into the folder.
• Gcf command is used to handle the current figure.
• Cd command is used to define path to the folder where the data has to be stored.

Code:

 %% Plotting the graph
    %plotting between Temperature and specific heat
    
    figure(1)
    plot(T,cp,'linewidth',3,'color','r')
    title('specific heat vs temperature',species_name);
    xlabel('Temperature[k]')
    ylabel('Specific Heat [kj/kg.k]')
    f = fullfile(species_name);
    cd(f);
    saveas(gcf,[pwd '/specific_heat vs temperature.jpg']);
    cd(file);
    
    %plotting between Temperature and Enthalphy
    
    figure(2)
    plot(T,h,'linewidth',3,'color','b')
    title('Enthalpy vs Temperature',species_name);
    xlabel('Temperature[K]')
    ylabel('Enthalpy[KJ]')
    f = fullfile(species_name);
    cd(f);
    saveas(gcf,[pwd '/enthalphy vs temperature.jpg']);
    cd(file);
    
    %plotting between Temperature and Entropy
    
    figure(3)
    plot(T,s,'LineWidth',3,'color','c');
    title('Entropy vs Temperature',species_name)
    xlabel('Temperature[K]')
    ylabel('Entropy[KJ/K]')
    f = fullfile(species_name);
    cd(f);
    saveas(gcf,[pwd '/entropy vs temperature.jpg']);
    cd(file);

Main Code:

%% NASA Thermodynamic Data and calculating the molecular weight
clear all 
close all
clc

diary Results_NASA_Thermo_Data.txt;
% Getting the header info
f1 = fopen('THERMO.data', 'r');

header = fgetl(f1);

% Getting the temperature data

global_temp = fgetl(f1);

x = strsplit(global_temp);

global_temp_low = str2double(x{2});
global_temp_mid = str2double(x{3});
global_temp_high = str2double(x{4});

fgetl(f1);
fgetl(f1);
fgetl(f1);

% getting the line
for i = 1:53
    line1 = fgetl(f1);
    species_name = line1(1:6);
    
    % getting local temp
    local_temp1 = line1(49:55);
    low_temp = str2double(local_temp1);
    
    local_temp2 = line1(58:65);
    mid_temp = str2double(local_temp2);
    
    local_temp3 = line1(68:75);
    high_temp = str2double(local_temp3);
    
    %% Extracting the main coefficients from the data
    
    % coefficient
    line2 = fgetl(f1);
    line3 = fgetl(f1);
    line4 = fgetl(f1);
    
    a = strfind(line2, 'E');
    b = strfind(line3, 'E');
    c = strfind(line4, 'E');
    
    % For low temp
    a1 = str2double(line2(1:a(1)+3));
    a2 = str2double(line2(a(1)+4:a(2)+3));
    a3 = str2double(line2(a(2)+4:a(3)+3));
    a4 = str2double(line2(a(3)+4:a(4)+3));
    a5 = str2double(line2(a(4)+4:a(5)+3));
    a6 = str2double(line3(1:b(2)+3));
    a7 = str2double(line3(b(1)+4:b(2)+3));
    a8 = str2double(line3(b(2)+4:b(3)+3));
    a9 = str2double(line3(b(3)+4:b(4)+3));
    a10 = str2double(line3(b(4)+4:b(5)+3));
    a11 = str2double(line4(1:c(1)+3));
    a12 = str2double(line4(c(1)+4:c(2)+3));
    a13 = str2double(line4(c(2)+4:c(3)+3));
    a14 = str2double(line4(c(3)+4:c(4)+3));
    
    % Creating a folder 'SPECIES' to save all the species_name
    species_name = strcat('Property ...',species_name);
   
    file = pwd;   
    
    % Applying universal gas constant
    R = 8.314;
    T = linspace(global_temp_low,global_temp_high,500);
    
    % Calculation speicfic heat
    [cp, h, s] = coeff(a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12,a13,a14,T);
    
    %% Plotting the graph
    %plotting between Temperature and specific heat
    
    figure(1)
    plot(T,cp,'linewidth',3,'color','r')
    title('specific heat vs temperature',species_name);
    xlabel('Temperature[k]')
    ylabel('Specific Heat [kj/kg.k]')
    f = fullfile(species_name);
    cd(f);
    saveas(gcf,[pwd '/specific_heat vs temperature.jpg']);
    cd(file);
    
    %plotting between Temperature and Enthalphy
    
    figure(2)
    plot(T,h,'linewidth',3,'color','b')
    title('Enthalpy vs Temperature',species_name);
    xlabel('Temperature[K]')
    ylabel('Enthalpy[KJ]')
    f = fullfile(species_name);
    cd(f);
    saveas(gcf,[pwd '/enthalphy vs temperature.jpg']);
    cd(file);
    
    %plotting between Temperature and Entropy
    
    figure(3)
    plot(T,s,'LineWidth',3,'color','c');
    title('Entropy vs Temperature',species_name)
    xlabel('Temperature[K]')
    ylabel('Entropy[KJ/K]')
    f = fullfile(species_name);
    cd(f);
    saveas(gcf,[pwd '/entropy vs temperature.jpg']);
    cd(file);
    
    
    
    
    molecular_weight_cal = mole_wt(species_name);
    k = sprintf('%s mol-wt = %f',species_name,molecular_weight_cal);
    disp(k);

end 
fclose(f1)

Plots:

Co2 :

N2:

O2:

Errors:

• There are basic syntactic error and some new command error which have been resolved immediately.

Saving the Command window Data:

• For saving the Command Window data I have used the command ‘dairy’ which will store in txt file in the same directory. For example (diary Results.txt).

Conclusion:

• The thermodynamic data has been parsed, and thermodynamic data has been calculated and plotted for each species.
• It has been clearly seen that the data of cp, h and s vary at different temperature range from low-temperature to high-temperature at each compound.
• Hence we have calculated the Molecular weight of different species and stored the data using file parsing technique.

References:

• https://en.wikipedia.org/wiki/Specific_heat_capacity
• https://en.wikipedia.org/wiki/Enthalpy
• https://en.wikipedia.org/wiki/Entropy