Project 1 - Parsing NASA thermodynamic data
Parsing NASA thermodynamic data
AIM:-
1. To write a Matlab code that extracts the 14 coefficients and calculates the enthalpy (H), entropy (S), and specific heats (Cp) for all the species in the data (THERMO.dat) file. The FIRST 7 coefficients are HIGH-temperature coefficients and the SECOND 7 coefficients are LOW-temperature coefficients. Only use the LOCAL TEMPERATURE RANGES which are unique for each species, not the GLOBAL TEMPERATURE RANGE.
2. To calculate the molecular weight of each species and display it in the command window.
3. Plot the specific heat, Enthalpy, and Entropy for the local temperature range, 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 these should be generated automatically.
4. The main program will take the "THERMO.dat" file as an input and generate Cp, H, and S plots, and write a detailed report to describe the logic used in the program.
Principle:-
File Parsing is how we read information from a text file or any type of file and how it can be used in any computer program. Parsing data helps us turn unstructured data/un-readable data into structured data. Data Parsers are also used to convert data from one type to another. In simple terms, parsing data can help turn unstructured data or unreadable data into more readable and valuable information. For example, let's assume we have an HTML file, and we would need to convert it into a PDF file. In this case, by parsing (reads) through the HTML file, extracts the necessary information and exports it to the PDF file. Similarly, there are several problems that a data parser could solve when dealing with massive amounts of data.
Enthalpy (H) -
Enthalpy is 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 chemical, biological, and physical measurements at a constant pressure, which is conveniently provided by the large ambient atmosphere. The pressure-volume term expresses the work required to establish the system's physical dimensions, i.e. to make room for it by displacing its surroundings. The pressure-volume term is very small for solids and liquids under common conditions and fairly small for gases. Therefore, enthalpy is a stand-in for energy in chemical systems; bonds, lattices, and other "energies" in chemistry are actually enthalpy different. As a state function, enthalpy depends only on the final configuration of internal energy, pressure, and volume, not on the path taken to achieve it.
Formula to apply:
Here R is the universal gas constant, T is the temperature and a1-a6 are coefficients from the data file.
Entropy (S) -
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 firrecognizedsed, to the microscopic description of nature in statistical physics, and the principles of information theory. It has found far-ranging applications in chemistry and physics, in biological systems and their relation to life, cosmology, economics, sociology, weather science, climate change, and information systems including the transmission of information in telecommunication. Entropy is central to the second law of thermodynamics, which states that the entropy of isolated systems left to spontaneous evolution cannot decrease with time, as they always arrive at a state of thermodynamic equilibrium, where the entropy is highest.
Formula to apply:
Here R is the universal gas constant, T is the temperature and a1-a7 are coefficients from the data file.
Specific Heat (Cp) -
In thermodynamics, the specific heat capacity (Cp) of a substance is the heat capacity of a sample of the substance divided by the mass of the sample. Specific heat is also sometimes referred to as massic heat capacity. Informally, it is the amount of heat that must be added to one unit of mass of the substance to cause an increase of one unit in temperature. The specific heat capacity SI unit is joule per kelvin per kilogram, J⋅kg−1⋅K−1.
Specific heat capacity often varies with temperature and is different for each state of matter. Liquid water has one of the highest specific heat capacities among common substances, about 4184 J⋅kg−1⋅K−1 at 20 °C; but that of ice, just below 0 °C, is only 2093 J⋅kg−1⋅K−1. The specific heat capacities of iron, granite, and hydrogen gas are about 449 J⋅kg−1⋅K−1, 790 J⋅kg−1⋅K−1, and 14300 J⋅kg−1⋅K−1, respectively. While the substance is undergoing a phase transition, such as melting or boiling, its specific heat capacity is technically infinite, because the heat goes into changing its state rather than raising its temperature.
Formula to apply:
Here R is the universal gas constant, T is the temperature and a1-a4 are coefficients from the data file.
Commands used -
1. fopen
Open the file, or obtain information about open files
Syntax
fileID = fopen(filename)
Description
fileID = fopen(filename) opens the file, filename, for binary read access, and returns an integer file identifier equal to or greater than 3. MATLAB® reserves file identifiers 0, 1, and 2 for standard input, standard output (the screen), and standard error, respectively.
If fopen cannot open the file, then fileID is -1.
2. fgetl
Read line from the file, removing newline characters
Syntax
tline = fgetl(fileID)
Description
tline = fgetl(fileID) returns the next line of the specified file, removing the newline characters.
- If the file is nonempty, then fgetl returns tline as a character vector.
- If the file is empty and contains only the end-of-file marker, then fgetl returns tline as a numeric value -1.
3. strsplit
Split string or character vector at the specified delimiter
Syntax
C = strsplit(str,delimiter)
Description
C = strsplit(str,delimiter) splits str at the delimiters specified by delimiter.
If str has consecutive delimiters, with no other characters between them, then strsplit treats them as one delimiter. For example, both strsplit('Hello,world',',') and strsplit('Hello,,,world',',') return the same output.
4. str2double
Convert strings to double precision values
Syntax
X = str2double(str)
Description
X = str2double(str) converts the text in str to double-precision values. str contains text that represents real or complex numeric values. str can be a character vector, a cell array of character vectors, or a string array. If str is a character vector or string scalar, then X is a numeric scalar. If str is a cell array of character vectors or a string array, then X is a numeric array that is the same size as str.
Text that represents a number can contain digits, a comma (thousands separator), a decimal point, a leading + or - sign, an e preceding power of 10 scale factor, and an i or a j for a complex unit. You cannot use a period as a thousand separator, or a comma as a decimal point.
If str2double cannot convert text to a number, then it returns a NaN value.
5. mkdir
Make new folder
Syntax
mkdir folderName
Description
mkdir folderName creates the folder folderName. If folderName exists, MATLAB® issues a warning. If the operation is not successful, mkdir throws an error to the Command Window.
6. cd
Change current folder
Syntax
cd newFolder
Description
cd newFolder changes the current folder to newFolder. Folder changes are global. Therefore, if you use cd within a function, the folder change persists after MATLAB® finishes executing the function.
7. saveas
Save figure to specific file format
Syntax
saveas(fig,filename)
Description
saveas(fig, filename) saves the figure or Simulink® block diagram specified by fig to file filename. Specify the file name as a character vector or string that includes a file extension, for example, 'myplot.jpg'. The file extension defines the file format. If you do not specify an extension, then saveas saves the figure to a FIG file. To save the current figure, specify fig as gcf.
Program:-
Two functions are being created for the main program to calculate thermodynamic properties like specific heat(Cp), enthalpy(H), and entropy(S), and another function is for calculating the molecular weight of the species in the THERMO.dat file.
A function name ‘thermodynamic_properties’ is created for Cp, H, and S as output. 3 different functions can also be created separately for Cp, H, and S.
% Function to calculate thermodynamic properties | specific heat(Cp) | enthalpy(H) | entropy(S)
function [Cp, H, S] = thermodynamic_properties(C1, C2, C3, C4, C5, C6, C7, C8, C9, C10, C11, C12, C13, C14, R, T, LTm)
% Condition of temperature to consider high-temperature and low-temperature coefficients
if T <= LTm
% Low Temperature Coefficients
Cp = (C8 + (C9.*T) + (C10.*T.^2) + (C11.*T.^3) + (C12.*T.^4)).*R;
H = (C8 + (C9.*(T./2)) + (C10.*(T.^2./3)) + (C11.*(T.^3./4)) + (C12.*(T.^4./5)) + (C13./T)).*(R.*T);
S = ((C8.*log(T)) + (C9.*T) + (C10.*T.^2./2) + (C11.*T.^3./3) + (C12.*T.^4./4) + C14).*R;
else
% High Temperature Coefficients
Cp = (C1 + (C2.*T) + (C3.*T.^2) + (C4.*T.^3) + (C5.*T.^4)).*R;
H = (C1 + (C2.*(T./2)) + (C3.*(T.^2./3)) + (C4.*(T.^3./4)) + (C5.*(T.^4./5)) + (C6./T)).*(R.*T);
S = ((C1.*log(T)) + (C2.*T) + (C3.*T.^2./2) + (C4.*T.^3./3) + (C5.*T.^4./4) + C7).*R;
end
end- In this function, an ‘if’ condition is given to consider low-temperature coefficients and high-temperature coefficients based on temperature value as the T value varies for every species. If the temperature value follows the condition which is if T is less than or equal to LTm(local temperature mid-value) then low-temperature coefficients that are from C8-C14 will be considered for calculating thermodynamic properties or else high-temperature coefficients that are from C1-C7 will be considered for calculation.
In the below function, molecular weight is calculated for all the species of the data file where the required output is molecular weight(M_W) and input is species.
% Function to calculate molecular weight
function [M_W] = molecular_weight(species_name)
% An array of all the atoms in every species and their atomic weight
atomic_names = ['H', 'C', 'O', 'N', 'Ar', 'S'];
atomic_weight = [1.00784, 12.0107, 15.999, 14.0067, 39.948, 32.065];
% Taking intial molecular weight as 0
M_W = 0;
% For loop to repeat the same steps for every species and atoms to get their molecular weight
% Length of the species or molecules
for i = 1:length(species_name)
% Length of the array for atomic names
for j = 1:length(atomic_names)
% Condition to get atomic weight of the same atom for every species
if species_name(i) == atomic_names(j)
% Molecular weight of species with single atomic bond of an atom
M_W = M_W + atomic_weight(j);
% Position of an atom from the array(atomic_names)
% atomic weight of the same position will be consider
position = j;
end
end
% Converting no. of atomic bonds of an atom into number
n = str2double(species_name(i));
% Condition for a species to follow for further calculation
if n>1
% Molecular weight of species with multiple atomic bonds of an atom
M_W = M_W + atomic_weight(position)*(n-1);
end
end
% Command to print molecular weight of a species
fprintf('Molecular Weight of %s = %f', species_name, M_W)
disp(' ')
end
- An array is created for all the atoms and atomic weights that are present in every species. Considering initial molecular weight(M_W) as 0. Nested for loop is used where the first loop is about the length of every species and the second loop is about the length of the array for atoms. One ‘if’ condition is given for the species name and the atoms to be the same, if it satisfies the condition then M_W is calculated for that particular species.
- Point out the position of the current atom, as it’ll be considered in case further steps are satisfied. If an atom has more than one bond then the no. of bonds is converted into numbers under variable ’n’. If the n value is greater than 1 then the molecular weight of that particular atom is calculated. Finally adding both the M_W we get the molecular weight of a molecule species.‘fprintf’ command is used to print the statement inside the command with the species name and M_W. Giving space after that will make the new statement start from the new line.
The main program is about extracting the 14 coefficients and calculating the enthalpy, entropy, and specific heats for all the species in the data file. Also, calculate the molecular weight of every species and plot the Cp, Enthalpy, and Entropy for the local temperature range specified for each species. At the end save the plots as images with appropriate names and save them in separate folders for each species with suitable names.
% File (THERMO.dat) Parsing
clear all
close all
clc
% Opening the data file and extracting header
f1 = fopen('THERMO.dat', 'r');
Head = fgetl(f1);
Head_1 = fgetl(f1);
% Leaving 3 comments line
for i = 1:3
fgetl(f1);
end
% Extracting 53 species Data
% Creating loop for all the species to repeat the same steps to get data
for j = 1:53
% Extracting local temperature for the species
species = fgetl(f1);
% Splitting the line with space to get local temperature values
LT = strsplit(species, ' ');
% Naming the temperature values from LT line and converting to numbers
species_name = LT{1};
LTl = str2double(LT{end - 3});
LTm = str2double(LT{end - 1});
LTh = str2double(LT{end - 2});
% Getting 1st line of coefficents for the species
species_l1 = fgetl(f1);
% seperating the coefficient string from E as a common factor
l = strfind(species_l1, 'E');
% Naming the coefficients after seperating the string from species_l1 line and converting into numbers
C1 = species_l1(1: l(1)+3);
C1 = str2double(C1);
C2 = species_l1(l(1)+4: l(2)+3);
C2 = str2double(C2);
C3 = species_l1(l(2)+4: l(3)+3);
C3 = str2double(C3);
C4 = species_l1(l(3)+4: l(4)+3);
C4 = str2double(C4);
C5 = species_l1(l(4)+4: l(5)+3);
C5 = str2double(C5);
% Getting 2nd line of coefficents for the species
species_l2 = fgetl(f1);
% seperating the coefficient string from E as a common factor
l1 = strfind(species_l2, 'E');
% Naming the coefficients after seperating the string from species_l2 line and converting into numbers
C6 = species_l2(1: l1(1)+3);
C6 = str2double(C6);
C7 = species_l2(l1(1)+4: l1(2)+3);
C7 = str2double(C7);
C8 = species_l2(l1(2)+4: l1(3)+3);
C8 = str2double(C8);
C9 = species_l2(l1(3)+4: l1(4)+3);
C9 = str2double(C9);
C10 = species_l2(l1(4)+4: l1(5)+3);
C10 = str2double(C10);
% Getting 3rd line of coefficents for the species
species_l3 = fgetl(f1);
% seperating the coefficient string from E as a common factor
l2 = strfind(species_l3, 'E');
% Naming the coefficients after seperating the string from species_l3 line and converting into numbers
C11 = species_l3(1: l2(1)+3);
C11 = str2double(C11);
C12 = species_l3(l2(1)+4: l2(2)+3);
C12 = str2double(C12);
C13 = species_l3(l2(2)+4: l2(3)+3);
C13 = str2double(C13);
C14 = species_l3(l2(3)+4: l2(4)+3);
C14 = str2double(C14);
% Universal gas constant (j.mol^-1.K^-1)
R = 8.314;
% Temperature array
T = linspace(LTl, LTh, 500);
% Evaluating thermodynamic properties | specific heat (Cp) | enthalpy (H) | entropy (S)
�unction is created to calculate the properties
[Cp, H, S] = thermodynamic_properties(C1, C2, C3, C4, C5, C6, C7, C8, C9, C10, C11, C12, C13, C14, R, T, LTm);
�lculating molecular weight
% Function is created to calculate molecular weight
[M_W] = molecular_weight(species_name);
% Making folders for each species with particular species name
mkdir(species_name)
% Changing the current folder to the species name
cd(species_name)
% Plot for specific heat (Cp)-(KJ/kg.K) vs temperature (K)
figure(1)
plot(T, Cp, 'LineWidth', 2, 'Color', 'r')
title('Specific Heat (Cp) vs temperature')
xlabel('Local Temperature Range in K')
ylabel('Specific Heat (Cp) in KJ/kg.K')
% Saving figure(1) in current folder
saveas(figure(1), 'specific heat(Cp).png')
% Plot for enthalpy (H)-(KJ/kg) vs temperature (K)
figure(2)
plot(T, H, 'LineWidth', 2, 'Color', 'b')
title('Enthalpy (H) vs Temperature')
xlabel('Local Temperature Range in K')
ylabel('Enthalpy (H) in KJ/kg')
% Saving figure(2) in current folder
saveas(figure(2), 'enthalpy(H).png')
% Plot for entropy (S)-(KJ/K) vs temperature (K)
figure(3)
plot(T, S, 'LineWidth', 2, 'Color', 'b')
title('Entropy (S) vs Temperature')
xlabel('Local Temperature Range in K')
ylabel('Entropy (S) in KJ/K')
% Saving figure(3) in current folder
saveas(figure(3), 'entropy(S).png')
% Changing current folder back to the original folder
cd ..
end- To extract the data from the THERMO.dat file ‘fopen’ command is used to open the file and read. Every time the ‘fgetl’ command is used we get the next line of the data file. At first, we the header of the file then skipped the comment lines of the file by using the ‘fgetl’ command inside a loop.
- In species line we need to extract local temperature values hence using the ‘strsplit’ command, using space to separate all the dates of that line. Saving the species name under species_name and defining the positions of the local temperature values as it varies for every species. So with the help of the positions, we can get the exact value for different species. Extract 14 coefficients after reading the lines and separating them from the position of ‘E’ using the ‘strfind’ command. Then convert them into a number using the ‘str2double’ command.
- Once we extract all the 14 coefficients, using the thermodynamic properties function we calculate specific heat, enthalpy, and entropy. Taking universal gas constant R-value 8.314 and creating an array of local temperatures for 500 values using the ‘linspace’ command for the function to operate. Same way of calculating the molecular weight of the species using the molecular_weight function.
- Now we need to generate a plot for T vs Cp, H, and S and save them in a specific folder with an appropriate name for each species and for the same using the ’mkdir’ command to create a folder with the species name and next change the current folder of the program to the folder of the species name using ‘cd’ command. Then save the figures into the current folder using the ‘saveas’ command. At last, change the current folder back to the original folder for the program to run smoothly using the ‘cd’ command again. To get all these processes done for every species and to get the desired output we create a loop for them of the length of 53 as we have a total of 53 species in the data file.
Results:-
- Screenshot of the window where all the plots are saved automatically for each species in their respective folders.
Figure 1. species folders
- Plots for O2
Figure 2. O2 - enthalpy(H) vs temperature
Figure 3. O2 - entropy(S) vs temperature
Figure 4. O2 - specific heat(Cp) vs temperature
- Plots for N2
Figure 5. N2 - enthalpy(H) vs temperature
Figure 6. N2 - entropy(S) vs temperature
Figure 7. N2 - specific heat(Cp)
- Plots for CO2
Figure 8. CO2 - enthalpy(H) vs temperature
Figure 9. CO2 - entropy(S) vs temperature
Figure 10. CO2 - specific heat(Cp) vs temperature
- The molecular weight of species
Figure 11. molecular weight (I)
Figure 12. molecular weight (II)
- Error faced while programming
The below error has been rectified using single quotes for the filename.
Figure 13. error (I)
The below error occurred due to a missing input in the function which was then fixed.
Figure 14. error (II)
Below is a typo error for the 'linspace' command which was then corrected.
Figure 15. error (III)
Reference:-
MathWorks - https://in.mathworks.com/help/matlab/ref/fopen.html?s_tid=doc_ta
Enthalpy - https://en.wikipedia.org/wiki/Enthalpy
Entropy - https://en.wikipedia.org/wiki/Entropy
Specific Heat - https://en.wikipedia.org/wiki/Specific_heat_capacity
Molecular Weight - https://en.wikipedia.org/wiki/Molecular_mass