View this page as an IPython notebook
# This block does the setup
from IPython.core.display import HTML
from IPython.display import display_html
from CoolProp.Plots.Plots import Ph,Ps
from CoolProp.CoolProp import Props, get_fluid_param_string
import CoolProp.CoolProp as CP
from math import sqrt
import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize
%matplotlib inline
labels = {'P' : '$p$',
'D' : r'$\rho$',
'C' : '$c_p$',
'O' : '$c_v$',
'H' : '$h$',
'S' : '$s$',
'V' : '$\mu$',
'L' : '$\lambda$'}
sat_labels = labels.copy()
sat_labels['I'] = '$\sigma$'
Fluid = 'n-Butane'
RPFluid = 'REFPROP-'+get_fluid_param_string(Fluid,'REFPROPName')
Tt = Props(Fluid,'Tmin')
Tc = Props(Fluid,'Tcrit')
rhoc = Props(Fluid,'Tcrit')
T = np.linspace(Tt+0.01,0.99*Tc,20)
h0 = Props('H','T',0.95*Tc,'Q',1,Fluid)
s0 = Props('S','T',0.95*Tc,'Q',1,Fluid)
h0_RP = Props('H','T',0.95*Tc,'Q',1,RPFluid)
s0_RP = Props('S','T',0.95*Tc,'Q',1,RPFluid)
from parse_bib import BibTeXerClass
BTC = BibTeXerClass()
EOSkey = CP.get_BibTeXKey(Fluid, "EOS")
CP0key = CP.get_BibTeXKey(Fluid, "CP0")
SURFACE_TENSIONkey = CP.get_BibTeXKey(Fluid, "SURFACE_TENSION")
VISCOSITYkey = CP.get_BibTeXKey(Fluid, "VISCOSITY")
CONDUCTIVITYkey = CP.get_BibTeXKey(Fluid, "CONDUCTIVITY")
ECS_LENNARD_JONESkey = CP.get_BibTeXKey(Fluid, "ECS_LENNARD_JONES")
ECS_FITSkey = CP.get_BibTeXKey(Fluid, "ECS_FITS")
BibInfo = ''
if EOSkey:
BibInfo += '<p><b>Equation of State</b>: ' + BTC.entry2HTML(EOSkey)
if CP0key:
BibInfo += '<p><b>Ideal-Gas Specific Heat</b>: ' + BTC.entry2HTML(CP0key)
if SURFACE_TENSIONkey:
BibInfo += '<p><b>Surface Tension</b>: ' + BTC.entry2HTML(SURFACE_TENSIONkey)
if VISCOSITYkey:
BibInfo += '<p><b>Viscosity</b>: ' + BTC.entry2HTML(VISCOSITYkey)
if CONDUCTIVITYkey:
BibInfo += '<p><b>Conductivity</b>: ' + BTC.entry2HTML(CONDUCTIVITYkey)
if ECS_LENNARD_JONESkey:
BibInfo += '<p><b>Lennard-Jones Parameters for ECS</b>: ' + BTC.entry2HTML(ECS_LENNARD_JONESkey)
if ECS_FITSkey:
BibInfo += '<p><b>ECS Correction Fit</b>: ' + BTC.entry2HTML(ECS_FITSkey)
display_html(HTML(BibInfo))
aliases = get_fluid_param_string(Fluid,'aliases')
if aliases:
display_html(HTML(', '.join(['"'+alias+'"' for alias in aliases.split(', ')])))
else:
print '"'+Fluid+'"'
params = dict(mm = CP.Props(Fluid,'molemass'),
Tt = CP.Props(Fluid,'Ttriple'),
pt = CP.Props(Fluid,'ptriple'),
Tc = CP.Props(Fluid,'Tcrit'),
pc = CP.Props(Fluid,'pcrit'),
rhoc = CP.Props(Fluid,'rhocrit'),
Tmin = CP.Props(Fluid,'Tmin'),
CAS = get_fluid_param_string(Fluid,'CAS'),
ASHRAE = get_fluid_param_string(Fluid,'ASHRAE34')
)
s = """<table border = "1">
<tr> <th>Parameter</th> <th>Value</th> </tr>
<tr > <td colspan="2"><center>Triple point</center></td> </tr>
<tr> <td>Triple Point Temp. [K]</td> <td>{Tt:0.3f}</td> </tr>
<tr> <td>Triple Point Press. [kPa]</td> <td>{pt:0.10g}</td> </tr>
<tr > <td colspan="2" ><center>Critical point</center></td> </tr>
<tr> <td>Critical Point Temp. [K]</td> <td>{Tc:0.3f}</td> </tr>
<tr> <td>Critical Point Press. [kPa]</td> <td>{pc:0.10g}</td> </tr>
<tr> <td>Critical Point Density. [kPa]</td> <td>{rhoc:0.10g}</td> </tr>
<tr> <td colspan="2"><center>Other Values</center></td> </tr>
<tr> <td>Mole Mass [kg/kmol]</td> <td>{mm:0.5f}</td> </tr>
<tr> <td>Minimum temperature [K]</td> <td>{Tmin:0.3f}</td> </tr>
<tr> <td>CAS number</td> <td>{CAS:s}</td> </tr>
<tr> <td>ASHRAE classification</td> <td>{ASHRAE:s}</td> </tr>
</table>""".format(**params)
HTML(s)
fig = plt.figure(figsize=(13,6))
ax1 = fig.add_subplot(121)
ax2 = fig.add_subplot(122)
for ax, Q in zip([ax1,ax2],[0,1]):
for key, label in sat_labels.iteritems():
CPval = Props(key, 'T', T, 'Q', Q, Fluid)
RPval = Props(key, 'T', T, 'Q', Q, RPFluid)
if key =='H':
CPval -= h0
RPval -= h0_RP
if key =='S':
CPval -= s0
RPval -= s0_RP
ax.semilogy(T/Tc, np.abs(CPval/RPval-1)*100, 'o', label=label)
if Q == 0:
ax.set_title('Saturated Liquid')
else:
ax.set_title('Saturated Vapor')
ax.set_ylim(1e-18,1000)
ax.set_xlabel('Reduced temperature T/Tc')
ax.set_ylabel('Absolute deviation [%]')
ax.legend(numpoints=1,loc='lower center',ncol = 5)
rho = np.linspace(1e-12,1.5*rhoc)
fig = plt.figure(figsize=(6,6))
ax = fig.add_subplot(111)
for key, label in labels.iteritems():
CPval = Props(key, 'T', Tc, 'D', rho, Fluid)
RPval = Props(key, 'T', Tc, 'D', rho, RPFluid)
if key =='H':
CPval -= h0
RPval -= h0_RP
if key =='S':
CPval -= s0
RPval -= s0_RP
ax.semilogy(rho/rhoc, np.abs(CPval/RPval-1)*100, 'o', label=label)
ax.set_ylim(1e-18,100)
ax.set_title('Critical isotherm Deviations from REFPROP 9.1')
ax.set_xlabel(r'Reduced density $\rho/\rho_c$')
ax.set_ylabel('Absolute deviation [%]')
ax.legend(numpoints=1,loc='lower center',ncol=4)
plt.show()
fig = plt.figure(figsize=(10,5))
ax1 = fig.add_subplot(121)
ax2 = fig.add_subplot(122)
Tmin = Props(Fluid,'Tmin')+3
pmin = Props('P','T',Tmin,'Q',0,Fluid)
pmax = Props(Fluid,'pcrit')*2
hmin = Props('H','T',Tmin,'Q',0,Fluid)
hmax = 2*Props('H','T',Tc-1,'Q',1,Fluid)-hmin
smin = Props('S','T',Tmin,'Q',0,Fluid)
smax = 2*Props('S','T',Tc-1,'Q',1,Fluid)-smin
Ph(Fluid, axis = ax1, Tmin = Tmin, Tmax = Tc-0.01)
Ps(Fluid, axis = ax2, Tmin = Tmin, Tmax = Tc-0.01)
for p in np.linspace(pmin,pmax,10):
for h in np.linspace(hmin,hmax):
_bad = False
try:
T = Props('T','H',h,'P',p,Fluid)
rho = Props('D','H',h,'P',p,Fluid)
hEOS = Props('H','T',T,'D',rho,Fluid)
except ValueError:
_bad = True
if _bad or abs(hEOS/h-1)>1e-6:
ax1.plot(h,p,'x',ms = 10)
else:
ax1.plot(h,p,'k.', ms = 1)
for p in np.linspace(pmin,pmax,10):
for s in np.linspace(smin,smax):
_bad = False
try:
T = Props('T','S',s,'P',p,Fluid)
rho = Props('D','S',s,'P',p,Fluid)
sEOS = Props('S','T',T,'D',rho,Fluid)
except ValueError:
_bad = True
if _bad or abs(sEOS/s-1)>1e-6:
ax2.plot(s,p,'x',ms = 10)
else:
ax2.plot(s,p,'k.', ms = 1)
plt.tight_layout()
