Custom Fluid Models (props="CUSTOM")

carbatpy normally uses REFPROP or TREND as the thermodynamic backend. With props="CUSTOM" you can plug in any Python callable as the property model — no external library required.

Typical use cases:

• Simplified analytic models (constant cp, ideal gas, incompressible liquid)

• Temperature- or pressure-dependent polynomials fitted to data

• Interpolators built from measurement data

• Testing or prototyping a cycle without a full EOS installation

Once defined, a custom fluid behaves exactly like any other carbatpy Fluid: set_state, print_state, val_dict, properties.enthalpy, set_state_v — all work unchanged.

The interface

You supply a callable with the signature

def my_props(values, given, fluid, composition, wanted):
    ...
    return [T, P, H, V, S, Q, U]   # SI, mass basis

Argument	Type	Description
values	list[float, float]	The two input values fixing the thermodynamic state
given	str	Which two properties are given, e.g. "TP", "PH", "TQ"
fluid	str	Fluid string from FluidModel.fluid — forwarded automatically
composition	list	Mole fractions from Fluid.composition — forwarded automatically
wanted	str	The property string requested by the caller (default _THERMO_STRING) — forwarded automatically

If wanted is unchanged and remains at its default _THERMO_STRING, the return value must be a list/array in the order:

Index	Property	Unit
0	Temperature T	K
1	Pressure P	Pa
2	Specific enthalpy H	J/kg
3	Specific volume V	m³/kg
4	Specific entropy S	J/(kg·K)
5	Vapour quality Q	— (use −1 for single-phase)
6	Specific internal energy U	J/kg

Transport properties (viscosity, thermal conductivity, …) may follow optionally in the following order:

Index	Property	Unit
7	Viscosity VIS	Pa*s
8	Thermal conductivity TCX	W/(m*K)
9	Prandtl number PRANDTL
10	Specific heat capacity CP	J/(kg*K)
11	Speed of Sound W	m/s
12	Molecular Mass M	g/mol

If wanted is customized, return the variables as defined.

The callable is passed to init_fluid via props_func=.

import numpy as np
import carbatpy as cb

Example 1 — Analytic model: incompressible liquid with constant cp

A simple approximation for liquid water:

• isobaric heat capacity: \(c_p = 4182\) J/(kg·K)

• density: \(\rho = 997\) kg/m³ (incompressible)

• reference state at \(T_\mathrm{ref} = 273.15\) K

All model parameters live inside the function — the callable is fully self-contained. The function supports both "TP" and "PH" input pairs; for "PH" the temperature is back-calculated from the enthalpy definition.

def simple_water(values, given, fluid, composition, wanted):
    """Incompressible liquid water with constant cp."""
    T_ref = 273.15   # K
    cp    = 4182.0   # J/(kg·K)
    rho   = 997.0    # kg/m³

    if given == "TP":
        T, P = values
    elif given == "PH":
        P, h_in = values
        T = h_in / cp + T_ref
    else:
        raise NotImplementedError(f"Input pair '{given}' not supported by this model.")

    h = cp * (T - T_ref)
    v = 1.0 / rho
    s = cp * np.log(T / T_ref)
    u = h - P * v
    q = -6667.0       # always single-phase
    return [T, P, h, v, s, q, u] #use default order of _THERMO_STRING

Pass it to init_fluid with props="CUSTOM" and props_func=:

water = cb.init_fluid("simple_water", [1.0], props="CUSTOM", props_func=simple_water)

water.set_state([350.0, 2e5], "TP")
water.print_state()

simple_water, composition: [1.0]
Temperature    : 350.0
Pressure    : 200000.0
spec_Enthalpy    : 321386.70000000007
spec_Volume    : 0.0010030090270812437
spec_Entropy    : 1036.7682334217266
quality_mass    : -6667.0
spec_internal_Energy    : 321186.09819458384
x    : None
y    : None
z    : None

The result is accessible as a dict, by named attribute, or as a plain array — exactly like REFPROP/TREND fluids:

water.val_dict

{'Temperature': 350.0,
 'Pressure': 200000.0,
 'spec_Enthalpy': 321386.70000000007,
 'spec_Volume': 0.0010030090270812437,
 'spec_Entropy': np.float64(1036.7682334217266),
 'quality_mass': -6667.0,
 'spec_internal_Energy': 321186.09819458384,
 'x': None,
 'y': None,
 'z': None}

