Challenge

Parsing NASA thermodynamic data

DDamodhar Jangam4 Dec 2020

AIM:-

To parse the NASA thermodynamic data file and then to evaluate thermodynamic properties such as Cp, H, and S of various gas species through MATLAB

OBJECTIVE:-

Writing a function that extracts the 14 co-efficient and calculation of the enthalpy, entropy, and specific heats for all the species in the data file.
Calculation the molecular weight of each species and display it in the command window.

Plot the Cp, Enthalpy, and Entropy for the local temperature range (low temperature: high temperature) specific for each species.
Save the plots as images with appropriate names and dump them in separate folders for each species with a suitable name.
All of these should be generated automatically.
The main program that you will be submitted will just take the "THERMO.dat" file as an input and generation of Cp, H, and S plots. Screenshot of the window where all the folders are to be saved.
Plots for O2, N2, and CO2 species.

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 the formula as below

where,

R = Universal Gas constant

T = Temperature

C_p = Specific heat at constant pressure (KJ/mol)

S = Entropy (KJ/mol-K)

H = Enthalpy (KJ/mol)

Specific heat at constant pressure (C_p)

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

Entropy (s)

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 (h)

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 to prevent the mass transfer, for processes at constant pressure, the heat absorbed or released equals the change in enthalpy.

CODE

Main code

clc
clearvars
close all
% converting given data into charts and finding their molecular weight


% input
R = 8.314; % Universal Gas constant f
f1 = fopen('THERMO.dat','r');
No_of_component = 53;
Coefficients = 14;

% reading section start now 
fgetl(f1);
A = fscanf(f1,'%f');
global_low_temp  = A(1);
global_mid_temp  = A(2);
global_high_temp = A(3);

% skipping comments
for i = 1 : 3
    fgetl(f1);
end
% creating new folder
mkdir('component'); 
% making for loop for reading components cofficient and local temperature
for i = 1 : No_of_component
    tline_0 = fgetl(f1);
    tline_1 =  strsplit(tline_0,' ');
    component(i,1) = tline_1(1);
    low_temp(i,1) =  str2double(tline_1{end-3});
    
    high_temp(i,1) =  str2double(tline_1{end-2});
    
    mid_temp(i,1) =  str2double(tline_1{end-1});
    
    tline_2 = fgetl(f1);
    
    E = strfind(tline_2,'E');
    
    
    for k = 1:4
        j = k + 1;
        if k < 2
            a(i,k) = str2double(tline_2(1:E(k)+3));
        end
        a(i,j) = str2double(tline_2(E(k)+4:E(j)+3));
    end
    
    
    tline_3 = fgetl(f1);
    
    E = strfind(tline_3,'E');

    for k = 1:4
        j = k + 1;
        if k < 2
            a(i,k+5) = str2double(tline_3(1:E(k)+3));
        end
        a(i,j+5) = str2double(tline_3(E(k)+4:E(j)+3));
    end
    
    tline_4 = fgetl(f1);
    
    E = strfind(tline_4,'E');

    for k = 1:3
        j = k + 1;
        if k < 2
            a(i,k+10) = str2double(tline_4(1:E(k)+3));
        end
        a(i,j+10) = str2double(tline_4(E(k)+4:E(j)+3));
    end
    
end
 %% creating new file to write coefficients and their respective molecular weights
     f1 = fopen('All NASA Data seperation','w');
     fprintf(f1,'  S.no)Components   | Molecular Weight  |   ');
     for j = 1 : Coefficients 
        fprintf(f1,'|      (a%.2d)      |',j);
     end


    fprintf(f1,'n');
    for i = 1 : No_of_component
        
        fprintf(f1,' %5.f)%7st    |t%8.5f    t|t',i,component{i,1},Molecular_Weight(component{i,1}));
        
        for j = 1 : Coefficients
            fprintf(f1,'| %15d |',a(i,j));
        end
            fprintf(f1,'n');
            
    end
    fclose(f1);
  %% plot and molecular weight calculation
