Challenge

Project 1 - Parsing NASA thermodynamic data

AIM : 1. Write a function that extracts the 14 co-efficients and calculates the enthalpy, entropy and specific heats for all the species in the data file.

2. Calculate the molecular weight of each species and display it in the command window.

3. Plot the Cp, Enthalpy and Entropy for the local temperature range (low temperature : high temperature) specific for each species.

4. Save the plots as images with appropriate names and dump them in separate folders for each species with suitable name. All these should be generated automatically.

5. Generate Cp, H and S plots.

THEORY :

FILE PARSING : File Parsing is any information that reads from a particular text file or any type of file and how that information can be used in a computer program.

FILE WRITING : It is the process of writing data to a particular file with the write mode ‘w’.

FILE READING : File reading is kind a reverse of file writing. Let us say, someone gives you simulation data in the form of an input file, how would we read that data and how would we take the content of that data and use it. It is written with the mode ‘r’.

NASA THERMODYNAMIC DATA : NASA came up with polynomials that can be used 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 but here we are dealing with only 53 species.

SPECIFIC HEAT : Specific heat is the amount of heat required to increase the temperature of 1kg of substance through 1 deg Celsius.

Cp = (a1 + a2T + a3T^2 + a4T^3 + a5T^4)*R

Function code for Specific heat :

% Specific heat calculation 

function Cp = Specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,Temp,R,global_mid_temp)
% Compairing with Global temperature
% If less than global temperature - minimum coefficients are taken
if Temp <= global_mid_temp
    Cp = R*(a8+(a9.*Temp)+(a10+(Temp).^2)+(a11.*(Temp).^3)+(a12.*(Temp).^4));
else

    Cp = R*(a1+(a2.*Temp)+(a3+(Temp).^2)+(a4.*(Temp).^3)+(a5.*(Temp).^4));
 % Else maximum coefficients are taken
end

ENTROPY : Between two state points, if number of reversible paths are there , then, the quantity (dQ/T)rev for each reversible path is same and this quantity (dQ/T)rev known as Entropy.

Since Entropy is a property, it does not depends on the path followed by the system as long as end points are same, Entropy is also same for reversible as well as Irreversible process.

S = (a1lnT + a2T + a3T^2/2 + a4T^3/3 + a5T^4/4 + a7)*R

Function code for Entropy :

% Calculation of Entropy function

function S = Entropy_calc(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,Temp,R,global_mid_temp)

if Temp <= global_mid_temp

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

ENTHALPY : Enthalpy is the sum of the internal energy. and pressure times volume. We cannot measure the enthalpy of a system, but we can look at changes in enthalpy.

H = (a1 + a2T/2 + a3T^2/3 + a4T^3/4 + a5T^4/5 + a6/T)*R*T

Function code for Enthalpy :

% Calculation of Enthalpy function

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

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

MOLECULAR WEIGHT : The molecular mass is the mass of a given molecule which is measured in dalton. Different molecules of the same compound may have different molecular masses because they contain different isotopes of an element.

Function code for Molecular Weight :

% Calculation of Molecular Weight function

function Molecularweight = Molecular_weight(Species_SP)

% Storing value
compound_name = (Species_SP);

% Creating an array of elements
elements_name = ['C','H','O','N','A'];

% Storing atomic mass value in the array
value_at = [12 1 15.994 14.0067 39.948];

% Assigning molecular mass as 0
Molecularweight = 0;

% Creating loop
for i = 1:length(compound_name)
    for j = 1:length(elements_name)

        % Compairing string value
        if strcmp(compound_name(i),elements_name(j))

            %If above comparision is true, value of molecular mass is
            %calculated
            Molecularweight = Molecularweight + value_at(j);

            % Storing the position for further use
            pos = j;
        end
    end

    % Converting the string to number
    k = str2num(compound_name(i));

    % Compairing if the value of k is greater than executed, the further
    % statements
    if (k>1)

        % Molecular mass calculation
        Molecularweight = Molecularweight + value_at(pos)*(k-1);
    end
end
fprintf('molecular weight %s :',Species_SP);
fprintf('%d',Molecularweight);
disp(' ');
end

MAIN CODE :

clear all
close all
clc

f1 = fopen('THERMO.dat','r')
tline_1 = fgetl(f1)
tline_2 = fgetl(f1)


R = 8.314;      % (Universal Gas Constant in KJ/Kmol-K)

for i = 1:3
    get = fgetl(f1)
end

for i = 1:53
tline_3 = fgetl(f1);
B = strsplit(tline_3,' ');
Species_SP = (B{1});

% Taking all the local temperatures

global_low_temp = str2num(B{end-3});
global_mid_temp = str2num(B{end-1});
global_high_temp = str2num(B{end-2});
Temp = linspace(global_low_temp,global_high_temp, 1000);
line_1 = fgetl(f1)

% Finding or extracting species 

a = findstr(line_1,'E')

a1 = line_1(1:a(1)+3)          %1
a1 = str2num(a1)
a2 = line_1(a(1)+4:a(2)+3)     %2
a2 = str2num(a2)
a3 = line_1(a(2)+4:a(3)+3)     %3
a3 = str2num(a3)
a4 = line_1(a(3)+4:a(4)+3)     %4
a4 = str2num(a4)
a5 = line_1(a(4)+4:a(5)+3)     %5
a5 = str2num(a5)