# Named attribute access via FluidState
print(f"Enthalpy : {water.properties.enthalpy:.1f} J/kg")
print(f"Entropy  : {water.properties.entropy:.2f} J/(kg·K)")
print(f"sp_volume: {water.properties.sp_volume:.6f} m³/kg")

Enthalpy : 321386.7 J/kg
Entropy  : 1036.77 J/(kg·K)
sp_volume: 0.001003 m³/kg

# Using a different input pair: PH
h_stored = water.properties.enthalpy
st2 = water.set_state([2e5, h_stored], "PH")
print(f"Round-trip T from PH: {st2[0]:.4f} K  (expected 350.0 K)")

Round-trip T from PH: 350.0000 K  (expected 350.0 K)

Example 2 — Temperature-polynomial model

When properties depend on the thermodynamic state, the independent variable changes depending on the input pair: for "TP" T is known directly; for "PH" it must first be recovered from H. This is handled naturally by branching on given.

Here we use a polynomial \(c_p(T) = a_0 + a_1 T\) (a linear fit). Enthalpy and entropy follow by integration:

\[h(T) = a_0 (T-T_\mathrm{ref}) + \tfrac{a_1}{2}(T^2 - T_\mathrm{ref}^2)\]

\[s(T) = a_0 \ln\!\left(\tfrac{T}{T_\mathrm{ref}}\right) + a_1 (T - T_\mathrm{ref})\]

Again, all coefficients live inside the function.

def poly_fluid(values, given, fluid, composition, wanted):
    """Liquid with a linearly temperature-dependent cp(T) = a0 + a1*T."""
    a0    = 1000.0   # J/(kg·K)
    a1    = 2.5      # J/(kg·K²)
    T_ref = 300.0    # K  (h=0, s=0 reference)
    rho   = 1200.0   # kg/m³ (incompressible)

    if given == "TP":
        T, P = values
    elif given == "PH":
        P, h_in = values
        # Invert h(T) = a0*(T-T_ref) + a1/2*(T²-T_ref²)  →  quadratic in T
        C = -(a0 * T_ref + 0.5 * a1 * T_ref**2 + h_in)
        T = (-a0 + np.sqrt(a0**2 - 4 * 0.5 * a1 * C)) / a1
    else:
        raise NotImplementedError(f"Input pair '{given}' not supported.")

    h = a0 * (T - T_ref) + 0.5 * a1 * (T**2 - T_ref**2)
    s = a0 * np.log(T / T_ref) + a1 * (T - T_ref)
    v = 1.0 / rho
    u = h - P * v
    return [T, P, h, v, s, -6667.0, u]

poly = cb.init_fluid("poly_fluid", [1.0], props="CUSTOM", props_func=poly_fluid)

poly.set_state([360.0, 3e5], "TP")
poly.print_state()

poly_fluid, composition: [1.0]
Temperature    : 360.0
Pressure    : 300000.0
spec_Enthalpy    : 109500.0
spec_Volume    : 0.0008333333333333334
spec_Entropy    : 332.3215567939546
quality_mass    : -6667.0
spec_internal_Energy    : 109250.0
x    : None
y    : None
z    : None

# Verify round-trip consistency: TP → PH → T
h_stored = poly.properties.enthalpy
st_ph = poly.set_state([3e5, h_stored], "PH")
print(f"T from TP : 360.000000 K")
print(f"T from PH : {st_ph[0]:.6f} K")

T from TP : 360.000000 K
T from PH : 360.000000 K

Example 3 — Tabulated / measurement data via a callable class

Interpolators from measurement data are expensive to build and should be constructed once, not on every set_state call. A class with __init__ and __call__ is the natural fit: the interpolators are stored as instance attributes and reused across calls.

carbatpy accepts any callable, so a class instance works just as well as a plain function.

from scipy.interpolate import RegularGridInterpolator

