[ANN] PyThermo.jl, an interface to the Thermo library for thermophysical properties

PyThermo may not be the fastest thermophysical property library in the land, but thanks to the tireless efforts of Caleb Bell (the author of the underlying Python library, Thermo), it’s quite thorough.

Data is available for a vast number of chemical species, with Unitful accessor functions for the most common properties. Additional properties can be queried via the underlying PyObject, and should be assumed to have units of J, kg, m, s, K unless otherwise specified in Thermo’s documentation.

Water:

julia> H2O = Species("water", T = 40u"°C")
Species(water, 313.1 K, 1.013e+05 Pa)

julia> density(H2O)
997.0326202881853 kg m^-3

julia> H2O.calculate(T = 1000, P = 3000) # K, Pa

julia> density(H2O)
0.006500222660286359 kg m^-3

Air:

julia> air = Mixture(["N2" => 0.78, "O2" => 0.21, "Ar" => 0.01], P = 1u"atm")
Mixture({N2: 0.78, O2: 0.21, Ar: 0.01}, 298.1 K, 1.013e+05 Pa)

julia> R_specific(air)
287.0055796266994 J kg^-1 K^-1

Vodka:

julia> vodka = Mixture(["water", "ethanol"], Vfls = [.6, .4]) # liquid volume frac
Mixture({water: 0.83, ethanol: 0.17}, 298.1 K, 1.013e+05 Pa)

julia> vodka.Pc # critical pressure
1.933817214110451e7

1,1,1,3,5,5,5-heptamethyltrisiloxane:

julia> siloxane = Species("1,1,1,3,5,5,5-heptamethyltrisiloxane", T=500u"°F")
Species(1,1,1,3,5,5,5-heptamethyltrisiloxane, 533.1 K, 1.013e+05 Pa)

julia> soundspeed(siloxane)
142.92248228212017 m s^-1

julia> siloxane.Cp
1853.573089850276
10 Likes

Nice. How does Thermo compare to CoolProp?

As a big hairball of a Python library, it’s a few times slower than CoolProp.jl (which is just a wrapper around the CoolProp C++ library), but includes more properties & many more species (see here for CoolProp’s, or here for a subset of Thermo’s).

using PyThermo, CoolProp

function ρ_TP(species::Species, T, P)
    species.calculate(T = T, P = P)
    density(species)
end

function ρ_TP(cp_handle, T, P)
    CoolProp.AbstractState_update(cp_handle,CoolProp.get_input_pair_index("PT_INPUTS"), P, T)
    CoolProp.AbstractState_keyed_output(cp_handle,CoolProp.get_param_index("D"))
end
julia> thermo_H2O = Species("water");

julia> CP_H2O = CoolProp.AbstractState_factory("HEOS", "Water");

julia> @time [ρ_TP(thermo_H2O, T, P) for T in LinRange(280, 330, 100), P in LinRange(8e4, 12e5, 100)];
  3.598842 seconds (448.88 k allocations: 22.063 MiB, 0.38% gc time, 1.86% compilation time)

julia> @time [ρ_TP(CP_H2O, T, P) for T in LinRange(280, 330, 100), P in LinRange(8e4, 12e5, 100)];
  0.538447 seconds (115.41 k allocations: 4.366 MiB, 6.34% compilation time)

Thermo also allows the CoolProp Python binding to be used for certain properties:

julia> thermo_H2O.Psat
17217.65709618894

julia> thermo_H2O.VaporPressure.set_user_methods("COOLPROP")

julia> thermo_H2O.Psat
17213.15165587802

Thanks. I will have to take a deeper look. With all the Symbolics discussion going on, it is beginning to feel that a pure julia solution would be nice.

Paulo

2 Likes

Yes, a similarly-featured Julia package would be great. There are a number of existing thermo-oriented Julia packages, but they’re all kinda pointing in different directions, with very different takes on API.

hi, im the main developer of ThermoModels.jl. basically im porting all the relevant functionality to OpenSAFT.jl (on Juliacon 2021, chat in Zulip). What i learned abut implementing EoS in julia is the following:

  • AD is a game changer, thermodynamic models are basically combinations of partial derivatives, being able to automatically differenciate though models allows to use wild methods (for example, using a virial aproximation as a first guess for gas volume).
  • its necesary to stablish a standard on how to specify a thermodynamic state. how do we say to a model, i want to calculate this with fixed enthalpy and temperature, or fixed temperature and on a saturation line. one of my attempts to do this is in ThermoState.jl using julia dispatch to select the correct specification. (i used this on MoistAir.jl and WaterIF97.jl ,as a proof of concept on how such system could be implemented). decoupling the state from the model is important, as allows to calculate multiple states in parallel. For example, in a rachford rice iteration, the liquid and volume root finders could be executed at the same time, or in global optimization methods for multiphase equilibria, we could just throw a bunch of processes on a cluster to calculate everithing
  • with general derivatives, the only things necessary to define a EoS is the base function (namely, gibbs or helmholtz energy) and a minimun, lower bound volume. this lower bound volume is used to calculate liquid volumes. this is specially important when defining SAFT equations, as we dont have the rich amount of correlations that critical models have, the idea is that the function and the lower bound volume should be enough for defining models, not caring about starting points for iterations, as those could in theory be derived by the information stored in the function and its behaviour.
  • data is necessary, thats where the Thermo python package hits the jackpot, it doesnt matter if you have the most advanced equations, if you don’t have data to use those equations.

On a related theme, have you seen GitHub - CalebBell/chemicals: chemicals: Chemical database of Chemical Engineering Design Library (ChEDL) , its basically the database of models + helper functions, a good approach could be to use Thermo or Chemicals to query the results and julia to calculate properties. for example, this issue numba_vectorized reports TypeError · Issue #28 · CalebBell/chemicals · GitHub is basically solved by julia broadcasting

1 Like

CoolProp seems to use the same general approach when calculating derivatives, Pure and Pseudo-Pure fluid properties — CoolProp 6.4.1 documentation, so in theory, we can just plug in the helmholtz model defined on coolprop and make it work on julia.

Yeah, I don’t think it makes much sense to try to duplicate Caleb’s massive property/coefficient databases in Julia - we can stick to nailing the computationally-intensive stuff and implementing clean, Julian interfaces. ThermoState seems like it could be a solid foundation for thermo efforts in Julia - have you done any benchmarks versus CoolProp?

1 Like

ThermoState.does not compete with CoolProp, instead proposes a common property interface api available to exchange between models. A dificulty on comparing packages is that we need the same equation of state to compare.a PC SAFT requires that requires an internal fixed point calculation on each property call, is about 100 times slower than a Peng Robinson eos, for example.

What could be done with a ThermoState object is to pass it between different models, for example, for speed, maybe a cubic suffices when simulating a compressor, but you need a more complete EoS when simulating a distillation tower. you just pass the corresponding state and the eos uses the appropiate equations to obtain all the needed properties.
the apropiate benchmark package is OpenSAFT, as it has component data and eos models. both packages have an IAPWS95 implementation, we can compare those as a starting point. from CoolProp docs, it seems that IAPWS95 is the standard when you calculate with water.

Pkg.add("OpenSAFT#database") #latest branch
using OpenSAFT
water = IAPWS95()
OpenSAFT.volume(water,101325,300)

on CoolProp (using CoolProp.jl)

Pkg.clone("https://github.com/vimalaad/CoolProp.jl.git")
Pkg.build("CoolProp") # to download the latest binaries
using CoolProp
PropsSI("D", "T", 300, "P", 101325, "Water")