line_2 = fgetl(f1)
b = findstr(line_2,'E')

a6 = line_2(1:b(1)+3)          %6
a6 = str2num(a6)
a7 = line_2(b(1)+4:b(2)+3)     %7
a7 = str2num(a7)
a8 = line_2(b(2)+4:b(3)+3)     %8
a8 = str2num(a8)
a9 = line_2(b(3)+4:b(4)+3)     %9
a9 = str2num(a9)
a10 = line_2(b(4)+4:b(5)+3)    %10
a10 = str2num(a10)

line_3 = fgetl(f1)
c = findstr(line_3,'E')

a11 = line_3(1:c(1)+3)         %11
a11 = str2num(a11)
a12 = line_3(c(1)+4:c(2)+3)    %12
a12 = str2num(a12)
a13 = line_3(c(2)+4:c(3)+3)    %13
a13 = str2num(a13)
a14 = line_3(c(3)+4:c(4)+3)    %14
a14 = str2num(a14)

% Calculate specific heat
Cp = Specific_heat(a1,a2,a3,a4,a5,a8,a9,a10,a11,a12,Temp,R,global_mid_temp);

% Calculate Entropy
S = Entropy_calc(a1,a2,a3,a4,a5,a7,a8,a9,a10,a11,a12,a14,Temp,R,global_mid_temp);

% Calculate Enthalpy
H = Enthalpy(a1,a2,a3,a4,a5,a6,a8,a9,a10,a11,a12,a13,Temp,R,global_mid_temp);

% Calculate Molecular weight of the Species
Molecularweight = Molecular_weight(Species_SP);



% Plotting(Cp,H,S,Temp,Species_SP)
% Saving a folder


mkdir(Species_SP)
cd(Species_SP)

% Plot 1, Specific heat(Cp) vs Temperature Range(Temp)

figure(1)
plot(Temp,Cp,'LineWidth',3,'Color','r');
xlabel('Temperature [K]');
ylabel('Specific Heat [KJ]');
title(sprintf('Varying Specific heat with Temperature range of %s',Species_SP));
grid on
saveas(figure(1),'Specific Heat.jpg');

% Plot 2, Enthalpy(H) vs Temperature range (Temp)

figure(2)
plot(Temp,H,'LineWidth',3,'Color','b');
xlabel('Temperature [K]');
ylabel('Enthalpy [KJ/mol]');
title(sprintf('Varying Enthalpy with Temperature range of %s',Species_SP));
grid on
saveas(figure(2),'Enthalpy.jpg');

% Plot 3, Entropy(S) vs Temperature range (Temp)

figure(3)
plot(Temp,S,'LineWidth',3,'Color','g');
xlabel('Temperature [K]');
ylabel('Entropy [KJ/mol-K]');
title(sprintf('Varying Entropy with Temperature range of %s',Species_SP));
grid on
saveas(figure(3),'Entropy.jpg');

cd ..

% Molecular weight

M = Molecular_weight(Species_SP);


end

EXPLANATION AND SYNTAX:

  1. First of all, we used clear all, close all and clc command to clear all the previous data.
  2. Loading the ‘THERMO.dat’ file through file reading mode.
  3. Reading and loading each lines step by step through ‘fgetl’ command.
  4. Skipping the three statement lines using for loop.
  5. Running the for loop again for all the 53 species.
  6. Taking and acknowledging all three local temperatures as asked in the question.
  7. Extracting species through ‘findstr’ command and finding all the coefficients while changing it from string to number through ‘str2num’ command.
  8. Calculating all the values for different species of specific heat, entropy and enthalpy through their respective function code.
  9. Plotting the graph for each and every species of Cp, S and H.
  10. And also calculating their respective molecular weight.

Fopen : It is used to open any file.

‘r’ : It is used for reading of any file that is opened.

Fgetl : It is used to extract lines from the file step by step.

Strsplit : Splitting of string at specified delimiter.

Str2num : It is used to convert any string into number.

Findstr : Find one string within another.

Linspace : To generate equally spaced vectors.

Fprintf : Write formatted data to file.

Mkdir : To make new directory.

Cd : To change working directory.

OUTPUT AND RESULTS :

1. Screenshot of all the folders have been attached.

2. Molecular weight of each and every species have been found.

3. Plot for O2, N2 and CO2 species have been attached.

O2 - ENTHALPY

O2 - ENTROPY

O2 - SPECIFIC HEAT

N2 - ENTHALPY

N2 - ENTROPY

N2 - SPECIFIC HEAT

CO2 - ENTHALPY

CO2 - ENTROPY

CO2 - SPECIFIC HEAT

ERRORS :

Not any major errors occured while programming except some minor.

CONCLUSION :

Using MATLAB, all the coefficients of the species have been found, also calculated the molecular weight for different species and have plotted all the graphs for varying enthalpy, entropy and specific heat against the local temperatures.

REFERNECES :

1. https://in.mathworks.com/help/matlab/ref/fprintf.html

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

Your portfolio, next

Ship work like this
in your own course.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.