class TabulatedFluid:
    """
    Property look-up via interpolation of a (T, P) property table.

    Parameters
    ----------
    T_tab : 1-D array  (K)
    P_tab : 1-D array  (Pa)
    H_tab : 2-D array  (J/kg),    shape (len(T_tab), len(P_tab))
    V_tab : 2-D array  (m³/kg),   same shape
    S_tab : 2-D array  (J/(kg·K)), same shape
    U_tab : 2-D array  (J/kg),    same shape
    """

    def __init__(self, T_tab, P_tab, H_tab, V_tab, S_tab, U_tab):
        kw = dict(bounds_error=True, method="linear")
        axes = (T_tab, P_tab)
        self._h = RegularGridInterpolator(axes, H_tab, **kw)
        self._v = RegularGridInterpolator(axes, V_tab, **kw)
        self._s = RegularGridInterpolator(axes, S_tab, **kw)
        self._u = RegularGridInterpolator(axes, U_tab, **kw)

    def __call__(self, values, given, fluid, composition, wanted):
        if given != "TP":
            raise NotImplementedError(f"Only 'TP' supported by this interpolator.")
        T, P = values
        pt = np.array([[T, P]])
        return [
            T, P,
            float(self._h(pt)),
            float(self._v(pt)),
            float(self._s(pt)),
            -6667.0,
            float(self._u(pt)),
        ]

# Build synthetic tabulated data — replace with real measurements
T_tab = np.linspace(280.0, 380.0, 20)   # K
P_tab = np.linspace(0.5e5, 5.0e5, 10)  # Pa
T_grid, P_grid = np.meshgrid(T_tab, P_tab, indexing="ij")

T_REF_TAB = 273.15
CP_TAB    = 4000.0
RHO_0_TAB = 1000.0

H_tab = CP_TAB * (T_grid - T_REF_TAB)
V_tab = 1.0 / (RHO_0_TAB - 0.3 * (T_grid - T_REF_TAB))  # slight T-dependence
S_tab = CP_TAB * np.log(T_grid / T_REF_TAB)
U_tab = H_tab - P_grid * V_tab

# Construct the callable — interpolators are built once here
tab_props = TabulatedFluid(T_tab, P_tab, H_tab, V_tab, S_tab, U_tab)

tab_fluid = cb.init_fluid("tabulated", [1.0], props="CUSTOM", props_func=tab_props)

tab_fluid.set_state([320.0, 2e5], "TP")
tab_fluid.print_state()

tabulated, composition: [1.0]
Temperature    : 320.0
Pressure    : 200000.0
spec_Enthalpy    : 187400.0000000001
spec_Volume    : 0.0010142559831652705
spec_Entropy    : 633.0694664672109
quality_mass    : -6667.0
spec_internal_Energy    : 187197.14880336705
x    : None
y    : None
z    : None

C:\Users\welp\AppData\Local\Temp\ipykernel_11608\3240967204.py:32: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.)
  float(self._h(pt)),
C:\Users\welp\AppData\Local\Temp\ipykernel_11608\3240967204.py:33: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.)
  float(self._v(pt)),
C:\Users\welp\AppData\Local\Temp\ipykernel_11608\3240967204.py:34: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.)
  float(self._s(pt)),
C:\Users\welp\AppData\Local\Temp\ipykernel_11608\3240967204.py:36: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.)
  float(self._u(pt)),

Example 4 — CoolProp backend via PropsSI

If CoolProp is installed, it can serve as the backend for a CUSTOM fluid. The callable maps carbatpy’s single-letter input codes (T, P, H, S, Q) to CoolProp’s property names and delegates every call to PropsSI.

Because carbatpy forwards the fluid name automatically via fluid, the same callable works for any CoolProp-supported substance — no model constants need to be hard-coded.

Be careful with mixtures! CoolProp might suggest a result even if it is not valid! See the Coolprop documentation for further information: https://coolprop.sourceforge.net/fluid_properties/Mixtures.html

from CoolProp.CoolProp import PropsSI

# Maps carbatpy/REFPROP property codes → CoolProp names.
# Covers both input-pair codes (single char) and wanted-string codes.
_CP_MAP = {
    'T':       'T',
    'P':       'P',
    'H':       'Hmass',
    'V':       'Dmass',          # density; inverted to specific volume below
    'S':       'Smass',
    'Q':       'Q',              # used in given (input pair)
    'QMASS':   'Q',              # used in wanted string
    'E':       'Umass',
    'VIS':     'viscosity',
    'TCX':     'conductivity',
    'PRANDTL': 'Prandtl',
    'CP':      'Cp0mass',
    'W':       'speed_of_sound',
    'M':       'molar_mass',
}
_CP_INVERT = {'V'}  # specific volume = 1 / density

