Activity models and mixing rules

Activity models

torch_flash.activity

Excess-Gibbs-energy activity models.

NRTL

NRTL(
    interaction: Tensor,
    nonrandomness: Tensor,
    *,
    trainable: bool = False,
)

Bases: Module

Non-random two-liquid excess Gibbs energy model.

Parameters:

Name	Type	Description	Default
interaction	Tensor	Matrix storing :math:g_{ij}-g_{jj} in J/mol.	required
nonrandomness	Tensor	Dimensionless :math:\alpha_{ij} matrix.	required
trainable	bool	Register the interaction-energy matrix as a trainable parameter.	False

References

Renon and Prausnitz, AIChE Journal 14 (1968), 135-144, doi:10.1002/aic.690140124.

Source code in src/torch_flash/activity/models.py

def __init__(
    self,
    interaction: Tensor,
    nonrandomness: Tensor,
    *,
    trainable: bool = False,
) -> None:
    super().__init__()
    if interaction.shape != nonrandomness.shape or interaction.ndim != 2:
        raise ValueError("NRTL parameter arrays must be equally sized matrices")
    value = interaction.clone()
    if trainable:
        self.interaction = nn.Parameter(value)
    else:
        self.register_buffer("interaction", value)
    self.register_buffer("nonrandomness", nonrandomness.clone())

excess_gibbs_rt

excess_gibbs_rt(
    temperature: Tensor, composition: Tensor
) -> Tensor

Evaluate dimensionless molar excess Gibbs energy.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required
composition	Tensor	Mole fractions with components on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless g^E/(RT) with the broadcast leading shape.

Source code in src/torch_flash/activity/models.py

