Chemical#
- class Chemical(ID, cache=None, *, search_ID=None, eos=None, phase_ref=None, CAS=None, default=False, phase=None, search_db=True, V=None, Cn=None, mu=None, Cp=None, rho=None, sigma=None, kappa=None, epsilon=None, Psat=None, Hvap=None, method=None, **data)[source]#
Creates a Chemical object which contains constant chemical properties, as well as thermodynamic and transport properties as a function of temperature and pressure.
- Parameters:
ID (str) –
- One of the following [-]:
Name, in IUPAC form or common form or an alias registered in PubChem
InChI name, prefixed by ‘InChI=1S/’ or ‘InChI=1/’
InChI key, prefixed by ‘InChIKey=’
PubChem CID, prefixed by ‘PubChem=’
SMILES (prefix with ‘SMILES=’ to ensure smiles parsing)
CAS number
cache (optional) – Whether or not to use cached chemicals and cache new chemicals.
search_ID (str, optional) – ID to search through database. Pass this key-word argument when you’d like to give a custom name to the chemical, but cannot find the chemical with that name.
eos (GCEOS subclass, optional) – Equation of state class for solving thermodynamic properties. Defaults to Peng-Robinson.
phase_ref ({'s', 'l', 'g'}, optional) – Reference phase of chemical.
CAS (str, optional) – CAS number of chemical.
phase ({'s', 'l' or 'g'}, optional) – Phase to set state of chemical.
search_db=True (bool, optional) – Whether to search the data base for chemical.
Cn (float or function(T), optional) – Molar heat capacity model [J/mol] as a function of temperature [K].
sigma (float or function(T), optional) – Surface tension model [N/m] as a function of temperature [K].
epsilon (float or function(T), optional) – Relative permittivity model [-] as a function of temperature [K].
Psat (float or function(T), optional) – Vapor pressure model [N/m] as a function of temperature [K].
Hvap (float or function(T), optional) – Heat of vaporization model [J/mol] as a function of temperature [K].
V (float or function(T, P), optional) – Molar volume model [m3/mol] as a function of temperature [K] and pressure [Pa].
mu (float or function(T, P), optional) – Dynamic viscosity model [Pa*s] as a function of temperature [K] and pressure [Pa].
kappa (float or function(T, P), optional) – Thermal conductivity model [W/m/K] as a function of temperature [K] and pressure [Pa].
Cp (float, optional) – Constant heat capacity model [J/g].
rho (float, optional) – Constant density model [kg/m3].
default=False (bool, optional) – Whether to default any missing chemical properties such as molar volume, heat capacity, surface tension, thermal conductivity, and molecular weight to that of water (on a weight basis).
Examples
Chemical objects contain the same thermodynamic models and data as thermo’s Chemical objects, but present an API more suitable for BioSTEAM’s needs. Chemical objects contain pure component properties:
>>> import thermosteam as tmo >>> # Initialize with an identifier >>> # (e.g. by name, CAS, InChI...) >>> Water = tmo.Chemical('Water') >>> # Water.show() >>> # Chemical: Water (phase_ref='l') >>> # [Names] CAS: 7732-18-5 >>> # InChI: H2O/h1H2 >>> # InChI_key: XLYOFNOQVPJJNP-U... >>> # common_name: water >>> # iupac_name: ('oxidane',) >>> # pubchemid: 962 >>> # smiles: O >>> # formula: H2O >>> # [Groups] Dortmund: <1H2O> >>> # UNIFAC: <1H2O> >>> # PSRK: <1H2O> >>> # NIST: <Empty> >>> # [Data] MW: 18.015 g/mol >>> # Tm: 273.15 K >>> # Tb: 373.12 K >>> # Tt: 273.16 K >>> # Tc: 647.1 K >>> # Pt: 611.65 Pa >>> # Pc: 2.2064e+07 Pa >>> # Vc: 5.5948e-05 m^3/mol >>> # Hf: -2.8582e+05 J/mol >>> # S0: 70 J/K/mol >>> # LHV: -44011 J/mol >>> # HHV: -0 J/mol >>> # Hfus: 6010 J/mol >>> # Sfus: None >>> # omega: 0.3443 >>> # dipole: 1.85 Debye >>> # similarity_variable: 0.16653 >>> # iscyclic_aliphatic: 0 >>> # combustion: {'H2O': 1.0}
All fields shown are accessible:
>>> Water.CAS '7732-18-5'
Functional group identifiers (e.g. Dortmund, UNIFAC, PSRK) allow for the estimation of activity coefficients through group contribution methods. In other words, these attributes define the functional groups for thermodynamic equilibrium calculations:
>>> Water.Dortmund <DortmundGroupCounts: 1H2O>
Temperature (in Kelvin) and pressure (in Pascal) dependent properties can be computed:
>>> # Vapor pressure (Pa) >>> Water.Psat(T=373.15) 101417.99 >>> # Surface tension (N/m) >>> Water.sigma(T=298.15) 0.0719 >>> # Molar volume (m^3/mol) >>> Water.V(phase='l', T=298.15, P=101325) 1.806...e-05
Note that the reference state of all chemicals is 25 degC and 1 atm:
>>> (Water.T_ref, Water.P_ref) (298.15, 101325.0) >>> # Enthalpy at reference conditions (J/mol; without excess energies) >>> Water.H(T=298.15, phase='l') 0.0
Constant pure component properties are also available:
>>> # Molecular weight (g/mol) >>> Water.MW 18.01528 >>> # Boiling point (K) >>> Water.Tb 373.124
Temperature dependent properties are managed by objects:
>>> # Water.Psat >>> # VaporPressure(CASRN="7732-18-5", Tb=373.124, Tc=647.14, Pc=22048320.0, omega=0.344, extrapolation="AntoineAB|DIPPR101_ABC", method="HEOS_FIT")
Phase dependent properties have attributes with model handles for each phase:
>>> Water.V <PhaseTPHandle(phase, T, P=None) -> V [m^3/mol]> >>> # Water.V.l >>> # VolumeLiquid(CASRN="7732-18-5", MW=18.01528, Tb=373.124, Tc=647.14, Pc=22048320.0, Vc=5.6-05, Zc=0.2294727, omega=0.344, dipole=1.85, Psat=VaporPressure(CASRN="7732-18-5", Tb=373.124, Tc=647.14, Pc=22048320.0, omega=0.344, extrapolation="AntoineAB|DIPPR101_ABC", method="WAGNER_MCGARRY"), eos=[PR(Tc=647.14, Pc=22048320.0, omega=0.344, T=298.15, P=101325.0)], extrapolation="constant", method="VDI_PPDS", method_P="COSTALD_COMPRESSED", tabular_extrapolation_permitted=True)
A new model can be added easily using add_method, for example:
>>> def User_antoine_model(T): ... return 10.0**(10.116 - 1687.537 / (T - 42.98)) >>> # Water.Psat.add_method(f=User_antoine_model, Tmin=273.20, Tmax=473.20) >>> # Water.Psat >>> # VaporPressure(CASRN="7732-18-5", Tb=373.124, Tc=647.14, Pc=22048320.0, omega=0.344, extrapolation="AntoineAB|DIPPR101_ABC", method="USER_METHOD")
The add_method method is a high level interface that even lets you create a constant model:
>>> # Water.Cn.l.add_method(75.31) >>> # Water.Cn('l', 350) 75.31
Note
Because no bounds were given, the model assumes it is valid across all temperatures and pressures.
Choose what model to use through the method attribute:
>>> # list(sorted(Water.Cn.l.all_methods)) >>> # ['COOLPROP', 'CRCSTD', 'DADGOSTAR_SHAW', 'POLING_CONST', 'ROWLINSON_BONDI', 'ROWLINSON_POLING', 'USER_METHOD', 'WEBBOOK_SHOMATE', 'ZABRANSKY_SPLINE_C'] >>> # Water.Cn.l.method = 'ZABRANSKY_SPLINE_C'
Note
For details on the available methods and the capabilities of objects like VaporPressure and VolumeLiquid, visit thermo’s documentation.
- mu(phase, T, P)#
Dynamic viscosity [Pa*s].
- kappa(phase, T, P)#
Thermal conductivity [W/m/K].
- V(phase, T, P)#
Molar volume [m^3/mol].
- Cn(phase, T)#
Molar heat capacity [J/mol/K].
- Psat(T)#
Vapor pressure [Pa].
- Hvap(T)#
Heat of vaporization [J/mol]
- sigma(T)#
Surface tension [N/m].
- epsilon(T)#
Relative permittivity [-]
- S(phase, T, P)#
Entropy [J/mol].
- H(phase, T)#
Enthalpy [J/mol].
- S_excess(T, P)#
Excess entropy [J/mol].
- H_excess(T, P)#
Excess enthalpy [J/mol].
- T_ref = 298.15#
[float] Reference temperature in Kelvin
- P_ref = 101325.0#
[float] Reference pressure in Pascal
- H_ref = 0.0#
[float] Reference enthalpy in J/mol
- chemical_cache = {}#
dict[str, Chemical] Cached chemicals
- cache = False#
[bool] Wheather or not to search cache by default
- classmethod new(ID, CAS, eos=None, phase_ref=None, phase=None, **data)[source]#
Create a new chemical from data without searching through the database, and load all possible models from given data.
- classmethod blank(ID, CAS=None, phase_ref=None, phase=None, formula=None, aliases=None, synonyms=None, free_energies=True, atoms=None, **data)[source]#
Return a new Chemical object without any thermodynamic models or data (unless provided).
- Parameters:
Examples
>>> from thermosteam import Chemical >>> Substance = Chemical.blank('Substance') >>> Substance.show() Chemical: Substance (phase_ref='l') [Names] CAS: Substance InChI: None InChI_key: None common_name: None iupac_name: None pubchemid: None smiles: None formula: None [Groups] Dortmund: <Empty> UNIFAC: <Empty> PSRK: <Empty> NIST: <Empty> [Data] MW: None Tm: None Tb: None Tt: None Tc: None Pt: None Pc: None Vc: None Hf: None S0: 0 J/K/mol LHV: None HHV: None Hfus: None Sfus: None omega: None dipole: None similarity_variable: None iscyclic_aliphatic: None combustion: None
- property charge#
Charge of chemical as described in the chemical formula.
- copy(ID, CAS=None, **data)[source]#
Return a copy of the chemical with a new ID.
Examples
>>> import thermosteam as tmo >>> Glucose = tmo.Chemical('Glucose') >>> Mannose = Glucose.copy('Mannose') >>> assert Mannose.MW == Glucose.MW
- rho(*args, **kwargs)[source]#
Return density given thermal condition [kg/m^3].
Examples
>>> import thermosteam as tmo >>> Water = tmo.Chemical('Water', cache=True) >>> Water.rho('l', 298.15, 101325) 997.064
- Cp(*args, **kwargs)[source]#
Return heat capacity given thermal condition [J/g/K].
Examples
>>> import thermosteam as tmo >>> Water = tmo.Chemical('Water', cache=True) >>> Water.Cp('l', 298.15, 101325) 4.181
- alpha(*args, **kwargs)[source]#
Return thermal diffusivity given thermal condition [m^2/s].
Examples
>>> import thermosteam as tmo >>> Water = tmo.Chemical('Water', cache=True) >>> Water.alpha('l', 298.15, 101325) 1.425...-07
- nu(*args, **kwargs)[source]#
Return kinematic viscosity given thermal condition [m^2/s].
Examples
>>> import thermosteam as tmo >>> Water = tmo.Chemical('Water', cache=True) >>> Water.nu('l', 298.15, 101325) 8.928...e-07
- Pr(*args, **kwargs)[source]#
Return the Prandtl number given thermal condition [-].
Examples
>>> import thermosteam as tmo >>> Water = tmo.Chemical('Water', cache=True) >>> Water.Pr('l', 298.15, 101325) 6.261
- get_property(name, units, *args, **kwargs)[source]#
Return property in requested units.
- Parameters:
Examples
>>> import thermosteam as tmo >>> Water = tmo.Chemical('Water', cache=True) >>> Water.get_property('sigma', 'N/m', 300.) # Surface tension 0.07168
>>> Water.get_property('rho', 'g/cm3', 'l', 300., 101325) # Density 0.9965
- property phase_ref#
{‘s’, ‘l’, ‘g’} Phase at 298 K and 101325 Pa.
- property ID#
[str] Identification of chemical.
- property CAS#
[str] CAS number of chemical.
- property aliases#
set[str] User-defined aliases.
- property synonyms#
set[str] User-defined aliases.
- property InChI#
[str] IUPAC International Chemical Identifier.
- property InChI_key#
[str] IUPAC International Chemical Identifier shorthand.
- property common_name#
[str] Common name identifier.
- property iupac_name#
[tuple[str]] Standard names as defined by IUPAC.
- property pubchemid#
[str] Chemical identifier as defined by PubMed.
- property smiles#
[str] Chemical SMILES formula.
- property formula#
[str] Chemical atomic formula.
- property Dortmund#
[DortmundGroupCounts] Dictionary-like object with functional group numerical identifiers as keys and the number of groups as values.
- property UNIFAC#
[UNIFACGroupCounts] Dictionary-like object with functional group numerical identifiers as keys and the number of groups as values.
- property PSRK#
[PSRKGroupCounts] Dictionary-like object with functional group numerical identifiers as keys and the number of groups as values.
- property NIST#
[NISTGroupCounts] Dictionary-like object with functional group numerical identifiers as keys and the number of groups as values.
- property eos#
[object] Instance for solving equations of state.
- property mu#
Dynamic viscosity [Pa*s].
- property kappa#
Thermal conductivity [W/m/K].
- property V#
Molar volume [m^3/mol].
- property Cn#
Molar heat capacity [J/mol/K].
- property Psat#
Vapor pressure [Pa].
- property Hvap#
Heat of vaporization [J/mol].
- property Hsub#
Heat of sublimation [J/mol].
- property Psub#
Sublimation pressure [Pa].
- property sigma#
Surface tension [N/m].
- property epsilon#
Relative permittivity [-].
- property S#
Entropy [J/mol].
- property H#
Enthalpy [J/mol].
- property S_excess#
Excess entropy [J/mol].
- property H_excess#
Excess enthalpy [J/mol].
- property MW#
Molecular weight [g/mol].
- property Tm#
Normal melting temperature [K].
- property Tb#
Normal boiling temperature [K].
- property Pt#
Triple point pressure [Pa].
- property Tt#
Triple point temperature [K].
- property Tc#
Critical point temperature [K].
- property Pc#
Critical point pressure [Pa].
- property Vc#
Critical point molar volume [m^3/mol].
- property Hfus#
Heat of fusion [J/mol].
- property Sfus#
Entropy of fusion [J/mol].
- property omega#
Acentric factor [-].
- property dipole#
Dipole moment [Debye].
- property similarity_variable#
Similarity variable [-].
- property iscyclic_aliphatic#
Whether the chemical is cyclic-aliphatic.
- property Hf#
Heat of formation at reference phase and temperature [J/mol].
- property LHV#
Lower heating value [J/mol].
- property HHV#
Higher heating value [J/mol].
- property combustion#
dict[str, int] Combustion reaction.
- property Stiel_Polar#
[float] Stiel Polar factor for computing surface tension.
- property Zc#
[float] Compressibility factor.
- property has_hydroxyl#
[bool] Whether or not chemical contains a hydroxyl functional group, as determined by the Dortmund/UNIFAC/PSRK functional groups.
- get_combustion_reaction(chemicals=None, conversion=1.0)[source]#
Return a Reaction object defining the combustion of this chemical. If no combustion stoichiometry is available, return None.
- get_phase(T=298.15, P=101325.0)[source]#
Return phase of chemical at given thermal condition.
Examples
>>> from thermosteam import Chemical >>> Water = Chemical('Water', cache=True) >>> Water.get_phase(T=400, P=101325) 'g'
- Tsat(P, Tguess=None, Tmin=None, Tmax=None, *, check_validity=True)[source]#
Return the saturated temperature (in Kelvin) given the pressure (in Pascal).
- reset(CAS, eos=None, phase_ref=None, smiles=None, InChI=None, InChI_key=None, pubchemid=None, iupac_name=None, common_name=None, formula=None, MW=None, Tm=None, Tb=None, Tc=None, Pc=None, Vc=None, omega=None, Tt=None, Pt=None, Hf=None, S0=None, LHV=None, combustion=None, HHV=None, Hfus=None, dipole=None, similarity_variable=None, iscyclic_aliphatic=None, aliases=None, synonyms=None, *, metadata=None, phase=None, free_energies=True)[source]#
Reset all chemical properties.
- reset_combustion_data(method='Stoichiometry')[source]#
Reset combustion data (LHV, HHV, and combustion attributes) based on the molecular formula and the heat of formation.
- get_key_property_names()[source]#
Return the attribute names of key properties required to model a process.
- default(properties=None)[source]#
Default all properties with the chemical properties of water. If no properties given, all essential chemical properties that are missing are defaulted. properties which are still missing are returned as set.
- Parameters:
properties (Iterable[str], optional) – Names of chemical properties to default.
- Returns:
missing_properties – Names of chemical properties that are still missing.
- Return type:
Examples
>>> from thermosteam import Chemical >>> Substance = Chemical.blank('Substance') >>> missing_properties = Substance.default() >>> sorted(missing_properties) ['Dortmund', 'Hvap', 'PSRK', 'Pc', 'Psat', 'Pt', 'Tb', 'Tc', 'Tm', 'Tt', 'UNIFAC', 'Vc', 'dipole', 'iscyclic_aliphatic', 'omega', 'similarity_variable']
Note that missing properties does not include essential properties volume, heat capacity, and conductivity.
- get_missing_properties(properties=None)[source]#
Return a list all missing thermodynamic properties.
Examples
>>> from thermosteam import Chemical >>> Substance = Chemical.blank('Substance', phase_ref='l') >>> sorted(Substance.get_missing_properties()) ['Cn', 'Dortmund', 'H', 'HHV', 'H_excess', 'Hf', 'Hfus', 'Hvap', 'LHV', 'MW', 'PSRK', 'Pc', 'Psat', 'Pt', 'S', 'S_excess', 'Sfus', 'Tb', 'Tc', 'Tm', 'Tt', 'UNIFAC', 'Vc', 'combustion', 'dipole', 'epsilon', 'iscyclic_aliphatic', 'kappa', 'mu', 'omega', 'sigma', 'similarity_variable']
- property locked_state#
[str] Constant phase of chemical.
- property N_solutes#
[int] Number of molecules formed when solvated.