for i = 1 : No_of_component
    % calling function
    W = Molecular_Weight(component{i,1});
    T = linspace(mid_temp(i,1),high_temp(i,1),1000);
    C_p = Constant_pressure(a(i,1),a(i,2),a(i,3),a(i,4),a(i,5),a(i,8),a(i,9),a(i,10),a(i,11),a(i,12),R,T,mid_temp(i,1));
    H = Enthalpy(a(i,1),a(i,2),a(i,3),a(i,4),a(i,5),a(i,6),a(i,8),a(i,9),a(i,10),a(i,11),a(i,12),a(i,13),R,T,mid_temp(i,1));
    S = Entropy(a(i,1),a(i,2),a(i,3),a(i,4),a(i,5),a(i,7),a(i,8),a(i,9),a(i,10),a(i,11),a(i,12),a(i,14),R,T,mid_temp(i,1));
    cd('D:mathlabweek4nasa data segerationcomponent')       
    % creating folder with component name
    mkdir(component{i,1});                                      
    cd(component{i,1});
    
    temp = fopen('coefficients.txt','w');
    % for printing making header
    fprintf(f1,'Molecular Weight ofn %.4s : ',component{i,1});   
    fprintf('Molecular Weight of n %6s : ',component{i,1});
    fprintf('%f',W);
    fprintf(f1,'%f',W);
    disp(' ')
    fprintf(f1,' r  Componet namet%6stn',component{i,1});
    % for printing 14 Coefficients
    for j = 1 : Coefficients
            fprintf(f1,'ta%.2d t - t%15dn',j,a(i,j));
    end
    fclose(temp);
    % save as Specific Heat.jpg
    figure(1)
    plot(T,C_p,'linewidth',2,'color','b');
    xlabel('Temperature(K)');
    ylabel('Specific Heat at constant pressure(KJ/mol-K)');
    title(sprintf('Specific heat at constant pressure v/s Temperature range of "%s"',component{i,1}));
    grid on
    saveas(figure(1),'Specific Heat.jpg');
    
    % save as Entropy.jpg
    figure(2)
    plot(T,S,'linewidth',2,'color','r');
    xlabel('Temperature(K)');
    ylabel('Entropy(KJ/mol-K)');
    title(sprintf('Entropy v/s Temperature range of "%s"',component{i,1}));
    grid on
    saveas(figure(2),'Entropy.jpg');
    
    % save as Enthalpy.jpg
    figure(3)
    plot(T,H,'linewidth',2,'color','g');
    xlabel('Temperature(K)');
    ylabel('Enthalpy(KJ/mol)');
    title(sprintf('Enthalpy v/s Temperature range of "%s"',component{i,1}));
    grid on
    saveas(figure(3),'Enthalpy.jpg');
    
    % come back to starting folder
    cd('D:mathlabweek4nasa data segeration'); 
    
end

STEPS:-

The program starts with Input the parameter R which is a universal constant and the open a file THERMO.dat file and reading this file using fopen command.
Then, we will read the first THERMO using fgetl command and then we read the 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 a number using str2num command to find low, high and mid temperatures.
Now, run the loop to find 53 components in the data file by assigning fgetl command and then splitting into the cell array to find the species line and then find the species name and then find the number of the local temperatures in the species line.
Then assign a local temperature range using linspace command. Then we run a loop to skip command lines the data file using fgetl command inside the loop.
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 the remaining 14 coefficients and convert them to numbers from strings.
Out of which we will find the first 7 Coefficients which are for Higher temperature and then the Second 7 Coefficients which are for Lower temperature.
Now we will specify Specific heat, Entropy, Enthalpy functions to find these parameters using these 14 coefficients for Lower and Higher Temperatures. then we will specify the molecular weight function to find the Molecular weight of the different species.
Now we will write the Molecular weight and fprintf all the molecular weights of the species. and then close the file using fclose command.
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 a different component. legend will use scanf command in the title command to create different species title.
Then we will use save as a 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 MATLAB.