def coolprop_props(values, given, fluid, composition, wanted):
    """Fluid properties via CoolProp.PropsSI, iterating through *wanted*.

    fluid must be a CoolProp-compatible name, e.g. 'Water', 'R134a', 'Propane'.
    Supported input pairs: any combination of T, P, H, S, Q.
    """
    name1 = _CP_MAP[given[0]]
    name2 = _CP_MAP[given[1]]
    v1, v2 = values

    result = []
    for code in wanted.split(";"):
        cp_name = _CP_MAP[code]
        val = PropsSI(cp_name, name1, v1, name2, v2, fluid)
        if code in _CP_INVERT:
            val = 1.0 / val
        if code == 'QMASS' and not 0.0 <= val <= 1.0:
            val = -6667.0
        result.append(val)
    return result

cp_water = cb.init_fluid("Water", [1.0], props="CUSTOM", props_func=coolprop_props)
cp_water.set_state([350.0, 2e5], "TP")
cp_water.print_state()

Water, composition: [1.0]
Temperature    : 350.0
Pressure    : 200000.0
spec_Enthalpy    : 321918.3569311111
spec_Volume    : 0.0010269339603223562
spec_Entropy    : 1037.888228707492
quality_mass    : -6667.0
spec_internal_Energy    : 321712.97013907303
x    : None
y    : None
z    : None

cp_water.set_state([350.0, 2e5], "TP", wanted=cb.CB_DEFAULTS["Fluid_Defaults"]["TRANS_STRING"])
cp_water.print_state()

Water, composition: [1.0]
Temperature    : 350.0
Pressure    : 200000.0
spec_Enthalpy    : 321918.3569311111
spec_Volume    : 0.0010269339603223562
spec_Entropy    : 1037.888228707492
quality_mass    : -6667.0
spec_internal_Energy    : 321712.97013907303
viscosity    : 0.00036849611088197846
thermal_conductivity    : 0.6649268689536709
Prandtl_number    : 2.32441423502578
isobaric_heat_capacity    : 1880.5502447228766
speed_of_sound    : 1555.1376386761594
molecular_mass    : 0.018015268
x    : None
y    : None
z    : None

cp_r134a = cb.init_fluid("R134a", [1.0], props="CUSTOM", props_func=coolprop_props)
cp_r134a.set_state([300.0, 6e5], "TP")
cp_r134a.print_state()

h_stored = cp_r134a.properties.enthalpy
st_ph = cp_r134a.set_state([6e5, h_stored], "PH")
print(f"\nRound-trip T from PH: {st_ph[0]:.4f} K  (expected 300.0 K)")

R134a, composition: [1.0]
Temperature    : 300.0
Pressure    : 600000.0
spec_Enthalpy    : 415849.56452239887
spec_Volume    : 0.03536749623274023
spec_Entropy    : 1735.2099107648696
quality_mass    : -6667.0
spec_internal_Energy    : 394629.06678279996
x    : None
y    : None
z    : None

Round-trip T from PH: 300.0000 K  (expected 300.0 K)

If you want to use a specific CoolProp-Fluid-Model, you can do so by defining the fluidname.

# Maps carbatpy/REFPROP property codes → CoolProp incompressible-fluid names.
# Note: CoolProp uses different single-letter codes for incompressible fluids:
#   'V' = dynamic viscosity, 'L' = thermal conductivity, 'C' = isobaric cp.
# PRANDTL and W, M require special treatment (see function body).
_CP_MAP_INCOMP = {
    'T':       'T',
    'P':       'P',
    'H':       'Hmass',
    'V':       'Dmass',   # density; inverted to specific volume (reuses _CP_INVERT)
    'S':       'Smass',
    'Q':       'Q',       # for input pairs only; Q as state is not valid for incompressible
    'QMASS':   None,      # always -6667.0 for incompressible fluids
    'E':       'Umass',
    'VIS':     'V',       # CoolProp incompressible: dynamic viscosity
    'TCX':     'L',       # CoolProp incompressible: thermal conductivity
    'CP':      'C',       # CoolProp incompressible: isobaric heat capacity
    'PRANDTL': None,      # computed as VIS * CP / TCX (see below)
    'W':       None,      # not available for incompressible fluids in CoolProp
    'M':       None,      # not available for incompressible fluids in CoolProp
}

_CP_INVERT = {'V'}  # specific volume = 1 / density

