[PYTHON] Run NASA CEA on Mac OS X
Operations required to execute NASA CEA with commands. Windows is relatively easy, but Mac is a bit annoying. I tried Mac book air, OSX yosemite For those who know what CEA is.
If you can read English, please refer to the PDF and youtube below. It's all written. http://spase.stanford.edu/CEA_files/cea_how_to_gfortran_snowleo.pdf https://www.youtube.com/watch?v=CGpxfxmqA-c
gfortran Install XCode and Command Line Tools to use GCC. Now you can compile Fortran with the command gfortran. http://qiita.com/ShinichiOkayama/items/962e55ac7f8920a5f007
CEA download and path
Go to this site at NASA Glenn Research Center. http://www.grc.nasa.gov/WWW/CEAWeb/ceaguiDownload-unix.htm Download the following file from here and move it to the folder you want to install.
- CEAgui JAR (CEAgui-jar.tar.Z)
- CEA+Fortran Package (CEA+Fortran.tar.Z)
- CEAexec Package (CEAexec-mac.tar.Z)
I created a folder called ~ / Applications / NASAcea.
Now start the terminal. Go to ~ / Applications / NASAcea on the terminal and
zcat CEA+Fortran.tar.Z | tar xvf -
zcat CEAexec-mac.tar.Z | tar xvf -
zcat CEAgui-jar.tar.Z | tar xvf -
And defrost. In addition, enter the following as file permission settings.
chmod a+x b1b2b3
chmod a+x syntax
chmod a+x FCEA2
chmod a+x runCEA.sh
Now, pass the path to the folder so that it can be executed.
echo 'export PATH=$HOME/Applications/NASAcea:$PATH' >> ~/.bash_profile
bash -login
Also used by other users
To make it available to other users, do the following. Not required for personal use.
sudo vi /etc/paths
So add the following to the last line of this file:
~/Applications/NASAcea
CEA file update
Update the executable file. Enter the following at the terminal.
gfortran cea2.f
mv a.out FCEA2
gfortran b1b2b3.f
mv a.out b1b2b3
gfortran syntax.f
mv a.out syntax
After that, enter the commands in order to update the library. (Instead of copying and pasting continuously, insert line by line). [Addition] FCEA2 is in the path, and when executing in another directory, it is necessary to put it in the directory where thrmo.inp, trans.inp, cea2.inp is executed in advance. After executing the following command, it is normal if the trans.lib, thermo.lib, (cea2.lib) file contains unreadable characters in the binary file. In the case of files from, I haven't said it well.
FCEA2
trans
FCEA2
thermo
FCEA2
cea2
That's it.
sample
As a test, let's show the performance of the Kerosin / LOX rocket engine. Create a text file called sample.inp. The extension .inp is important. The contents are as follows.
sample.inp
problem o/f=2.5,
rocket equilibrium frozen nfz=2
p,bar=100,
pi/p=100,
react
fuel=RP-1 wt=100 t,k=298.15
oxid=O2(L) wt=100 t,k=90.17
output
plot o/f isp aeat cf t
end
Now go to the folder with sample.inp on the terminal and
FCEA2
sample
This second line contains the name part of the inp file. This is the result because a file called sample.out appears in the same folder.
If you don't know how to create a .inp file for CEA, write what you want to do with the online version of CEA and click on it. Finally, the .inp file will appear along with the result file, so refer to that. https://cearun.grc.nasa.gov/
reference
There is an example near the end of the paper NASA-RP1311. http://www.frad.t.u-tokyo.ac.jp/public/cea/doc/xRP-1311P2.pdf For example, rockets are in EXAMPLE 8 ~ 14. The part from problem to end after the comment # should be an .inp file.
Scripting
The results produced by CEA are old-fashioned and difficult to summarize, so it is necessary to summarize them. Try it with a python script. Some people are making CEA interfaces in Matlab https://github.com/PurdueH2Lab/MatlabCEA Rather than creating a decent interface, hitting it directly with a script would be a straightforward approach. For example, let's calculate the performance of a liquid oxygen / liquid hydrogen rocket engine. You can use it by changing the file name or .inp file part as appropriate. Put the following .py file in the same folder as FCEA2 and execute .py from the terminal. I moved to the FCEA2 folder (~ / Applications / NASAcea) with ipython
run LH2LO2.py
Is running.
LH2LO2.py
# -*- coding: utf-8 -*-
import sys
reload(sys)
#Change the default character code.
sys.setdefaultencoding('utf-8')
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
from matplotlib.font_manager import FontProperties
from scipy import optimize
import csv
import os
import subprocess
# for Mac
font_path = '/Library/Fonts/Osaka.ttf'
font_prop = matplotlib.font_manager.FontProperties(fname=font_path)
matplotlib.rcParams['font.family'] = font_prop.get_name()
def cea (set_of, set_chamber):
file_name = '00temp' #Temporarily created file
# ---- .inp file creation----
input_file = file_name + '.inp'
str = """problem o/f=%.1f,
rocket equilibrium frozen nfz=2
p,bar = %d
pi/p= %d
react
fuel=H2(L) wt=100 t,k=20.27
oxid=O2(L) wt=100 t,k=90.17
output
plot p pi/p o/f isp ivac aeat cf t
end
""" % (set_of, set_chamber, set_chamber)
f = open(input_file, 'w') #Open in write mode
f.write(str) #Write the argument string to a file
f.close() #Close file
# ----FCEA2 command----
#The location of FCEA2~/Applications/NASAcea/FCEA2 is used, but rewrite as appropriate
cmd = '~/Applications/NASAcea/FCEA2'
p = subprocess.Popen(cmd, shell=True, stdin=subprocess.PIPE, stdout=subprocess.PIPE)
p.communicate(file_name)
# ---- .plt file read----
output_file = file_name + '.plt'
with open(output_file, 'r') as f:
reader = csv.reader(f, delimiter=' ')
header = next(reader) #Skip header
index = 0
for row in reader:
while row.count("") > 0:
row.remove("")
# print row #List display line by line
if index == 0:
chamber = float(row[0])
elif index == 2:
of = float(row[2])
isp = float(row[3])
ivac = float(row[4])
aeat = float(row[5])
cf = float(row[6])
index += 1
data = [chamber, of, isp, ivac, aeat, cf]
#File deletion
os.remove(input_file)
os.remove(output_file)
return data
if __name__ == '__main__': #main (starts here)
data_file = '00output.csv'
plt.ion()
plt.figure()
# ----Parameter specification----
of_array = np.linspace(3.0, 6.0, 7)
# chamber_array = np.linspace(5, 40, 30)
chamber_array = np.linspace(10, 200, 30)
data_array = np.zeros((len(chamber_array), 6))
save_array = np.zeros([1,6]) #Stock location for data to write to CSV
# ----Loop the command----
#In this example, the combustion chamber pressure is on the horizontal axis, O/I want to PLOT with F on the vertical axis, so O/Loop in the order of F → pressure
for j in of_array:
k = 0
for i in chamber_array:
data = cea(j, i)
print data
data_array[k] = np.array(data)
k+=1
g = 9.80665
press = data_array[:,0] / g #Combustion chamber pressure MPa
of = data_array[:,1] # O/F
Isp = data_array[:,2] / g #Maritime specific impulse sec
Ivac = data_array[:,3] / g #Specific impulse in vacuum sec
Cf = data_array[:,4] #Thrust coefficient Cf
# ---- PLOT ----For example, PLOT the combustion chamber pressure and the specific impulse at sea.
plt.plot(press, Isp, label='O/F=%.1f' % (j))
# ---- ----
save_array = np.vstack([save_array, data_array])
# ----Save results in CSV----
header = u'chamber[MPa]\tof[-]\tisp[sec]\tivac[sec]\taeat[-]\tcf[-]'
np.savetxt(data_file, save_array, fmt='%.3f', delimiter='\t', header = header)
# ----PLOT settings----
plt.xlim(0,20)
plt.xlabel('Pressure (MPa)')
plt.ylabel('Isp (sec)')
plt.grid()
plt.title('LOX/LH2 100% equilibrium frozen nfz=2')
plt.legend(loc='best', fontsize=12)
I get this result.
Recommended Posts