Substance objects¶
Substance objects are the heart of ReactorD library. Each reagent, product and inert of the reactive system must be individually represented by a Substance object.
For this tutorial, import ReactorD, numpy and matplotlib.pyplot
[1]:
import matplotlib.pyplot as plt
import numpy as np
import reactord as rd
Definition of substance objects¶
A fully defined Substance contains a lot of information. Specific attributes definition will be required for the reactors, described in each reactor documentation. For example, an adiabatic reactor will require that substances define a heat capacity function, on the other hand, when using isothermic reactors this won’t be necessary. Substance has the from_thermo_data_base alternative construction method. Substances objects can be saved as pickle files with the method to_pickle. substances
can be loaded from a pickle file with the method from_pickle.
The minimum necessary information for each substance will be specified on each mixture and each reactor. The complete signature for the Substance objects is:
rd.Substance(
name: str,
molecular_weight: float = None,
normal_boiling_point: float = None,
normal_melting_point: float = None,
critical_temperature: float = None,
critical_pressure: float = None,
acentric_factor: float = None,
formation_enthalpy: float = None,
formation_enthalpy_ig: float = None,
vaporization_enthalpy: Callable = None,
sublimation_enthalpy: Callable = None,
volume_solid: Callable = None,
volume_liquid: Callable = None,
heat_capacity_solid: Callable = None,
heat_capacity_liquid: Callable = None,
heat_capacity_gas: Callable = None,
thermal_conductivity_liquid: Callable = None,
thermal_conductivity_gas: Callable = None,
viscosity_liquid: Callable = None,
viscosity_gas: Callable = None,
heat_capacity_solid_dt_integral: Callable = None,
heat_capacity_liquid_dt_integral: Callable = None,
heat_capacity_gas_dt_integral: Callable = None,
vectorize_functions: bool = False
)
The Temperature and pressure dependent properties must be specified as Python functions of temperature and pressure.
Properties calculated by integration must be specified as Python functions of two temperatures and pressure to return the definite integral between temperature1 and temperature2.
Instantiate a substance¶
For this tutorial, a water Substance object is defined with some example properties.
The data for defining the water is obtained from:
Don W. Green; Robert H. Perry. Perry’s Chemical Engineers’ Handbook, Eighth Edition (McGraw-Hill: New York, Chicago, San Francisco, Lisbon, London, Madrid, Mexico City, Milan, New Delhi, San Juan, Seoul, Singapore, Sydney, Toronto, 2008, 1997, 1984, 1973, 1963, 1950, 1941, 1934).
[2]:
# Defining functions for molar liquid volume and liquid heat capacity.
def water_volume_liquid_tp(temperature, pressure):
"""User defined function for liquid molar volume. [m^3/mol]
Even if pressure is not used in the molar volume function it must be
specified as an argument.
Parameters
----------
temperature : float
Temperature [K]
pressure : float
Pressure [Pa]
Returns
-------
float
Molar volume of the mixture [m^3/mol]
"""
c1 = -13.851
c2 = 0.64038
c3 = -0.00191
c4 = 1.8211e-06
t = temperature
molar_vol = 1 / (c1 + c2 * t + c3 * t**2 + c4 * t**3) / 1000
return molar_vol
def water_heat_capacity_liquid_t(temperature, pressure):
"""User defined function for liquid heat capacity. [J/mol/K]
Parameters
----------
temperature : float
Temperature [K]
Returns
-------
float
Liquid heat capacity [J/mol/K]
"""
c1 = 276370
c2 = -2090.1
c3 = 8.125
c4 = -0.014116
c5 = 9.3701e-06
t = temperature
cp_liquid = (c1 + c2 * t + c3 * t**2 + c4 * t**3 + c5 * t**4) / 1000
return cp_liquid
Now we instantiate a water Substance object in the variable water_byhand. It’s important and mandatory to give a name to your substance. You can use whatever string you want, so use something easy to remember. The name of the substance will be important to define border condition for reactors and to read the simulation results of the reactors.
[3]:
# Defining water Substance by hand
water_byhand = rd.Substance(
name="water",
molecular_weight=18.015,
volume_liquid=water_volume_liquid_tp,
heat_capacity_liquid=water_heat_capacity_liquid_t,
)
# Calling substance methods
print(water_byhand.molecular_weight)
print(water_byhand.volume_liquid(298.15, 101325))
print(water_byhand.heat_capacity_liquid(298.15, 101325))
18.015
1.7999364032438635e-05
75.38420366723898
Broadcasting¶
It’s important that the temperature-pressure functions of the substances follows the broadcasting rules of Numpy. Reactors’ behavior ruled by differential equations will evaluate the substances function with Numpy arrays and will expect an array as a return.
[4]:
temperatures = np.array([300, 400, 500, 600])
pressures = np.array([100_000, 200_000, 300_000, 400_000])
water_byhand.volume_liquid(temperatures, pressures)
[4]:
array([1.80074082e-05, 1.87788490e-05, 1.77064797e-05, 1.31346326e-05])
Defining a constant property¶
For example, if you want to define a constant molar volume you can use the numpy.full() function as:
[5]:
def volume(temperature, pressure):
v = 1.8e-5
return np.full(np.size(temperature), v)
This function will return a constant volume value for all temperatures and pressures.
[6]:
temperatures = np.array([300, 400, 500, 600])
pressures = np.array([100_000, 200_000, 300_000, 400_000])
volume(temperatures, pressures)
[6]:
array([1.8e-05, 1.8e-05, 1.8e-05, 1.8e-05])
Also to assure the broadcasting rules in the functions, you can turn to True the vectorize_functions parameter in the Substance object instantiation. For example:
[7]:
def volume_not_broadcast(temperature, pressure):
return 1.8e-5
The function volume_not_broadcast will not follow the broadcasting rules and will return a unique value for an array evaluation.
[8]:
temperatures = np.array([300, 400, 500, 600])
pressures = np.array([100_000, 200_000, 300_000, 400_000])
volume_not_broadcast(temperatures, pressures)
[8]:
1.8e-05
We can still define a substance with this volume function, but enabling the vectorize_functions parameter.
[9]:
example_substance = rd.Substance(
name="example",
volume_liquid=volume_not_broadcast,
vectorize_functions=True
)
temperatures = np.array([300, 400, 500, 600])
pressures = np.array([100_000, 200_000, 300_000, 400_000])
example_substance.volume_liquid(temperatures, pressures)
[9]:
array([1.8e-05, 1.8e-05, 1.8e-05, 1.8e-05])
The function is vectorized via numpy.vectorize(), leading to a worse calculation performance.
Storing a substance¶
For convinence the object subject create by hand could be save into a file for its posterior use. This file y generate thanks to the `dill <https://github.com/uqfoundation/dill>`__ library.
[10]:
# The object is save
water_byhand.to_pickle("water_file")
# The substance object is called
water_from_pickle = rd.Substance.from_pickle("water_file")
print(water_from_pickle.volume_liquid(298.15, 101325))
1.7999364032438635e-05
Alternative Substance constructor¶
For convenience, an alternative Substance constructor is provided:
substance_object = rd.Substance.from_thermo_database(name, identification)
Where name is a string chosen that will be assigned to the Substance object and identification is a string that may be either the name or the CAS number of the substance. This alternative method uses the Caleb’s Thermo library:
Caleb Bell and Contributors (2016-2021). Thermo: Chemical properties component of Chemical Engineering Design Library (ChEDL) https://github.com/CalebBell/thermo.
[11]:
water_thermo = rd.Substance.from_thermo_database("substance1", "water")
print(water_thermo.molecular_weight)
print(water_thermo.volume_liquid(298.15, 101325))
print(water_thermo.heat_capacity_liquid(298.15, 101325))
18.01528
1.8069338439592966e-05
75.31465144297992
[12]:
# Comparisson
temperature = np.linspace(298.15, 370, 10)
water_vol_byhand = water_byhand.volume_liquid(temperature, 101325)
water_vol_thermo = water_thermo.volume_liquid(temperature, 101325)
plt.plot(temperature, water_vol_byhand)
plt.plot(temperature, water_vol_thermo)
plt.title("Water liquid molar volume ")
plt.legend(["water_byhand", "water_thermo"])
plt.ylabel(r"Molar volume [$m^{3}$/mol]")
plt.xlabel("Temperature [K]")
[12]:
Text(0.5, 0, 'Temperature [K]')
Finally, a list of all the attributes and methods that Substance objects may have is shown below.
IMPORTANT: some data may be not available for a given substance!
[13]:
print(
f"name: {water_thermo.name} \n"
f"molecular weight: {water_thermo.molecular_weight} g/mol \n"
f"normal boiling point: {water_thermo.normal_boiling_point} K \n"
f"normal melting point: {water_thermo.normal_melting_point} K \n"
f"critical temperature: {water_thermo.critical_temperature} K \n"
f"critical pressure: {water_thermo.critical_pressure} Pa \n"
f"acentric factor: {water_thermo.acentric_factor} \n"
f"formation enthalpy: {water_thermo.formation_enthalpy} J/mol \n"
f"formation enthalpy ideal gas: {water_thermo.formation_enthalpy_ig} J/mol"
"\n"
f"formation Gibbs energy ideal gas: {water_thermo.formation_enthalpy_ig}"
" J/mol \n"
f"vaporization enthalpy: {water_thermo.vaporization_enthalpy(373.15)}"
" J/mol \n"
f"sublimation enthalpy: {water_thermo.sublimation_enthalpy(273.15)} J/mol"
"\n"
f"fusion enthalpy: {water_thermo.fusion_enthalpy(298.15)} J/mol \n"
f"volume solid: {water_thermo.volume_solid(273.15, 101325)} m^3/mol \n"
f"volume liquid: {water_thermo.volume_liquid(298.15, 101325)} m^3/mol \n"
f"heat capacity solid: {water_thermo.heat_capacity_solid(373.15, 101325)} "
"J/K/mol \n"
f"heat capacity liquid: {water_thermo.heat_capacity_liquid(373, 101325)}"
" J/K/mol \n"
f"heat capacity gas: {water_thermo.heat_capacity_gas(373.15, 101325)} "
"J/K/mol \n"
"thermal conductivity liquid: "
f"{water_thermo.thermal_conductivity_liquid(298.15, 101325)}"
" W/m/K \n"
"thermal conductivity gas: "
f"{water_thermo.thermal_conductivity_gas(373.15, 101325)} W/m/K \n"
f"viscosity liquid: {water_thermo.viscosity_liquid(298.15, 101325)} Pa/s"
"\n"
f"viscosity gas: {water_thermo.viscosity_gas(400, 101325)} Pa/s \n"
"cp integral solid: "
f"{water_thermo.heat_capacity_solid_dt_integral(250, 280, 101325)} "
"J/mol \n"
"cp integral liquid liq: "
f"{water_thermo.heat_capacity_liquid_dt_integral(300, 330, 101325)} "
"J/mol \n"
"cp integral gas: "
f"{water_thermo.heat_capacity_gas_dt_integral(380, 410, 101325)} "
"J/mol"
)
name: substance1
molecular weight: 18.01528 g/mol
normal boiling point: 373.124 K
normal melting point: 273.15 K
critical temperature: 647.14 K
critical pressure: 22048320.0 Pa
acentric factor: 0.344
formation enthalpy: -285825.0 J/mol
formation enthalpy ideal gas: -241822.0 J/mol
formation Gibbs energy ideal gas: -241822.0 J/mol
vaporization enthalpy: 40798.29512500727 J/mol
sublimation enthalpy: 50742.6591346729 J/mol
fusion enthalpy: 7127.725817752558 J/mol
volume solid: 1.5989414456471014e-05 m^3/mol
volume liquid: 1.8069338439592966e-05 m^3/mol
heat capacity solid: 30.969944935783165 J/K/mol
heat capacity liquid: 76.05035935651769 J/K/mol
heat capacity gas: 34.03105522722508 J/K/mol
thermal conductivity liquid: 0.606299909107096 W/m/K
thermal conductivity gas: 0.023865413505450152 W/m/K
viscosity liquid: 0.0009125307951858123 Pa/s
viscosity gas: 1.3587844730035658e-05 Pa/s
cp integral solid: 681.6705747748953 J/mol
cp integral liquid liq: 2258.000649097829 J/mol
cp integral gas: 1026.395658287748 J/mol