def excess_gibbs_rt(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Evaluate dimensionless molar excess Gibbs energy.

    Parameters
    ----------
    temperature
        Temperature in K.
    composition
        Mole fractions with components on the final axis.

    Returns
    -------
    Tensor
        Dimensionless ``g^E/(RT)`` with the broadcast leading shape.
    """
    x = normalize_composition(composition)
    tau = self.interaction / (R * temperature[..., None, None])
    weights = torch.exp(-self.nonrandomness * tau)
    denominator = torch.einsum("...k,...ki->...i", x, weights)
    numerator = torch.einsum("...j,...ji,...ji->...i", x, tau, weights)
    return torch.sum(x * numerator / denominator, dim=-1)

log_activity_coefficients

log_activity_coefficients(
    temperature: Tensor, composition: Tensor
) -> Tensor

Return component log activity coefficients.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required
composition	Tensor	Mole fractions with components on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless ln(gamma_i) obtained from derivatives of the extensive excess Gibbs energy.

Source code in src/torch_flash/activity/models.py

def log_activity_coefficients(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Return component log activity coefficients.

    Parameters
    ----------
    temperature
        Temperature in K.
    composition
        Mole fractions with components on the final axis.

    Returns
    -------
    Tensor
        Dimensionless ``ln(gamma_i)`` obtained from derivatives of the
        extensive excess Gibbs energy.
    """
    x = normalize_composition(composition)

    def extensive(moles: Tensor) -> Tensor:
        total = moles.sum(dim=-1)
        fractions = moles / total[..., None]
        return (total * self.excess_gibbs_rt(temperature, fractions)).sum()

    result: Tensor = torch.func.grad(extensive)(x)
    return result

AnchoredHuronVidalNRTL

AnchoredHuronVidalNRTL(
    tau_at_lower_temperature: Tensor,
    tau_at_upper_temperature: Tensor,
    nonrandomness: Tensor,
    covolumes: Tensor,
    *,
    lower_temperature: Tensor | float,
    upper_temperature: Tensor | float,
    trainable_nonrandomness: bool = False,
    nonrandomness_bounds: tuple[float, float] = (0.05, 0.6),
)

Bases: Module

Fit-ready HV-NRTL model anchored at two temperatures.

Directly regressing :math:A_{ij} and :math:B_{ij} in :math:\tau_{ij}=A_{ij}/T+B_{ij} is poorly scaled and can be strongly correlated over a finite experimental temperature interval. This exactly equivalent parameterization stores trainable values of :math:\tau_{ij} at two user-selected temperatures and interpolates linearly in :math:1/T.

The optional trainable non-randomness matrix is kept symmetric and within explicit open bounds by a sigmoid transform. Diagonal interactions remain exactly zero. Call :meth:freeze after fitting to obtain the standard :class:HuronVidalNRTL representation used by parameter databases.

Parameters:

Name	Type	Description	Default
tau_at_lower_temperature	Tensor	Dimensionless interaction matrices at the two anchors.	required
tau_at_upper_temperature	Tensor	Dimensionless interaction matrices at the two anchors.	required
nonrandomness	Tensor	Fixed or initially bounded-symmetric nonrandomness matrix.	required
covolumes	Tensor	Positive component cubic covolumes in m3/mol.	required
lower_temperature	Tensor | float	Positive increasing anchor temperatures in K.	required
upper_temperature	Tensor | float	Positive increasing anchor temperatures in K.	required
trainable_nonrandomness	bool	Fit the symmetric nonrandomness matrix through a bounded transform.	False
nonrandomness_bounds	tuple[float, float]	Open lower/upper bounds for fitted off-diagonal nonrandomness.	(0.05, 0.6)

Source code in src/torch_flash/activity/models.py

def __init__(
    self,
    tau_at_lower_temperature: Tensor,
    tau_at_upper_temperature: Tensor,
    nonrandomness: Tensor,
    covolumes: Tensor,
    *,
    lower_temperature: Tensor | float,
    upper_temperature: Tensor | float,
    trainable_nonrandomness: bool = False,
    nonrandomness_bounds: tuple[float, float] = (0.05, 0.60),
) -> None:
    super().__init__()
    if (
        tau_at_lower_temperature.ndim != 2
        or tau_at_lower_temperature.shape[0] != tau_at_lower_temperature.shape[1]
        or tau_at_upper_temperature.shape != tau_at_lower_temperature.shape
        or nonrandomness.shape != tau_at_lower_temperature.shape
    ):
        raise ValueError("anchored HV-NRTL arrays must be equally sized square matrices")
    if covolumes.shape != (tau_at_lower_temperature.shape[0],):
        raise ValueError("one positive anchored HV-NRTL covolume is required per component")
    if bool((covolumes <= 0.0).any()):
        raise ValueError("anchored HV-NRTL covolumes must be positive")

    lower = torch.as_tensor(
        lower_temperature,
        dtype=tau_at_lower_temperature.dtype,
        device=tau_at_lower_temperature.device,
    )
    upper = torch.as_tensor(
        upper_temperature,
        dtype=tau_at_lower_temperature.dtype,
        device=tau_at_lower_temperature.device,
    )
    if (
        lower.ndim != 0
        or upper.ndim != 0
        or not bool(
            torch.isfinite(lower) & torch.isfinite(upper) & (lower > 0.0) & (upper > lower)
        )
    ):
        raise ValueError(
            "anchored HV-NRTL temperatures must be finite, positive, and increasing"
        )

    self.raw_tau_at_lower_temperature = nn.Parameter(tau_at_lower_temperature.clone())
    self.raw_tau_at_upper_temperature = nn.Parameter(tau_at_upper_temperature.clone())
    self.register_buffer("lower_temperature", lower.clone())
    self.register_buffer("upper_temperature", upper.clone())
    self.register_buffer("covolumes", covolumes.clone())
    size = tau_at_lower_temperature.shape[0]
    mask = torch.ones(
        (size, size),
        dtype=tau_at_lower_temperature.dtype,
        device=tau_at_lower_temperature.device,
    ) - torch.eye(
        size,
        dtype=tau_at_lower_temperature.dtype,
        device=tau_at_lower_temperature.device,
    )
    self.register_buffer("_off_diagonal", mask)

    if trainable_nonrandomness:
        lower_bound, upper_bound = nonrandomness_bounds
        if (
            not 0.0 <= lower_bound < upper_bound
            or not torch.isfinite(torch.tensor([lower_bound, upper_bound])).all()
        ):
            raise ValueError(
                "anchored HV-NRTL nonrandomness bounds must be finite and increasing"
            )
        off_diagonal_values = nonrandomness[mask.bool()]
        if not bool(
            (off_diagonal_values > lower_bound).all()
            & (off_diagonal_values < upper_bound).all()
        ):
            raise ValueError(
                "initial trainable HV-NRTL nonrandomness must lie inside its bounds"
            )
        if not torch.allclose(nonrandomness, nonrandomness.mT):
            raise ValueError("trainable HV-NRTL nonrandomness must be symmetric")
        scaled = (nonrandomness - lower_bound) / (upper_bound - lower_bound)
        scaled = torch.where(mask.bool(), scaled, torch.full_like(scaled, 0.5))
        raw_nonrandomness = torch.log(scaled) - torch.log1p(-scaled)
        self.raw_nonrandomness = nn.Parameter(raw_nonrandomness)
        self.register_buffer(
            "_fixed_nonrandomness",
            torch.empty(
                0,
                dtype=nonrandomness.dtype,
                device=nonrandomness.device,
            ),
        )
    else:
        self.register_parameter("raw_nonrandomness", None)
        self.register_buffer("_fixed_nonrandomness", nonrandomness.clone())
    self.nonrandomness_bounds = (float(nonrandomness_bounds[0]), float(nonrandomness_bounds[1]))

tau_at_lower_temperature property

tau_at_lower_temperature: Tensor

Return the lower-anchor interaction matrix with an exact zero diagonal.

tau_at_upper_temperature property

tau_at_upper_temperature: Tensor

Return the upper-anchor interaction matrix with an exact zero diagonal.

nonrandomness property

nonrandomness: Tensor

Return the fixed or bounded-symmetric non-randomness matrix.

energy_over_r property

energy_over_r: Tensor

Return :math:A_{ij} in kelvin for tau = A/T + B.

temperature_coefficient property

temperature_coefficient: Tensor

Return dimensionless :math:B_{ij} for tau = A/T + B.

tau_matrix

tau_matrix(temperature: Tensor) -> Tensor

Evaluate interactions by interpolation in inverse temperature.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required

Returns:

Type	Description
Tensor	Dimensionless interaction matrices with trailing component axes.

Raises:

Type	Description
ValueError	If any temperature is nonpositive.

Source code in src/torch_flash/activity/models.py

def tau_matrix(self, temperature: Tensor) -> Tensor:
    """Evaluate interactions by interpolation in inverse temperature.

    Parameters
    ----------
    temperature
        Positive temperature in K.

    Returns
    -------
    Tensor
        Dimensionless interaction matrices with trailing component axes.

    Raises
    ------
    ValueError
        If any temperature is nonpositive.
    """
    if bool((temperature <= 0.0).any()):
        raise ValueError("HV-NRTL temperature must be positive")
    inverse_fraction = (temperature.reciprocal() - self.upper_temperature.reciprocal()) / (
        self.lower_temperature.reciprocal() - self.upper_temperature.reciprocal()
    )
    return (
        inverse_fraction[..., None, None] * self.tau_at_lower_temperature
        + (1.0 - inverse_fraction[..., None, None]) * self.tau_at_upper_temperature
    )

excess_gibbs_rt

excess_gibbs_rt(
    temperature: Tensor, composition: Tensor
) -> Tensor

Evaluate anchored HV dimensionless excess Gibbs energy.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required
composition	Tensor	Mole fractions on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless molar g^E/(RT).

Source code in src/torch_flash/activity/models.py

def excess_gibbs_rt(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Evaluate anchored HV dimensionless excess Gibbs energy.

    Parameters
    ----------
    temperature
        Positive temperature in K.
    composition
        Mole fractions on the final axis.

    Returns
    -------
    Tensor
        Dimensionless molar ``g^E/(RT)``.
    """
    x = normalize_composition(composition)
    return _huron_vidal_excess_gibbs_rt(
        x,
        self.tau_matrix(temperature),
        self.nonrandomness,
        self.covolumes,
    )

log_activity_coefficients

log_activity_coefficients(
    temperature: Tensor, composition: Tensor
) -> Tensor

Return anchored HV component log activity coefficients.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required
composition	Tensor	Mole fractions on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless ln(gamma_i).

Source code in src/torch_flash/activity/models.py

def log_activity_coefficients(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Return anchored HV component log activity coefficients.

    Parameters
    ----------
    temperature
        Positive temperature in K.
    composition
        Mole fractions on the final axis.

    Returns
    -------
    Tensor
        Dimensionless ``ln(gamma_i)``.
    """
    x = normalize_composition(composition)

    def extensive(moles: Tensor) -> Tensor:
        total = moles.sum(dim=-1)
        fractions = moles / total[..., None]
        return (total * self.excess_gibbs_rt(temperature, fractions)).sum()

    result: Tensor = torch.func.grad(extensive)(x)
    return result

freeze

freeze() -> HuronVidalNRTL

Convert fitted anchor parameters to a detached standard HV model.

Returns:

Type	Description
HuronVidalNRTL	Prediction/storage model containing detached A, B, nonrandomness, and covolume tensors.

Source code in src/torch_flash/activity/models.py

def freeze(self) -> HuronVidalNRTL:
    """Convert fitted anchor parameters to a detached standard HV model.

    Returns
    -------
    HuronVidalNRTL
        Prediction/storage model containing detached ``A``, ``B``,
        nonrandomness, and covolume tensors.
    """
    return HuronVidalNRTL(
        self.energy_over_r.detach(),
        self.temperature_coefficient.detach(),
        self.nonrandomness.detach(),
        self.covolumes.detach(),
    )

HuronVidalNRTL

HuronVidalNRTL(
    energy_over_r: Tensor,
    temperature_coefficient: Tensor,
    nonrandomness: Tensor,
    covolumes: Tensor,
    *,
    trainable: bool = False,
)

Bases: Module

Covolume-weighted NRTL model used by the original HV formulation.

The dimensionless interaction energy is

.. math::

\tau_{ji}(T) = A_{ji}/T + B_{ji},

where energy_over_r stores :math:A in kelvin and temperature_coefficient stores the dimensionless linear coefficient :math:B. The excess Gibbs energy is the modified NRTL expression from Huron and Vidal (1979), with :math:x_j b_j replacing :math:x_j in the local-composition sums. Consequently its parameters must be fitted for the HV infinite-pressure reference; ordinary low-pressure NRTL parameters are not interchangeable.

Reference: Huron and Vidal, Fluid Phase Equilibria 3 (1979), 255-271, doi:10.1016/0378-3812(79)80001-1.

Parameters:

Name	Type	Description	Default
energy_over_r	Tensor	Matrix A in K.	required
temperature_coefficient	Tensor	Dimensionless matrix B.	required
nonrandomness	Tensor	Dimensionless NRTL nonrandomness matrix.	required
covolumes	Tensor	Positive component cubic covolumes in m3/mol.	required
trainable	bool	Register A and B as trainable parameters.	False

Source code in src/torch_flash/activity/models.py

def __init__(
    self,
    energy_over_r: Tensor,
    temperature_coefficient: Tensor,
    nonrandomness: Tensor,
    covolumes: Tensor,
    *,
    trainable: bool = False,
) -> None:
    super().__init__()
    if (
        energy_over_r.ndim != 2
        or energy_over_r.shape[0] != energy_over_r.shape[1]
        or temperature_coefficient.shape != energy_over_r.shape
        or nonrandomness.shape != energy_over_r.shape
    ):
        raise ValueError("HV-NRTL interaction arrays must be equally sized square matrices")
    if covolumes.shape != (energy_over_r.shape[0],):
        raise ValueError("one positive HV-NRTL covolume is required per component")
    if bool((covolumes <= 0.0).any()):
        raise ValueError("HV-NRTL covolumes must be positive")
    if trainable:
        self.energy_over_r = nn.Parameter(energy_over_r.clone())
        self.temperature_coefficient = nn.Parameter(temperature_coefficient.clone())
    else:
        self.register_buffer("energy_over_r", energy_over_r.clone())
        self.register_buffer("temperature_coefficient", temperature_coefficient.clone())
    self.register_buffer("nonrandomness", nonrandomness.clone())
    self.register_buffer("covolumes", covolumes.clone())

tau_matrix

tau_matrix(temperature: Tensor) -> Tensor

Evaluate the dimensionless HV-NRTL interaction matrix.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required

Returns:

Type	Description
Tensor	tau(T) = A/T + B with trailing component-pair axes.

Source code in src/torch_flash/activity/models.py

def tau_matrix(self, temperature: Tensor) -> Tensor:
    """Evaluate the dimensionless HV-NRTL interaction matrix.

    Parameters
    ----------
    temperature
        Temperature in K.

    Returns
    -------
    Tensor
        ``tau(T) = A/T + B`` with trailing component-pair axes.
    """
    return self.energy_over_r / temperature[..., None, None] + self.temperature_coefficient

excess_gibbs_rt

excess_gibbs_rt(
    temperature: Tensor, composition: Tensor
) -> Tensor

Evaluate covolume-weighted HV excess Gibbs energy.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required
composition	Tensor	Mole fractions on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless molar g^E/(RT).

Source code in src/torch_flash/activity/models.py

def excess_gibbs_rt(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Evaluate covolume-weighted HV excess Gibbs energy.

    Parameters
    ----------
    temperature
        Temperature in K.
    composition
        Mole fractions on the final axis.

    Returns
    -------
    Tensor
        Dimensionless molar ``g^E/(RT)``.
    """
    x = normalize_composition(composition)
    tau = self.tau_matrix(temperature)
    return _huron_vidal_excess_gibbs_rt(
        x,
        tau,
        self.nonrandomness,
        self.covolumes,
    )

log_activity_coefficients

log_activity_coefficients(
    temperature: Tensor, composition: Tensor
) -> Tensor

Return HV-NRTL component log activity coefficients.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required
composition	Tensor	Mole fractions on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless ln(gamma_i) from the extensive HV excess Gibbs energy derivative.

Source code in src/torch_flash/activity/models.py

def log_activity_coefficients(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Return HV-NRTL component log activity coefficients.

    Parameters
    ----------
    temperature
        Temperature in K.
    composition
        Mole fractions on the final axis.

    Returns
    -------
    Tensor
        Dimensionless ``ln(gamma_i)`` from the extensive HV excess Gibbs
        energy derivative.
    """
    x = normalize_composition(composition)

    def extensive(moles: Tensor) -> Tensor:
        total = moles.sum(dim=-1)
        fractions = moles / total[..., None]
        return (total * self.excess_gibbs_rt(temperature, fractions)).sum()

    result: Tensor = torch.func.grad(extensive)(x)
    return result

Wilson

Wilson(
    interaction: Tensor,
    molar_volumes: Tensor,
    *,
    trainable: bool = False,
)

Bases: Module

Wilson excess Gibbs model with energy and molar-volume parameters.

Parameters:

Name	Type	Description	Default
interaction	Tensor	Square energy-difference matrix in J/mol.	required
molar_volumes	Tensor	Positive pure-liquid molar volumes in m3/mol.	required
trainable	bool	Register interaction as a trainable PyTorch parameter.	False

References

Wilson, Journal of the American Chemical Society 86 (1964), 127-130, doi:10.1021/ja01056a002.

Source code in src/torch_flash/activity/models.py

def __init__(
    self,
    interaction: Tensor,
    molar_volumes: Tensor,
    *,
    trainable: bool = False,
) -> None:
    super().__init__()
    if interaction.ndim != 2 or interaction.shape[0] != interaction.shape[1]:
        raise ValueError("Wilson interaction must be a square matrix")
    if molar_volumes.shape != (interaction.shape[0],):
        raise ValueError("one Wilson molar volume is required per component")
    if trainable:
        self.interaction = nn.Parameter(interaction.clone())
    else:
        self.register_buffer("interaction", interaction.clone())
    self.register_buffer("molar_volumes", molar_volumes.clone())

lambda_matrix

lambda_matrix(temperature: Tensor) -> Tensor

Evaluate the temperature-dependent Wilson Lambda matrix.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required

Returns:

Type	Description
Tensor	Dimensionless Lambda matrices with trailing component-pair axes.

Source code in src/torch_flash/activity/models.py

def lambda_matrix(self, temperature: Tensor) -> Tensor:
    """Evaluate the temperature-dependent Wilson Lambda matrix.

    Parameters
    ----------
    temperature
        Temperature in K.

    Returns
    -------
    Tensor
        Dimensionless Lambda matrices with trailing component-pair axes.
    """
    volume_ratio = self.molar_volumes[None, :] / self.molar_volumes[:, None]
    return volume_ratio * torch.exp(-self.interaction / (R * temperature[..., None, None]))

excess_gibbs_rt

excess_gibbs_rt(
    temperature: Tensor, composition: Tensor
) -> Tensor

Evaluate dimensionless molar excess Gibbs energy.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required
composition	Tensor	Mole fractions on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless g^E/(RT).

Source code in src/torch_flash/activity/models.py

def excess_gibbs_rt(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Evaluate dimensionless molar excess Gibbs energy.

    Parameters
    ----------
    temperature
        Temperature in K.
    composition
        Mole fractions on the final axis.

    Returns
    -------
    Tensor
        Dimensionless ``g^E/(RT)``.
    """
    x = normalize_composition(composition)
    lambda_matrix = self.lambda_matrix(temperature)
    sums = torch.einsum("...j,...ij->...i", x, lambda_matrix)
    return -torch.sum(x * torch.log(sums), dim=-1)

log_activity_coefficients

log_activity_coefficients(
    temperature: Tensor, composition: Tensor
) -> Tensor

Evaluate analytical Wilson log activity coefficients.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Temperature in K.	required
composition	Tensor	Mole fractions on the final axis.	required

Returns:

Type	Description
Tensor	Dimensionless ln(gamma_i).

Source code in src/torch_flash/activity/models.py

def log_activity_coefficients(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Evaluate analytical Wilson log activity coefficients.

    Parameters
    ----------
    temperature
        Temperature in K.
    composition
        Mole fractions on the final axis.

    Returns
    -------
    Tensor
        Dimensionless ``ln(gamma_i)``.
    """
    x = normalize_composition(composition)
    lambda_matrix = self.lambda_matrix(temperature)
    row_sums = torch.einsum("...j,...ij->...i", x, lambda_matrix)
    correction = torch.einsum("...j,...ji,...j->...i", x, lambda_matrix, row_sums.reciprocal())
    return 1.0 - torch.log(row_sums) - correction

UNIFAC

UNIFAC(
    group_counts: Tensor,
    relative_volume: Tensor,
    relative_surface_area: Tensor,
    subgroup_main_indices: Tensor,
    interaction: Tensor,
    *,
    coordination_number: float = 10.0,
    trainable: bool = False,
)

Bases: Module

Original UNIFAC excess-Gibbs-energy model.

Parameters:

Name	Type	Description	Default
group_counts	Tensor	Component-by-subgroup occurrence matrix :math:\nu_k^{(i)}.	required
relative_volume	Tensor	Selected subgroup :math:R_k and :math:Q_k values.	required
relative_surface_area	Tensor	Selected subgroup :math:R_k and :math:Q_k values.	required
subgroup_main_indices	Tensor	Zero-based index of the selected main group for each subgroup.	required
interaction	Tensor	Directed selected-main-group matrix :math:a_{mn} in kelvin. Every required off-diagonal value must be supplied.	required
coordination_number	float	Lattice coordination number :math:z; original UNIFAC uses 10.	10.0
trainable	bool	Store the main-group interaction matrix as an nn.Parameter.	False

Notes

The model is vectorized over all leading dimensions of temperature and composition. PyTorch autodiff can therefore differentiate activities, excess Gibbs energy, and temperature/composition derivatives on CPU or GPU. Group assignment is discrete and is intentionally not trainable.

Source code in src/torch_flash/activity/unifac.py

def __init__(
    self,
    group_counts: Tensor,
    relative_volume: Tensor,
    relative_surface_area: Tensor,
    subgroup_main_indices: Tensor,
    interaction: Tensor,
    *,
    coordination_number: float = 10.0,
    trainable: bool = False,
) -> None:
    super().__init__()
    if group_counts.ndim != 2 or not group_counts.shape[0] or not group_counts.shape[1]:
        raise ValueError("UNIFAC group counts must be a nonempty component-by-group matrix")
    group_shape = (group_counts.shape[1],)
    if relative_volume.shape != group_shape or relative_surface_area.shape != group_shape:
        raise ValueError("UNIFAC requires one R and Q value per selected subgroup")
    if subgroup_main_indices.shape != group_shape:
        raise ValueError("UNIFAC requires one main-group index per selected subgroup")
    if subgroup_main_indices.dtype != torch.long:
        raise ValueError("UNIFAC main-group indices must have torch.long dtype")
    if interaction.ndim != 2 or interaction.shape[0] != interaction.shape[1]:
        raise ValueError("UNIFAC interactions must be a square directed matrix")
    if not interaction.shape[0]:
        raise ValueError("UNIFAC requires at least one selected main group")
    tensors = (group_counts, relative_volume, relative_surface_area, interaction)
    if not all(bool(torch.isfinite(value).all()) for value in tensors):
        raise ValueError("UNIFAC numerical inputs must be finite")
    if bool((group_counts < 0.0).any()) or bool((group_counts.sum(dim=-1) <= 0.0).any()):
        raise ValueError("each UNIFAC component needs nonnegative, nonempty group counts")
    if bool((relative_volume <= 0.0).any()) or bool((relative_surface_area < 0.0).any()):
        raise ValueError("UNIFAC R values must be positive and Q values nonnegative")
    if bool((subgroup_main_indices < 0).any()) or bool(
        (subgroup_main_indices >= interaction.shape[0]).any()
    ):
        raise ValueError("UNIFAC main-group indices are outside the interaction matrix")
    if not torch.allclose(
        torch.diagonal(interaction),
        torch.zeros(interaction.shape[0], dtype=interaction.dtype, device=interaction.device),
    ):
        raise ValueError("UNIFAC same-main-group interactions must be zero")
    if coordination_number <= 0.0:
        raise ValueError("UNIFAC coordination number must be positive")

    molecular_volume = group_counts @ relative_volume
    molecular_surface = group_counts @ relative_surface_area
    if bool((molecular_volume <= 0.0).any()) or bool((molecular_surface <= 0.0).any()):
        raise ValueError("each UNIFAC component must have positive molecular R and Q")

    self.coordination_number = float(coordination_number)
    self.subgroup_keys = ()
    self.main_group_ids = ()
    self.register_buffer("group_counts", group_counts.clone())
    self.register_buffer("relative_volume", relative_volume.clone())
    self.register_buffer("relative_surface_area", relative_surface_area.clone())
    self.register_buffer("subgroup_main_indices", subgroup_main_indices.clone())
    if trainable:
        self.interaction = nn.Parameter(interaction.clone())
    else:
        self.register_buffer("interaction", interaction.clone())

molecular_relative_volume property

molecular_relative_volume: Tensor

Return component molecular relative-volume parameters.

Returns:

Type	Description
Tensor	Dimensionless r_i obtained from subgroup counts and R_k.

molecular_relative_surface_area property

molecular_relative_surface_area: Tensor

Return component molecular relative-surface parameters.

Returns:

Type	Description
Tensor	Dimensionless q_i obtained from subgroup counts and Q_k.

interaction_factors

interaction_factors(temperature: Tensor) -> Tensor

Evaluate directed main-group interaction factors.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required

Returns:

Type	Description
Tensor	Dimensionless :math:\Psi_{mn}=\exp(-a_{mn}/T) with trailing main-group axes.

Source code in src/torch_flash/activity/unifac.py

def interaction_factors(self, temperature: Tensor) -> Tensor:
    r"""Evaluate directed main-group interaction factors.

    Parameters
    ----------
    temperature
        Positive temperature in K.

    Returns
    -------
    Tensor
        Dimensionless :math:`\Psi_{mn}=\exp(-a_{mn}/T)` with trailing
        main-group axes.
    """
    if not torch.compiler.is_compiling() and bool((temperature <= 0.0).any()):
        raise ValueError("temperature must be positive")
    interaction = self.interaction - torch.diag_embed(torch.diagonal(self.interaction))
    return torch.exp(-interaction / temperature[..., None, None])

combinatorial_log_activity_coefficients

combinatorial_log_activity_coefficients(
    composition: Tensor,
) -> Tensor

Evaluate the original UNIFAC combinatorial contribution.

Parameters:

Name	Type	Description	Default
composition	Tensor	Mole fractions on the final component axis.	required

Returns:

Type	Description
Tensor	Component contributions to dimensionless ln(gamma_i).

Source code in src/torch_flash/activity/unifac.py

def combinatorial_log_activity_coefficients(
    self,
    composition: Tensor,
) -> Tensor:
    r"""Evaluate the original UNIFAC combinatorial contribution.

    Parameters
    ----------
    composition
        Mole fractions on the final component axis.

    Returns
    -------
    Tensor
        Component contributions to dimensionless ``ln(gamma_i)``.
    """
    x = normalize_composition(composition)
    if x.shape[-1] != self.group_counts.shape[0]:
        raise ValueError("UNIFAC composition size must match the component group matrix")
    r = self.molecular_relative_volume
    q = self.molecular_relative_surface_area
    r_bar = torch.sum(x * r, dim=-1, keepdim=True)
    q_bar = torch.sum(x * q, dim=-1, keepdim=True)
    phi_over_x = r / r_bar
    theta_over_phi = (q / r) * (r_bar / q_bar)
    return (
        torch.log(phi_over_x)
        + 1.0
        - phi_over_x
        + 0.5
        * self.coordination_number
        * q
        * (torch.log(theta_over_phi) - 1.0 + theta_over_phi.reciprocal())
    )

residual_log_activity_coefficients

residual_log_activity_coefficients(
    temperature: Tensor, composition: Tensor
) -> Tensor

Evaluate the original UNIFAC residual group contribution.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required
composition	Tensor	Mole fractions on the final component axis.	required

Returns:

Type	Description
Tensor	Component residual contributions to ln(gamma_i).

Source code in src/torch_flash/activity/unifac.py

def residual_log_activity_coefficients(
    self,
    temperature: Tensor,
    composition: Tensor,
) -> Tensor:
    r"""Evaluate the original UNIFAC residual group contribution.

    Parameters
    ----------
    temperature
        Positive temperature in K.
    composition
        Mole fractions on the final component axis.

    Returns
    -------
    Tensor
        Component residual contributions to ``ln(gamma_i)``.
    """
    temperature, x = self._state(temperature, composition)
    group_abundance = torch.einsum("...i,ik->...k", x, self.group_counts)
    mixture_surface = self.relative_surface_area * group_abundance
    mixture_surface = mixture_surface / mixture_surface.sum(dim=-1, keepdim=True)

    main_factors = self.interaction_factors(temperature)
    subgroup_factors = main_factors[
        ...,
        self.subgroup_main_indices[:, None],
        self.subgroup_main_indices[None, :],
    ]
    mixture_log_group = self._group_log_activity_coefficients(
        mixture_surface,
        subgroup_factors,
    )

    pure_surface = self.group_counts * self.relative_surface_area
    pure_surface = pure_surface / pure_surface.sum(dim=-1, keepdim=True)
    batch_shape = x.shape[:-1]
    pure_surface = pure_surface.expand((*batch_shape, *pure_surface.shape))
    pure_log_group = self._group_log_activity_coefficients(
        pure_surface,
        subgroup_factors[..., None, :, :],
    )
    return torch.sum(
        self.group_counts * (mixture_log_group[..., None, :] - pure_log_group),
        dim=-1,
    )

log_activity_coefficients

log_activity_coefficients(
    temperature: Tensor, composition: Tensor
) -> Tensor

Return original-UNIFAC component log activity coefficients.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required
composition	Tensor	Mole fractions on the final component axis.	required

Returns:

Type	Description
Tensor	Sum of combinatorial and residual ln(gamma_i) contributions.

Source code in src/torch_flash/activity/unifac.py

def log_activity_coefficients(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Return original-UNIFAC component log activity coefficients.

    Parameters
    ----------
    temperature
        Positive temperature in K.
    composition
        Mole fractions on the final component axis.

    Returns
    -------
    Tensor
        Sum of combinatorial and residual ``ln(gamma_i)`` contributions.
    """
    _, x = self._state(temperature, composition)
    return self.combinatorial_log_activity_coefficients(x) + (
        self.residual_log_activity_coefficients(temperature, x)
    )

excess_gibbs_rt

excess_gibbs_rt(
    temperature: Tensor, composition: Tensor
) -> Tensor

Return dimensionless original-UNIFAC excess Gibbs energy.

Parameters:

Name	Type	Description	Default
temperature	Tensor	Positive temperature in K.	required
composition	Tensor	Mole fractions on the final component axis.	required

Returns:

Type	Description
Tensor	sum_i x_i ln(gamma_i) = g^E/(RT).

Source code in src/torch_flash/activity/unifac.py

def excess_gibbs_rt(self, temperature: Tensor, composition: Tensor) -> Tensor:
    """Return dimensionless original-UNIFAC excess Gibbs energy.

    Parameters
    ----------
    temperature
        Positive temperature in K.
    composition
        Mole fractions on the final component axis.

    Returns
    -------
    Tensor
        ``sum_i x_i ln(gamma_i) = g^E/(RT)``.
    """
    _, x = self._state(temperature, composition)
    return torch.sum(x * self.log_activity_coefficients(temperature, x), dim=-1)

activity_model

activity_model(
    parameter_set: ParameterSource,
    names: tuple[str, ...] | None = None,
    *,
    dtype: dtype | None = None,
    device: device | str | None = None,
    trainable: bool = False,
    covolumes: Tensor | None = None,
    group_assignments: Sequence[GroupAssignment]
    | None = None,
) -> ActivityModel

Construct NRTL, HV-NRTL, Wilson, or original UNIFAC from YAML.

Parameters:

Name	Type	Description	Default
parameter_set	ParameterSource	Bundled identifier, YAML path, mapping, or loaded activity parameter set.	required
names	tuple[str, ...] | None	Optional requested component order. Omitting it uses the stored order.	None
dtype	dtype | None	Tensor placement, defaulting to runtime configuration.	None
device	dtype | None	Tensor placement, defaulting to runtime configuration.	None
trainable	bool	Register supported interaction tensors as PyTorch parameters.	False
covolumes	Tensor | None	Optional HV-NRTL component covolumes in m3/mol.	None
group_assignments	Sequence[GroupAssignment] | None	Optional explicit UNIFAC subgroup counts, one mapping per component.	None

Returns:

Type	Description
ActivityModel	Model matching the parameter document's exact identity.

Raises:

Type	Description
KeyError	If requested components or UNIFAC groups are unsupported.
ValueError	If component order or model-specific overrides are inconsistent.
ParameterDatabaseError	If model identity, units, or parameter records are malformed.

Notes

The requested component order may differ from the stored order; vectors and matrices are permuted consistently. For custom fitting workflows, the model classes continue to accept explicit tensors directly. UNIFAC uses bundled fragmentations selected by names or explicit group_assignments.

Source code in src/torch_flash/activity/named.py

def activity_model(
    parameter_set: ParameterSource,
    names: tuple[str, ...] | None = None,
    *,
    dtype: torch.dtype | None = None,
    device: torch.device | str | None = None,
    trainable: bool = False,
    covolumes: Tensor | None = None,
    group_assignments: Sequence[GroupAssignment] | None = None,
) -> ActivityModel:
    """Construct NRTL, HV-NRTL, Wilson, or original UNIFAC from YAML.

    Parameters
    ----------
    parameter_set
        Bundled identifier, YAML path, mapping, or loaded activity parameter
        set.
    names
        Optional requested component order. Omitting it uses the stored order.
    dtype, device
        Tensor placement, defaulting to runtime configuration.
    trainable
        Register supported interaction tensors as PyTorch parameters.
    covolumes
        Optional HV-NRTL component covolumes in m3/mol.
    group_assignments
        Optional explicit UNIFAC subgroup counts, one mapping per component.

    Returns
    -------
    ActivityModel
        Model matching the parameter document's exact identity.

    Raises
    ------
    KeyError
        If requested components or UNIFAC groups are unsupported.
    ValueError
        If component order or model-specific overrides are inconsistent.
    ParameterDatabaseError
        If model identity, units, or parameter records are malformed.

    Notes
    -----
    The requested component order may differ from the stored order; vectors
    and matrices are permuted consistently. For custom fitting workflows, the
    model classes continue to accept explicit tensors directly. UNIFAC uses
    bundled fragmentations selected by ``names`` or explicit
    ``group_assignments``.
    """
    dtype, device = resolve_tensor_options(dtype, device)
    loaded = load_model_parameters(parameter_set)
    if loaded.model_kind != "activity":
        raise ParameterDatabaseError(
            f"{loaded.identifier!r} is {loaded.model_kind!r}, not 'activity'"
        )
    parameters = loaded.parameters
    model_name = loaded.model.strip().lower().replace("_", "-")
    if model_name in ("unifac", "original-unifac"):
        return _unifac_from_loaded(
            loaded,
            names=names,
            group_assignments=group_assignments,
            dtype=dtype,
            device=device,
            trainable=trainable,
        )
    if group_assignments is not None:
        raise ValueError("group_assignments are only valid for UNIFAC")
    stored_names = parameters.get("components")
    if not isinstance(stored_names, Sequence) or isinstance(stored_names, str):
        raise ParameterDatabaseError(f"{loaded.identifier!r} requires a components list")
    canonical_stored = tuple(
        canonical_component_name(str(name), strict=False) for name in stored_names
    )
    selected = (
        canonical_stored
        if names is None
        else tuple(canonical_component_name(name, strict=False) for name in names)
    )
    if not selected or len(set(selected)) != len(selected):
        raise ValueError("activity-model component names must be non-empty and unique")
    try:
        order = tuple(canonical_stored.index(name) for name in selected)
    except ValueError as exc:
        raise KeyError(
            f"{loaded.identifier} has no activity parameters for one or more of {selected!r}"
        ) from exc

    if model_name == "nrtl":
        return NRTL(
            _numeric_matrix(parameters, "interaction", order, dtype=dtype, device=device),
            _numeric_matrix(parameters, "nonrandomness", order, dtype=dtype, device=device),
            trainable=trainable,
        )
    if model_name == "wilson":
        return Wilson(
            _numeric_matrix(parameters, "interaction", order, dtype=dtype, device=device),
            _numeric_vector(parameters, "molar_volumes", order, dtype=dtype, device=device),
            trainable=trainable,
        )
    if model_name not in ("hv-nrtl", "huron-vidal-nrtl"):
        raise ParameterDatabaseError(
            f"{loaded.identifier!r} has unsupported activity model {loaded.model!r}"
        )

    if covolumes is None and "covolumes" in parameters:
        covolumes = _numeric_vector(
            parameters,
            "covolumes",
            order,
            dtype=dtype,
            device=device,
        )
    if covolumes is None:
        cubic_source = parameters.get("covolume_cubic_parameter_set")
        if not isinstance(cubic_source, str):
            raise ParameterDatabaseError(
                f"{loaded.identifier!r} requires covolumes or 'covolume_cubic_parameter_set'"
            )
        cubic = load_model_parameters(cubic_source)
        omega_b = cubic.parameters.get("omega_b")
        if cubic.model_kind != "cubic" or not isinstance(omega_b, int | float):
            raise ParameterDatabaseError(f"{cubic.identifier!r} cannot provide cubic covolumes")
        components = component_set(selected, dtype=dtype, device=device)
        covolumes = (
            float(omega_b) * R * components.critical_temperature / components.critical_pressure
        )
    else:
        covolumes = covolumes.to(dtype=dtype, device=device)
    if covolumes.shape != (len(selected),):
        raise ValueError("one activity-model covolume is required per selected component")
    return HuronVidalNRTL(
        _numeric_matrix(parameters, "energy_over_r", order, dtype=dtype, device=device),
        _numeric_matrix(
            parameters,
            "temperature_coefficient",
            order,
            dtype=dtype,
            device=device,
        ),
        _numeric_matrix(parameters, "nonrandomness", order, dtype=dtype, device=device),
        covolumes,
        trainable=trainable,
    )

unifac_groups_from_identifiers

unifac_groups_from_identifiers(
    identifiers: Sequence[str],
    *,
    identifier_type: str = "smiles",
) -> tuple[dict[int, float], ...]

Fragment molecules with the optional MIT-licensed ugropy package.

Parameters:

Name	Type	Description	Default
identifiers	Sequence[str]	Molecular SMILES strings or names.	required
identifier_type	str	"smiles" or "name".	'smiles'

Returns:

Type	Description
tuple	One mapping of published original-UNIFAC subgroup number to occurrence count per identifier.

Raises:

Type	Description
ImportError	If the optional groups dependency is not installed.
ValueError	If the identifier type is unsupported or fragmentation is ambiguous or unavailable.

Notes

SMILES inputs are recommended for reproducibility. Name lookup delegates to PubChem and therefore requires network access and may change outside torch-flash. Fragmentation is discrete; inspect every returned map before using it in regression or safety-critical calculations.

Source code in src/torch_flash/activity/unifac.py

def unifac_groups_from_identifiers(
    identifiers: Sequence[str],
    *,
    identifier_type: str = "smiles",
) -> tuple[dict[int, float], ...]:
    """Fragment molecules with the optional MIT-licensed ``ugropy`` package.

    Parameters
    ----------
    identifiers
        Molecular SMILES strings or names.
    identifier_type
        ``"smiles"`` or ``"name"``.

    Returns
    -------
    tuple
        One mapping of published original-UNIFAC subgroup number to occurrence
        count per identifier.

    Raises
    ------
    ImportError
        If the optional ``groups`` dependency is not installed.
    ValueError
        If the identifier type is unsupported or fragmentation is ambiguous or
        unavailable.

    Notes
    -----
    SMILES inputs are recommended for reproducibility.  Name lookup delegates
    to PubChem and therefore requires network access and may change outside
    ``torch-flash``.  Fragmentation is discrete; inspect every returned map
    before using it in regression or safety-critical calculations.
    """
    try:
        from ugropy import unifac as ugropy_unifac
    except ImportError as exc:
        raise ImportError(
            "automatic UNIFAC fragmentation requires 'ugropy'; install torch-flash[groups]"
        ) from exc
    if identifier_type not in ("name", "smiles"):
        raise ValueError("identifier_type must be 'name' or 'smiles'")
    assignments: list[dict[int, float]] = []
    for identifier in identifiers:
        result = ugropy_unifac.get_groups(identifier, identifier_type)
        if isinstance(result, list):
            if len(result) != 1:
                raise ValueError(
                    f"ugropy returned {len(result)} fragmentations for {identifier!r}; "
                    "select a group assignment explicitly"
                )
            result = result[0]
        raw = getattr(result, "subgroups_num", None)
        if not isinstance(raw, Mapping) or not raw:
            raise ValueError(f"ugropy could not assign original-UNIFAC groups to {identifier!r}")
        assignments.append({int(key): float(value) for key, value in raw.items()})
    return tuple(assignments)

unifac_model

unifac_model(
    parameter_set: ParameterSource = "unifac-original",
    names: Sequence[str] | None = None,
    *,
    group_assignments: Sequence[GroupAssignment]
    | None = None,
    dtype: dtype | None = None,
    device: device | str | None = None,
    trainable: bool = False,
) -> UNIFAC

Construct original UNIFAC from cached YAML or explicit parameters.

Parameters:

Name	Type	Description	Default
parameter_set	ParameterSource	Original-UNIFAC parameter source.	'unifac-original'
names	Sequence[str] | None	Optional component names selecting bundled audited fragmentations.	None
group_assignments	Sequence[GroupAssignment] | None	Optional explicit subgroup counts, one mapping per component.	None
dtype	dtype | None	Tensor placement.	None
device	dtype | None	Tensor placement.	None
trainable	bool	Register the directed main-group interaction matrix as a trainable parameter.	False

Returns:

Type	Description
UNIFAC	Original-UNIFAC model containing only the groups required by the selected components.

Raises:

Type	Description
KeyError	If a component or subgroup is absent.
ValueError	If names and assignments are inconsistent.
ParameterDatabaseError	If the parameter set is not original UNIFAC or lacks required directed group interactions.

Notes

names select bundled, audited component fragmentations. Arbitrary molecules use group_assignments, where each mapping may be keyed by subgroup key, published subgroup number, or an unambiguous subgroup name. The optional :func:unifac_groups_from_identifiers adapter can generate these mappings with ugropy.

Source code in src/torch_flash/activity/unifac.py

def unifac_model(
    parameter_set: ParameterSource = "unifac-original",
    names: Sequence[str] | None = None,
    *,
    group_assignments: Sequence[GroupAssignment] | None = None,
    dtype: torch.dtype | None = None,
    device: torch.device | str | None = None,
    trainable: bool = False,
) -> UNIFAC:
    """Construct original UNIFAC from cached YAML or explicit parameters.

    Parameters
    ----------
    parameter_set
        Original-UNIFAC parameter source.
    names
        Optional component names selecting bundled audited fragmentations.
    group_assignments
        Optional explicit subgroup counts, one mapping per component.
    dtype, device
        Tensor placement.
    trainable
        Register the directed main-group interaction matrix as a trainable
        parameter.

    Returns
    -------
    UNIFAC
        Original-UNIFAC model containing only the groups required by the
        selected components.

    Raises
    ------
    KeyError
        If a component or subgroup is absent.
    ValueError
        If names and assignments are inconsistent.
    ParameterDatabaseError
        If the parameter set is not original UNIFAC or lacks required directed
        group interactions.

    Notes
    -----
    ``names`` select bundled, audited component fragmentations.  Arbitrary
    molecules use ``group_assignments``, where each mapping may be keyed by
    subgroup key, published subgroup number, or an unambiguous subgroup name.
    The optional :func:`unifac_groups_from_identifiers` adapter can generate
    these mappings with ``ugropy``.
    """
    dtype, device = resolve_tensor_options(dtype, device)
    loaded = load_model_parameters(parameter_set)
    if loaded.model_kind != "activity":
        raise ParameterDatabaseError(
            f"{loaded.identifier!r} is {loaded.model_kind!r}, not 'activity'"
        )
    if loaded.model.strip().lower().replace("_", "-") not in ("unifac", "original-unifac"):
        raise ParameterDatabaseError(
            f"{loaded.identifier!r} is {loaded.model!r}, not original UNIFAC"
        )
    return _unifac_from_loaded(
        loaded,
        names=names,
        group_assignments=group_assignments,
        dtype=dtype,
        device=device,
        trainable=trainable,
    )

Cubic mixing rules

torch_flash.mixing

Mixing rules for equations of state.

HuronVidalMixing

HuronVidalMixing(
    activity_model: ActivityModel,
    *,
    delta1: float,
    delta2: float,
)

Bases: Module

Original infinite-pressure Huron-Vidal mixing rule.

The implementation follows Michelsen and Mollerup, chapter 7, Eq. 7:

a/(bRT) = sum(x_i a_i/(b_i RT)) - gE/(Delta RT).

Primary source: Huron and Vidal (1979), doi:10.1016/0378-3812(79)80001-1. The implementation convention follows Michelsen and Mollerup, Thermodynamic Models, 2nd ed. (2007), chapter 7, ISBN 978-87-989961-3-2.

Parameters:

Name	Type	Description	Default
activity_model	ActivityModel	Excess-Gibbs model compatible with the HV infinite-pressure reference.	required
delta1	float	Distinct generalized-cubic denominator-shape constants.	required
delta2	float	Distinct generalized-cubic denominator-shape constants.	required

Source code in src/torch_flash/mixing/rules.py

def __init__(
    self,
    activity_model: ActivityModel,
    *,
    delta1: float,
    delta2: float,
) -> None:
    super().__init__()
    if delta1 == delta2:
        raise ValueError("cubic delta parameters must be distinct")
    self.activity_model = activity_model  # type: ignore[assignment]
    self.delta1 = float(delta1)
    self.delta2 = float(delta2)
    ratio = (1.0 + delta2) / (1.0 + delta1)
    self.delta = float(torch.log(torch.tensor(ratio)) / (delta2 - delta1))

forward

forward(
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]

Return mixture a and b parameters.

Source code in src/torch_flash/mixing/rules.py

def forward(
    self,
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]:
    """Return mixture ``a`` and ``b`` parameters."""
    x = normalize_composition(composition)
    bm = torch.sum(x * pure_b, dim=-1)
    pure_alpha = pure_a / (pure_b * R * temperature[..., None])
    ge_rt = self.activity_model.excess_gibbs_rt(temperature, x)
    alpha = torch.sum(x * pure_alpha, dim=-1) - ge_rt / self.delta
    return alpha * bm * R * temperature, bm

PPR78Mixing

PPR78Mixing(
    group_fractions: Tensor,
    group_a: Tensor,
    group_b: Tensor,
    *,
    reference_temperature: float = 298.15,
    trainable: bool = False,
    parameter_set: str | None = None,
)

Bases: Module

Predictive PR78 group-contribution attraction mixing.

This is Eq. (5) of Jaubert and Mutelet (2004), doi:10.1016/j.fluid.2004.06.059. For component group fractions \(\alpha_{ik}\), the method evaluates

\[ k_{ij}(T) = \frac{ -\frac12\sum_{kl}\Delta\alpha_{ij,k}\Delta\alpha_{ij,l} A_{kl}(T_r/T)^{B_{kl}/A_{kl}-1} -(\sqrt{a_i}/b_i-\sqrt{a_j}/b_j)^2 }{ 2\sqrt{a_i a_j}/(b_i b_j) }. \]

group_a and group_b are symmetric pressure-valued group interaction matrices. Only their unique off-diagonal entries are stored, avoiding redundant degrees of freedom during fitting. The paper derives the correlation with the conventional linear covolume rule, which this class preserves.

Parameters:

Name	Type	Description	Default
group_fractions	Tensor	Normalized component-by-group fraction matrix.	required
group_a	Tensor	Symmetric pressure-valued universal interaction matrices in Pa.	required
group_b	Tensor	Symmetric pressure-valued universal interaction matrices in Pa.	required
reference_temperature	float	Positive correlation reference temperature in K.	298.15
trainable	bool	Register unique off-diagonal group interactions as trainable parameters.	False
parameter_set	str | None	Optional source parameter-set identifier.	None

Source code in src/torch_flash/mixing/rules.py

def __init__(
    self,
    group_fractions: Tensor,
    group_a: Tensor,
    group_b: Tensor,
    *,
    reference_temperature: float = 298.15,
    trainable: bool = False,
    parameter_set: str | None = None,
) -> None:
    super().__init__()
    if group_fractions.ndim != 2:
        raise ValueError("PPR78 group_fractions must be a component-by-group matrix")
    component_count, group_count = group_fractions.shape
    if component_count == 0 or group_count == 0:
        raise ValueError("PPR78 requires at least one component and one group")
    if group_a.shape != (group_count, group_count) or group_b.shape != group_a.shape:
        raise ValueError("PPR78 group A and B must be square matrices matching the groups")
    if not bool(
        torch.isfinite(group_fractions).all()
        & torch.isfinite(group_a).all()
        & torch.isfinite(group_b).all()
    ):
        raise ValueError("PPR78 fractions and group interactions must be finite")
    if bool((group_fractions < 0.0).any()):
        raise ValueError("PPR78 group fractions cannot be negative")
    row_sums = group_fractions.sum(dim=-1)
    if not torch.allclose(row_sums, torch.ones_like(row_sums)):
        raise ValueError("each PPR78 component group-fraction row must sum to one")
    if not torch.allclose(group_a, group_a.mT) or not torch.allclose(group_b, group_b.mT):
        raise ValueError("PPR78 group A and B matrices must be symmetric")
    if bool(torch.diagonal(group_a).count_nonzero()) or bool(
        torch.diagonal(group_b).count_nonzero()
    ):
        raise ValueError("PPR78 group A and B matrices must have zero diagonals")
    if not torch.isfinite(torch.tensor(reference_temperature)) or reference_temperature <= 0.0:
        raise ValueError("PPR78 reference_temperature must be finite and positive")
    undefined = (group_a == 0.0) & (group_b != 0.0)
    if bool(undefined.any()):
        raise ValueError("PPR78 group B must be zero wherever group A is zero")

    indices = torch.triu_indices(
        group_count,
        group_count,
        offset=1,
        device=group_a.device,
    )
    self.register_buffer("group_fractions", group_fractions.clone())
    self.register_buffer("pair_indices", indices)
    self.register_buffer(
        "reference_temperature",
        group_a.new_tensor(reference_temperature),
    )
    if trainable:
        self.raw_group_a = nn.Parameter(group_a[indices[0], indices[1]].clone())
        self.raw_group_b = nn.Parameter(group_b[indices[0], indices[1]].clone())
    else:
        self.register_buffer(
            "raw_group_a",
            group_a[indices[0], indices[1]].clone(),
        )
        self.register_buffer(
            "raw_group_b",
            group_b[indices[0], indices[1]].clone(),
        )
    self.parameter_set = parameter_set

ncomponents property

ncomponents: int

Number of component decompositions.

ngroups property

ngroups: int

Number of structural groups.

group_a property

group_a: Tensor

Return the symmetric :math:A_{kl} matrix in pascals.

group_b property

group_b: Tensor

Return the symmetric :math:B_{kl} matrix in pascals.

group_interaction_energy

group_interaction_energy(temperature: Tensor) -> Tensor

Return :math:A_{kl}(298.15/T)^{B_{kl}/A_{kl}-1} in pascals.

Source code in src/torch_flash/mixing/rules.py

def group_interaction_energy(self, temperature: Tensor) -> Tensor:
    """Return :math:`A_{kl}(298.15/T)^{B_{kl}/A_{kl}-1}` in pascals."""
    if bool((~torch.isfinite(temperature) | (temperature <= 0.0)).any()):
        raise ValueError("temperature must be finite and positive")
    group_a = self.group_a
    group_b = self.group_b
    is_zero = group_a == 0.0
    safe_a = torch.where(is_zero, torch.ones_like(group_a), group_a)
    exponent = group_b / safe_a - 1.0
    ratio = self.reference_temperature / temperature[..., None, None]
    value = group_a * ratio.pow(exponent)
    return torch.where(is_zero, torch.zeros_like(value), value)

kij

kij(
    temperature: Tensor, pure_a: Tensor, pure_b: Tensor
) -> Tensor

Evaluate the symmetric component BIP matrix.

pure_a and pure_b must be the PR78 pure-component parameters at temperature. This dependence is essential: PPR78 is not a standalone A + B/T law.

Source code in src/torch_flash/mixing/rules.py

def kij(
    self,
    temperature: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> Tensor:
    """Evaluate the symmetric component BIP matrix.

    ``pure_a`` and ``pure_b`` must be the PR78 pure-component parameters
    at ``temperature``. This dependence is essential: PPR78 is not a
    standalone ``A + B/T`` law.
    """
    if pure_a.shape[-1] != self.ncomponents or pure_b.shape[-1] != self.ncomponents:
        raise ValueError("PPR78 pure parameters must have one value per component")
    if not bool(
        torch.isfinite(pure_a).all()
        & torch.isfinite(pure_b).all()
        & (pure_a > 0.0).all()
        & (pure_b > 0.0).all()
    ):
        raise ValueError("PPR78 pure a and b parameters must be finite and positive")

    differences = self.group_fractions[:, None, :] - self.group_fractions[None, :, :]
    group_energy = self.group_interaction_energy(temperature)
    group_term = -0.5 * torch.einsum(
        "ijg,...gh,ijh->...ij",
        differences,
        group_energy,
        differences,
    )
    delta = torch.sqrt(pure_a) / pure_b
    pure_term = (delta[..., :, None] - delta[..., None, :]).square()
    denominator = (
        2.0
        * torch.sqrt(pure_a[..., :, None] * pure_a[..., None, :])
        / (pure_b[..., :, None] * pure_b[..., None, :])
    )
    result = (group_term - pure_term) / denominator
    return result - torch.diag_embed(torch.diagonal(result, dim1=-2, dim2=-1))

cross_a

cross_a(
    temperature: Tensor, pure_a: Tensor, pure_b: Tensor
) -> Tensor

Return PPR78 unlike attraction parameters.

Source code in src/torch_flash/mixing/rules.py

def cross_a(
    self,
    temperature: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> Tensor:
    """Return PPR78 unlike attraction parameters."""
    return torch.sqrt(pure_a[..., :, None] * pure_a[..., None, :]) * (
        1.0 - self.kij(temperature, pure_a, pure_b)
    )

partial_b

partial_b(composition: Tensor, pure_b: Tensor) -> Tensor

Return partial covolumes for the PPR78 linear covolume rule.

Source code in src/torch_flash/mixing/rules.py

def partial_b(self, composition: Tensor, pure_b: Tensor) -> Tensor:
    """Return partial covolumes for the PPR78 linear covolume rule."""
    x = normalize_composition(composition)
    return pure_b + torch.zeros_like(x)

forward

forward(
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]

Return PPR78 mixture a and linear mixture b.

Source code in src/torch_flash/mixing/rules.py

def forward(
    self,
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]:
    """Return PPR78 mixture ``a`` and linear mixture ``b``."""
    x = normalize_composition(composition)
    aij = self.cross_a(temperature, pure_a, pure_b)
    am = torch.einsum("...i,...ij,...j->...", x, aij, x)
    bm = torch.sum(x * pure_b, dim=-1)
    return am, bm

QuadraticMixing

QuadraticMixing(
    kij: Tensor,
    lij: Tensor | None = None,
    *,
    trainable: bool = False,
    trainable_lij: bool = False,
)

Bases: Module

van der Waals one-fluid mixing with kij and lij parameters.

The unlike attraction and covolume parameters are

aij = sqrt(ai*aj)*(1-kij) and bij = 0.5*(bi+bj)*(1-lij), respectively.

lij=0 recovers the conventional linear covolume rule exactly and uses a dedicated O(N) evaluation path. This convention is also exposed by ThermoPack's get_lij/set_lij cubic interface.

Parameters:

Name	Type	Description	Default
kij	Tensor	Symmetric dimensionless attraction-interaction matrix.	required
lij	Tensor | None	Optional symmetric dimensionless cross-covolume interactions.	None
trainable	bool	Register attraction or covolume interactions as trainable parameters.	False
trainable_lij	bool	Register attraction or covolume interactions as trainable parameters.	False

Source code in src/torch_flash/mixing/rules.py

def __init__(
    self,
    kij: Tensor,
    lij: Tensor | None = None,
    *,
    trainable: bool = False,
    trainable_lij: bool = False,
) -> None:
    super().__init__()
    if kij.ndim != 2 or kij.shape[0] != kij.shape[1]:
        raise ValueError("kij must be a square matrix")
    if not bool(torch.isfinite(kij).all()):
        raise ValueError("kij must contain only finite values")
    if not torch.allclose(kij, kij.mT):
        raise ValueError("kij must be symmetric")
    if lij is None:
        lij = torch.zeros_like(kij)
    elif lij.shape != kij.shape:
        raise ValueError("lij must have the same square shape as kij")
    elif not bool(torch.isfinite(lij).all()):
        raise ValueError("lij must contain only finite values")
    elif not torch.allclose(lij, lij.mT):
        raise ValueError("lij must be symmetric")
    if trainable:
        self.raw_kij = nn.Parameter(kij.clone())
    else:
        self.register_buffer("raw_kij", kij.clone())
    if trainable_lij:
        self.raw_lij = nn.Parameter(lij.clone())
    else:
        self.register_buffer("raw_lij", lij.clone())
    self._linear_covolume = not trainable_lij and not bool(torch.count_nonzero(lij))

kij property

kij: Tensor

Symmetrized interaction matrix with an exactly zero diagonal.

lij property

lij: Tensor

Symmetrized covolume interaction matrix with zero diagonal.

forward

forward(
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]

Return mixture a and b parameters.

Source code in src/torch_flash/mixing/rules.py

def forward(
    self,
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]:
    """Return mixture ``a`` and ``b`` parameters."""
    del temperature
    x = normalize_composition(composition)
    aij = self.cross_a(pure_a)
    am = torch.einsum("...i,...ij,...j->...", x, aij, x)
    if self._linear_covolume:
        bm = torch.sum(x * pure_b, dim=-1)
    else:
        bm = torch.einsum("...i,...ij,...j->...", x, self.cross_b(pure_b), x)
    return am, bm

cross_a

cross_a(pure_a: Tensor) -> Tensor

Return the matrix of unlike energy parameters aij.

Source code in src/torch_flash/mixing/rules.py

def cross_a(self, pure_a: Tensor) -> Tensor:
    """Return the matrix of unlike energy parameters ``aij``."""
    return torch.sqrt(pure_a[..., :, None] * pure_a[..., None, :]) * (1.0 - self.kij)

cross_b

cross_b(pure_b: Tensor) -> Tensor

Return unlike covolumes bij = (bi+bj)*(1-lij)/2.

Source code in src/torch_flash/mixing/rules.py

def cross_b(self, pure_b: Tensor) -> Tensor:
    """Return unlike covolumes ``bij = (bi+bj)*(1-lij)/2``."""
    return 0.5 * (pure_b[..., :, None] + pure_b[..., None, :]) * (1.0 - self.lij)

partial_b

partial_b(composition: Tensor, pure_b: Tensor) -> Tensor

Return partial molar covolumes d(n*b_mix)/d(n_i).

Source code in src/torch_flash/mixing/rules.py

def partial_b(self, composition: Tensor, pure_b: Tensor) -> Tensor:
    """Return partial molar covolumes ``d(n*b_mix)/d(n_i)``."""
    x = normalize_composition(composition)
    if self._linear_covolume:
        return pure_b + torch.zeros_like(x)
    bij = self.cross_b(pure_b)
    bm = torch.einsum("...i,...ij,...j->...", x, bij, x)
    return 2.0 * torch.einsum("...j,...ij->...i", x, bij) - bm[..., None]

TemperatureDependentQuadraticMixing

TemperatureDependentQuadraticMixing(
    a: Tensor,
    b: Tensor,
    lij: Tensor | None = None,
    *,
    trainable: bool = False,
    trainable_lij: bool = False,
)

Bases: Module

van der Waals mixing with kij(T) = Aij + Bij/T.

Aij is dimensionless and Bij is in kelvin. This form is used by Yan et al. for petroleum CPA calculations involving light hydrocarbons and water (Fluid Phase Equilibria 276, 2009, 75-85, doi:10.1016/j.fluid.2008.10.007).

Parameters:

Name	Type	Description	Default
a	Tensor	Symmetric dimensionless Aij matrix.	required
b	Tensor	Symmetric Bij matrix in K.	required
lij	Tensor | None	Optional symmetric dimensionless cross-covolume interactions.	None
trainable	bool	Register attraction or covolume interactions as trainable parameters.	False
trainable_lij	bool	Register attraction or covolume interactions as trainable parameters.	False

Source code in src/torch_flash/mixing/rules.py

def __init__(
    self,
    a: Tensor,
    b: Tensor,
    lij: Tensor | None = None,
    *,
    trainable: bool = False,
    trainable_lij: bool = False,
) -> None:
    super().__init__()
    if a.ndim != 2 or a.shape[0] != a.shape[1]:
        raise ValueError("temperature-dependent kij A must be a square matrix")
    if b.shape != a.shape:
        raise ValueError("temperature-dependent kij A and B matrices must have equal shapes")
    if not bool(torch.isfinite(a).all()) or not bool(torch.isfinite(b).all()):
        raise ValueError("temperature-dependent kij A and B matrices must be finite")
    if not torch.allclose(a, a.mT) or not torch.allclose(b, b.mT):
        raise ValueError("temperature-dependent kij A and B matrices must be symmetric")
    if lij is None:
        lij = torch.zeros_like(a)
    elif lij.shape != a.shape:
        raise ValueError("lij must have the same square shape as temperature-dependent kij")
    elif not bool(torch.isfinite(lij).all()):
        raise ValueError("lij must contain only finite values")
    elif not torch.allclose(lij, lij.mT):
        raise ValueError("lij must be symmetric")
    if trainable:
        self.raw_a = nn.Parameter(a.clone())
        self.raw_b = nn.Parameter(b.clone())
    else:
        self.register_buffer("raw_a", a.clone())
        self.register_buffer("raw_b", b.clone())
    if trainable_lij:
        self.raw_lij = nn.Parameter(lij.clone())
    else:
        self.register_buffer("raw_lij", lij.clone())
    self._linear_covolume = not trainable_lij and not bool(torch.count_nonzero(lij))

a property

a: Tensor

Return the symmetric dimensionless Aij matrix.

b property

b: Tensor

Return the symmetric Bij matrix in kelvin.

lij property

lij: Tensor

Return the symmetric dimensionless covolume interaction matrix.

kij

kij(temperature: Tensor) -> Tensor

Evaluate the interaction matrix at temperature in kelvin.

Source code in src/torch_flash/mixing/rules.py

def kij(self, temperature: Tensor) -> Tensor:
    """Evaluate the interaction matrix at temperature in kelvin."""
    if bool((~torch.isfinite(temperature) | (temperature <= 0.0)).any()):
        raise ValueError("temperature must be finite and positive")
    return self.a + self.b / temperature[..., None, None]

forward

forward(
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]

Return mixture a and b parameters.

Source code in src/torch_flash/mixing/rules.py

def forward(
    self,
    temperature: Tensor,
    composition: Tensor,
    pure_a: Tensor,
    pure_b: Tensor,
) -> tuple[Tensor, Tensor]:
    """Return mixture ``a`` and ``b`` parameters."""
    x = normalize_composition(composition)
    aij = torch.sqrt(pure_a[..., :, None] * pure_a[..., None, :]) * (
        1.0 - self.kij(temperature)
    )
    am = torch.einsum("...i,...ij,...j->...", x, aij, x)
    if self._linear_covolume:
        bm = torch.sum(x * pure_b, dim=-1)
    else:
        bm = torch.einsum("...i,...ij,...j->...", x, self.cross_b(pure_b), x)
    return am, bm

cross_a

cross_a(temperature: Tensor, pure_a: Tensor) -> Tensor

Return temperature-dependent unlike energy parameters.

Source code in src/torch_flash/mixing/rules.py

def cross_a(self, temperature: Tensor, pure_a: Tensor) -> Tensor:
    """Return temperature-dependent unlike energy parameters."""
    return torch.sqrt(pure_a[..., :, None] * pure_a[..., None, :]) * (
        1.0 - self.kij(temperature)
    )

cross_b

cross_b(pure_b: Tensor) -> Tensor

Return unlike covolumes bij = (bi+bj)*(1-lij)/2.

Source code in src/torch_flash/mixing/rules.py

def cross_b(self, pure_b: Tensor) -> Tensor:
    """Return unlike covolumes ``bij = (bi+bj)*(1-lij)/2``."""
    return 0.5 * (pure_b[..., :, None] + pure_b[..., None, :]) * (1.0 - self.lij)

partial_b

partial_b(composition: Tensor, pure_b: Tensor) -> Tensor

Return partial molar covolumes d(n*b_mix)/d(n_i).

Source code in src/torch_flash/mixing/rules.py

def partial_b(self, composition: Tensor, pure_b: Tensor) -> Tensor:
    """Return partial molar covolumes ``d(n*b_mix)/d(n_i)``."""
    x = normalize_composition(composition)
    if self._linear_covolume:
        return pure_b + torch.zeros_like(x)
    bij = self.cross_b(pure_b)
    bm = torch.einsum("...i,...ij,...j->...", x, bij, x)
    return 2.0 * torch.einsum("...j,...ij->...i", x, bij) - bm[..., None]