matplotlib

Previous topic

cis-2-Butene

Next topic

n-Butane

This Page

m-XyleneΒΆ

View this page as an IPython notebook

In [1]:
# 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 = 'm-Xylene'
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)

Literature References

In [2]:
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))

Equation of State: Zhou, Yong; Wu, Jiangtao; Lemmon, Eric W., 2012, Thermodynamic Properties of o-Xylene, m-Xylene, p-Xylene, and Ethylbenzene, J. Phys. Chem. Ref. Data, 41:023103-1 -- 023103-26

Aliases

In [3]:
aliases = get_fluid_param_string(Fluid,'aliases')
if aliases:
    display_html(HTML(', '.join(['"'+alias+'"' for alias in aliases.split(', ')])))
else:
    print '"'+Fluid+'"'
"mXylene", "m-xylene", "M-XYLENE"

Fluid Constants

In [4]:
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)
Out[4]:
Parameter Value
Triple point
Triple Point Temp. [K] 225.300
Triple Point Press. [kPa] 0.003123267493
Critical point
Critical Point Temp. [K] 616.890
Critical Point Press. [kPa] 3534.6
Critical Point Density. [kPa] 282.929725
Other Values
Mole Mass [kg/kmol] 106.16500
Minimum temperature [K] 225.300
CAS number 108-38-3
ASHRAE classification UNKNOWN

Saturation Property Deviations

In [5]:
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)

Along the critical isotherm where \(T=T_c\)

In [6]:
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()
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-6-719d680df4f4> in <module>()
      6 for key, label in labels.iteritems():
      7     CPval = Props(key, 'T', Tc, 'D', rho, Fluid)
----> 8     RPval = Props(key, 'T', Tc, 'D', rho, RPFluid)
      9     if key =='H':
     10         CPval -= h0

C:\Python27\lib\site-packages\CoolProp\CoolProp.pyd in CoolProp.CoolProp.Props (build\temp.win32-2.7\Release\pyrex\CoolProp\CoolProp.cpp:21147)()

C:\Python27\lib\site-packages\CoolProp\CoolProp.pyd in CoolProp.CoolProp.Props (build\temp.win32-2.7\Release\pyrex\CoolProp\CoolProp.cpp:19968)()

ValueError: CoolProp error: [TDFLSH error  4] one or more inputs are out of range:    pressure > 2 x upper limit, P =  405.57 MPa, Pmax =  200.00 MPa  :: inputs were:"P",'T',6.1688999999999999e+02,'D',8.8756622448979601e+02,"REFPROP-MXYLENE"

Check of p,h and p,s as inputs (X: Failure .: Success)

In [7]:
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()