Function to find Specific Heat at constant pressure

a function call with Constant_pressure with the following parameters

function Cp = Constant_pressure(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,R,t,mt)
if t > mt
    Cp = R*(a1+(a2*t)+(a3*t.^2)+(a4*t.^3)+(a5*t.^4));
else
    Cp = R*(a8+(a9*t)+(a10*t.^2)+(a11*t.^3)+(a12*t.^4));
end

STEPS :-

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

Function to find Entropy

function    S = Entropy(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,R,t,mt)
 if t > mt
    S = R*((a1*(log(t)))+(a2*t)+((a3*(t.^2))/2)+((a4*(t.^3))/3)+((a5*(t.^4))/4)+a7);
else
    S = R*((a8*(log(t)))+(a9*t)+((a10*(t.^2))/2)+((a11*(t.^3))/3)+((a12*(t.^4))/4)+a14);
end

STEPS:-

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

Function to find Enthalpy

function   H = Enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,R,t,mt)
if t > mt
    H = R*t.*(a1+((a2*t)./2)+((a3*(t.^2))./3)+((a4*(t.^3))./4)+((a5*(t.^4))./5)+(a6./t));
else
    H = R*t.*(a8+((a9*t)./2)+((a10*(t.^2))./3)+((a11*(t.^3))./4)+((a12*(t.^4))./5)+(a13./t));
end

STEPS:-

Now we will find the Enthalpy by creating the Enthalpy function and assigning its related formula coefficient from 14 coefficients, Gas constant, Temperature, and mid-temperature.
Now we will create a condition using an if loop and make the temperature less than and equal to mid-temperature to specify a local lower temperature coefficient in the formula.
Then we will use the else condition to specify the 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(component)

    % Defining Atomic name of the species
    atomic_name = ['H','S','C','O','N','Ar'];

    % Defining Atomic masses of main species
    atomic_weight = [1.00797 32.065 12.011 15.9994 14.0067 39.948];
    
    W=0;  % Molecular weight starting range
  
    for k=1:length(component)  
     
        for j=1:length(atomic_name)          
            if strcmp(component(k), atomic_name(j))  
                W = W + atomic_weight(j);
                position = j;
            end
        end
        n = str2num(component(k)); 
        if n>1
            W = W + atomic_weight(position)*(n-1);
        end
    end    

end

STEPS:-

Firstly we create a Molecular_weight function and assigning species input to that.
Then define the main element of all 53 species in an array atomic_name that are:- 'H','S', 'C', 'O',' N', 'Ar'.
Then by taking the above main elements, we will create their atomic_wieght array from any source like the internet that is:- 1.008 32.065 12.011 15.999 14.007 39.948.
Then we will assign molecular starting range W To zero and then create a for loop condition which runs to the length of the species and another nested for loop is created which runs to the length of the string atomic_name.
Then we will create a condition using the if condition and use strcmp command to compare the component and atomic_name strings to calculate molecular masses if, by comparison, both strings found the same then we will assign its molecular weight.
Now, we will string to the number using str2num command assign to n. Then we run a condition of if loop to find whether the species contains more than one element if there is more than 1 element that is n>1 then we will increment the atomic weight by the multiplication of position and n-1 to the atomic_weight array.

FINAL RESULTS:-

PLOTS:-

Plots of component O₂

_{Plots of component N2}

_{Plots of component CO2}

_{Screenshot of the Directory where Different Component is saved}

_{The output of the Molecular weight of the Different Component}

References:-

https://www.mathworks.com/matlabcentral/answers/317006-create-new-folders-by-mkdir

https://www.mathworks.com/matlabcentral/answers/17580-how-to-define-a-path-in-saveas-command

https://www.mathworks.com/matlabcentral/answers/106320-make-a-loop-to-calculate-molecular-weight

Parsing NASA thermodynamic data

Ship work like this
in your own course.

Talk to a career counsellor

Ship work like thisin your own course.

Talk to a career counsellor

Ship work like this
in your own course.