Project 1 - Parsing NASA thermodynamic data

Aim: Write a function that extracts the 14 coefficients and calculates the enthalpy, entropy, and specific heats for all the species in the data file. You will be reading this file NASA thermodynamic data(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.

Also use the 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. 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 a suitable name. All of these should be generated automatically.
5. The main program that you will submit will just take the "THERMO.dat" file as an input and generate Cp, H, and S plots.

DESCRIPTION

Functions for file parsing are as follows:

FGETL

• Read line from the file, discard the newline character.
• TLINE = FGETL(FID) returns the next line of a file associated with file identifier FID as a MATLAB character vector. The line terminator is NOT included.
• Use FGETS to get the next line with the line terminator INCLUDED. If just an end-of-file is encountered, -1 is returned. If an error occurs while reading from the file, FGETL returns an empty character vector.
• FGETL is intended for use with files that contain newline characters. Given a file with no newline characters, FGETL may take a long time to execute.

FGETS

• Read line from file, keeping the newline character.
• TLINE = FGETS(FID) returns the next line of a file associated with file identifier FID as a MATLAB character vector. The line terminator is included. Use FGETL to get the next line WITHOUT the line terminator. If just an end-of-file is encountered then -1 is returned.
• If an error occurs while reading from the file, FGETS returns an empty character vector. Use FERROR to determine the nature of the error.
• TLINE = FGETS(FID, NCHAR) returns at most NCHAR characters of the next line. No additional characters are read after the line terminator(s) or an end-of-file.
• MATLAB reads characters using the encoding scheme associated with the file.

FOPEN

• Open file.
• FID = FOPEN(FILENAME) opens the file FILENAME for read access. FILENAME is the name of the file to be opened.
• FILENAME can be a MATLABPATH relative partial pathname. If the file is not found in the current working directory, FOPEN searches for it on the MATLAB search path. On UNIX systems, FILENAME may also start with a "~/" or a "~username/”, which FOPEN expands to the current user's home directory or the specified user's home directory, respectively.
• FID is a scalar MATLAB integer-valued double, called a file identifier.

FCLOSE

• fclose('all')closes all open files.
• status = fclose(___)returns a status of 0 when the close operation is successful.

FERROR

• File I/O error information.
• ferror(fileID,'clear') clears the error indicator for the specified file. Specify the literal clear as a character vector or a string scalar.

FPRINTF

• fprintf(fileID,formatSpec,A1,...,An)applies theformatSpec to all elements of arrays A1,...An in column order, and writes the data to a text file. fprintf uses the encoding scheme specified in the call to

FREAD

• Read data from binary file
• A= fread(fileID) reads data from an open binary file into column vector A and positions the file pointer at the end-of-file marker. The binary file is indicated by the file identifier, fileID. Use fopen to open the file and obtain the fileID When you finish reading, close the file by calling fclose(fileID).
• A= fread(fileID,sizeA) reads file data into an array, A, with dimensions, sizeA, and positions the file pointer after the last value read. fread populates A in column order.

FSCANF

• Read data from text file.
• A = fscanf(fileID,formatSpec)reads data from an open text file into column vector A and interprets values in the file according to the format specified by formatSpec. The fscanffunction reapplies the format throughout the entire file and positions the file pointer at the end-of-file marker. If fscanfcannot match formatSpec to the data, it reads only the portion that matches and stops processing.
• The text file is indicated by the file identifier,fileID. Use fopento open the file, specify the character encoding, and obtain thefileID When you finish reading, close the file by calling fclose(fileID).

FWRITE

• Write data to binary file.
• fwrite(fileID,A)write the elements of array A as 8-bit unsigned integers to a binary file in column order. The binary file is indicated by the file identifier, fileID. Use fopen to open the file and obtain the fileID When you finish reading, close the file by calling fclose(fileID).
• fwrite(fileID,A,precision)writes the values in A in the form and size described by precision.

OBJECTIVE

• To generate plots of Cp v/s temperature, h v/s temperature, and s v/s temperature.

Figure: Flow of the project.

CODE

% NASA Thermodynamic data File Parsing

clear all

close all

clc

% Opening the Data file named THERMO.dat.

f1 = fopen('THERMO.dat','r');%opening the file

Header= fgetl(f1);% Reading header

Global_temperature = fgetl(f1);% Reading Global Temperature

A = strsplit(Global_temperature,' ');

global_low_temperature =str2num(A{2}); % Reading the Global Low Temperature

global_mid_temperature =str2num(A{3}); % Reading the Global Mid Temperature

global_high_temperature =str2num(A{4}); % Reading the Global High Temperature

% Skipping 3 lines

fgetl(f1);

fgetl(f1);

fgetl(f1);

for i = 1:53%initialization of for loop

line= fgetl(f1); %reading the line

C=strsplit(line,' ');

species=C(1); %Naming the First variable as 'Species'

Q = strfind(line,'G');

local_low_temperature = line(Q+3:55);

local_high_temperature = line(Q+13:65);

local_mid_temperature = line(Q+23:75);

lt = str2double(local_low_temperature); %Local Low Temperature

ht = str2double(local_high_temperature); %Local Mid Temperature

mt = str2double(local_mid_temperature); %Local High Temperature

% Reading the coefficient data of 1st line

coeff_line1= fgetl(f1);

X = strfind(coeff_line1,'E');

a1 = coeff_line1(1:X(1)+3);

a1 = str2double(a1); %storing 1st coefficient

a2 = coeff_line1(X(1)+4:X(2)+3);

a2 = str2double(a2); %storing 2nd coefficient

a3 = coeff_line1(X(2)+4:X(3)+3);

a3 = str2double(a3); %storing 3rd coefficient

a4 = coeff_line1(X(3)+4:X(4)+3);

a4 = str2double(a4); %storing 4th coefficient

a5 = coeff_line1(X(4)+4:X(5)+3);

a5 = str2double(a5); %storing 5th coefficient

cfline2 = fgetl(f1);% Reading coefficient of 2nd line

G = strfind(cfline2,'E');

a6 = cfline2(1:G(1)+3);

a6 = str2double(a6); %storing 6th coefficient

a7 = cfline2(G(1)+4:G(2)+3);

a7 = str2double(a7); %storing 7th coefficient

a8 = cfline2(G(2)+4:G(3)+3);

a8 = str2double(a8); %storing 8th coefficient

a9 = cfline2(G(3)+4:G(4)+3);

a9 = str2double(a9); %storing 9th coefficient

a10 = cfline2(G(4)+4:G(5)+3);

a10 = str2double(a10); %storing 10th coefficient

cfline3 = fgetl(f1);% Reading coefficient data of 3rd line

I = strfind(cfline3,'E');

a11 = cfline3(1:I(1)+3);

a11 = str2double(a11); %storing 11th coefficient

a12 = cfline3(I(1)+4:I(2)+3);

a12 = str2double(a12); %storing 12th coefficient

a13 = cfline3(I(2)+4:I(3)+3);

a13 = str2double(a13); %storing 13th coefficient

a14 = cfline3(I(3)+4:I(4)+3);

a14=str2double(a14); %storing 14th coefficient

R = 8.314 ;% Universal Gas Constant

T=linspace(lt,ht,500);% Defining Temperature Range low to high

% Calculation of Properties

if T>=mt

Cp=(a1+a2.*(T)+a3.*(T.^2)+a4.*(T.^3)+a5.*(T.^4))*R;% Specific Heat

else

Cp=(a8+a9.*(T)+a10.*(T.^2)+a11.*(T.^3)+a12.*(T.^4))*R;

end

if T>=mt

S=(a1.*log(T)+a2.*T+a3.*(T.^2)./2+a4.*(T.^3)./3+a5.*(T.^4)./4+a7)*R;% Entropy

else

S=(a8.*log(T)+a9.*T+a10.*(T.^2)./2+a11.*(T.^3)./3+a12.*(T.^4)./4+a14)*R;

end

if T>=mt

H=(a1.*T+a2.*(T.^2)./2+a3.*(T.^3)./3+a4.*(T.^4)./4+a5.*(T.^5)./5+a6)*R;% Enthapy

else

H=(a8.*T+a9.*(T.^2)./2+a10.*(T.^3)./3+a11.*(T.^4)./4+a12.*(T.^5)./5+a13)*R;

end

% Saving the files in the folder

spec=char(species);

fd = sprintf('%s',spec);

cur_dir=pwd;

mkdir(spec);

cd(spec);

%% Plotting

% Plot 1st for Specific heat Cp (vs) Temperature T

figure(1);

plot(T,Cp,'linewidth',3,'color','r');

xlabel('Temperature in K');

ylabel('Specific heat in kJ');

title('Specific heat vs Temperature');

% Plot 2nd for Enthalpy H (vs) Temperature T

figure(2);

plot(T,H,'linewidth',2,'color','b');

xlabel('Temperature in K');

ylabel('Enthalpy in kJ/(mol)');

title('Enthalpy vs Temperature');

%Plot 3rd for Entropy (S) vs Temperature (T)

figure(3);

plot(T,S,'linewidth',2,'color','g');

xlabel('Temperature in K');

ylabel('Entropy in kJ/(mol K)');

title('Entropy vs Temperature');

%% Saving the generated plots

saveas(figure(1), 'cp.png');

saveas(figure(2), 'h.png');

saveas(figure(3), 's.png');

cd(cur_dir);

%% Molecular Weight Calculation

MW= molecular_weight(species)

f2 = fopen('molecular weight','w')

fprintf('%s \n %s \n %f \n','species','Molecular weight ',MW)

fclose(f2)

end

fclose(f1);

Function for calculating Molecular Weight

% Function_for_Molecular_Weight

function M= molecular_weight(species)

atoms={'H','C','N','O','AR'};

AtomicWeight=[1 12 14 16 40];

M=0;

for i=1:length(species)

for j=1:length(AtomicWeight)

if (strcmp(species(i),atoms(j))==1)

M=M+AtomicWeight(j);

g=j;

end

end

% Adding the weights of the repeating atoms

n=str2double(species(i));

if n>1

p=(n-1);

M=M+AtomicWeight(p)*(n-1);

end

end

end

EXPLANATION OF CODE

• File ‘THERMO.dat’ is opened using fopen command.
• Read the file by using fgetl command and found the coefficients.
• String functions are used to store the values of coefficients.
• Values of coefficients are used to calculate the following equations
• Values of Cp , H, and S are used to plot graphs by using plot command.
• Pictures in png format are saved by using saveas command.
• Function to find the molecular weight is written and then called in the main program ton find the molecular weight.
• Molecular weight is printed by using fprint command

OUTPUT

• O₂

N₂

CO₂

Figure: Saved Files.

CONCLUSIONS

• Objectives are achieved successfully.

REFERENCES

https://in.mathworks.com/help/documentation-center.html