def coolprop_incomp_props(values, given, fluid, composition, wanted):
    """Fluid properties for incompressible fluids category via CoolProp.PropsSI.

    Fluid must be of category incompressible fluids (see: https://coolprop.org/fluid_properties/Incompressibles.html).
    fluid must be a CoolProp-compatible name WITHOUT the composition in brackets, e.g. "INCOMP::MEG".
    composition must be defined in correct units specified by CoolProp, see CoolProp documentation.

    Supported input pairs: any combination of T, P, H, S.
    (Q is not allowed as input for incompressible fluids in CoolProp.)
    """
    name1 = _CP_MAP_INCOMP[given[0]]
    name2 = _CP_MAP_INCOMP[given[1]]
    v1, v2 = values

    fluid_extended = fluid + str(composition)

    result = []
    for code in wanted.split(";"):
        if code == 'QMASS':
            result.append(-6667.0)
        elif code == 'PRANDTL':
            mu   = PropsSI('V', name1, v1, name2, v2, fluid_extended)
            cp_v = PropsSI('C', name1, v1, name2, v2, fluid_extended)
            lamb = PropsSI('L', name1, v1, name2, v2, fluid_extended)
            result.append(mu * cp_v / lamb)
        elif code in ('W', 'M'):
            result.append(-1000.)
        else:
            cp_name = _CP_MAP_INCOMP[code]
            val = PropsSI(cp_name, name1, v1, name2, v2, fluid_extended)
            if code in _CP_INVERT:
                val = 1.0 / val
            result.append(val)
    return result

cp_MEG = cb.init_fluid("INCOMP::MEG", [0.5], props="CUSTOM", props_func=coolprop_incomp_props)
cp_MEG.set_state([273.15, 1e5], "TP")
cp_MEG.print_state()

cp_MEG.set_state([273.15, 1e5], "TP", wanted=cb.CB_DEFAULTS["Fluid_Defaults"]["TRANS_STRING"])
cp_MEG.print_state()

INCOMP::MEG, composition: [0.5]
Temperature    : 273.15
Pressure    : 100000.0
spec_Enthalpy    : -65160.29676241186
spec_Volume    : 0.0009305583742099314
spec_Entropy    : -230.17235690828966
quality_mass    : -6667.0
spec_internal_Energy    : -65253.352599832855
x    : None
y    : None
z    : None

INCOMP::MEG, composition: [0.5]
Temperature    : 273.15
Pressure    : 100000.0
spec_Enthalpy    : -65160.29676241186
spec_Volume    : 0.0009305583742099314
spec_Entropy    : -230.17235690828966
quality_mass    : -6667.0
spec_internal_Energy    : -65253.352599832855
viscosity    : 0.007929773683300138
thermal_conductivity    : 0.3768161213114663
Prandtl_number    : 67.4018018553448
isobaric_heat_capacity    : 3202.8764702355397
speed_of_sound    : -1000.0
molecular_mass    : -1000.0
x    : None
y    : None
z    : None

Returning a FluidState object

As with any carbatpy fluid, passing output="FluidState" gives named attribute access:

state_obj = water.set_state([330.0, 1e5], "TP", output="FluidState")

print(f"Type     : {type(state_obj)}")
print(f"Enthalpy : {state_obj.enthalpy:.1f} J/kg")
print(f"Entropy  : {state_obj.entropy:.3f} J/(kg·K)")
print(f"\nFull dict:\n{state_obj.state_to_dict()}")

Type     : <class 'carbatpy.models.cb_fluids.fluid_props.FluidState'>
Enthalpy : 237746.7 J/kg
Entropy  : 790.697 J/(kg·K)

Full dict:
{'Temperature': 330.0, 'Pressure': 100000.0, 'spec_Enthalpy': 237746.7000000001, 'spec_Volume': 0.0010030090270812437, 'spec_Entropy': np.float64(790.6972623258188), 'quality_mass': -6667.0, 'spec_internal_Energy': 237646.399097292, 'x': None, 'y': None, 'z': None}

Summary

Step	Code
Define your function	def my_props(values, given, fluid, composition, wanted): ... return [T, P, H, V, S, Q, U]
Create the fluid	fluid = cb.init_fluid("name", [1.0], props="CUSTOM", props_func=my_props)
Use it	fluid.set_state([T, P], "TP") — identical to REFPROP/TREND

Tips

• All units must be SI mass-based: K, Pa, J/kg, m³/kg, J/(kg·K)

• fluid, composition, and wanted are always forwarded — declare them even if unused

• Use Q = -6667.0 to indicate single-phase

• Raise NotImplementedError for input pairs your function does not support

• Keep all model constants inside the function body so the callable is self-contained

• When the callable needs to hold expensive state (e.g. interpolators built from measurement data), use a class with __init__ and __call__