"""Cubic EoS mixing rules implementations module."""
from abc import ABC, abstractmethod
import numpy as np
from yaeos.core import GeModel
from yaeos.lib import yaeos_c
[docs]
class CubicMixRule(ABC):
"""Cubic mix rule abstract class."""
[docs]
@abstractmethod
def set_mixrule(self, ar_model_id: int) -> None:
"""Set mix rule abstract method.
Changes the default mixing rule of the given cubic EoS model.
Parameters
----------
ar_model_id : int
ID of the cubic EoS model
Raises
------
NotImplementedError
Abstract error, this method must be implemented in the subclass
"""
raise NotImplementedError
@abstractmethod
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
pass
@abstractmethod
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
pass
[docs]
class QMR(CubicMixRule):
"""Quadratic mixing rule.
Parameters
----------
kij : array_like
kij binary interaction parameters matrix
lij : array_like
lij binary interaction parameters matrix
Attributes
----------
kij : array_like
kij binary interaction parameters matrix
lij : array_like
lij binary interaction parameters matrix
Example
-------
.. code-block:: python
from yaeos import QMR, SoaveRedlichKwong
kij = [[0.0, 0.1], [0.1, 0.0]]
lij = [[0.0, 0.02], [0.02, 0.0]]
mixrule = QMR(kij, lij) # Quadratic mixing rule instance
tc = [305.32, 469.7] # critical temperature [K]
pc = [48.72, 33.7] # critical pressure [bar]
w = [0.0995, 0.152] # acentric factor
model = SoaveRedlichKwong(tc, pc, w, mixrule)
"""
name = "QMR"
def __init__(self, kij, lij) -> None:
self.kij = np.array(kij, order="F")
self.lij = np.array(lij, order="F")
[docs]
def set_mixrule(self, ar_model_id: int) -> None:
"""Set quadratic mix rule method.
Parameters
----------
ar_model_id : int
ID of the cubic EoS model
"""
yaeos_c.set_qmr(ar_model_id, self.kij, self.lij)
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
fcode = ""
kij_c = ""
lij_c = ""
for i in range(len(self.kij)):
kij_c += f"kij({i + 1}, :) = ["
lij_c += f"lij({i + 1}, :) = ["
for j in range(len(self.kij)):
if j < len(self.kij) - 1:
kij_c += f"{self.kij[i][j]}_pr, "
lij_c += f"{self.lij[i][j]}_pr, "
else:
kij_c += f"{self.kij[i][j]}_pr]\n"
lij_c += f"{self.lij[i][j]}_pr]\n"
fcode += kij_c + "\n"
fcode += lij_c + "\n"
fcode += "mixrule = QMR(k=kij, l=lij)\n\n"
return fcode
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
fcode = (
"type(QMR) :: mixrule"
"\n"
"real(pr) :: kij(nc, nc), lij(nc, nc)\n\n"
)
return fcode
[docs]
class QMRTD(CubicMixRule):
r"""Quadratic mixing rule, with temperature dependence.
..math::
k_{ij}(T) = k_{ij}^{\infty} + k_{ij}^0 \exp{\left(-T/T^{ref}\right)}
Parameters
----------
kij_0 : matrix_like
kij_0 binary interaction parameters matrix
kij_inf : matrix_like
kij_inf binary interaction parameters matrix
t_ref: matrix_like
Reference temperature
lij : matrix_like
lij binary interaction parameters matrix
Attributes
----------
kij_0 : matrix_like
kij_0 binary interaction parameters matrix
kij_inf : matrix_like
kij_inf binary interaction parameters matrix
t_ref: matrix_like
Reference temperature
lij : matrix_like
lij binary interaction parameters matrix
Example
-------
.. code-block:: python
from yaeos import QMRTD, SoaveRedlichKwong
kij_0 = [[0.0, 0.1], [0.1, 0.0]]
kij_inf = [[0.0, 0.1], [0.1, 0.0]]
Tref = [[0.0, 390], [390, 0.0]]
lij = [[0.0, 0.02], [0.02, 0.0]]
# Quadratic mixing rule instance
mixrule = QMRTD(kij_0, kij_inf, Tref, lij)
tc = [305.32, 469.7] # critical temperature [K]
pc = [48.72, 33.7] # critical pressure [bar]
w = [0.0995, 0.152] # acentric factor
model = SoaveRedlichKwong(tc, pc, w, mixrule)
"""
name = "QMRTD"
def __init__(self, kij_0, kij_inf, t_ref, lij) -> None:
self.kij_0 = np.array(kij_0, order="F")
self.kij_inf = np.array(kij_inf, order="F")
self.t_ref = np.array(t_ref, order="F")
self.lij = np.array(lij, order="F")
[docs]
def set_mixrule(self, ar_model_id: int) -> None:
"""Set mix rule method."""
yaeos_c.set_qmrtd(
ar_model_id,
kij_0=self.kij_0,
kij_inf=self.kij_inf,
t_star=self.t_ref,
lij=self.lij,
)
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
fcode = ""
kij0_c = ""
kijinf_c = ""
lij_c = ""
t_ref_c = ""
for i in range(len(self.kij_0)):
kij0_c += f"kij_0({i + 1}, :) = ["
kijinf_c += f"kij_inf({i + 1}, :) = ["
lij_c += f"lij({i + 1}, :) = ["
t_ref_c += f"t_ref({i + 1}, :) = ["
for j in range(len(self.kij_0)):
if j < len(self.kij_0) - 1:
kij0_c += f"{self.kij_0[i][j]}_pr, "
kijinf_c += f"{self.kij_inf[i][j]}_pr, "
lij_c += f"{self.lij[i][j]}_pr, "
t_ref_c += f"{self.t_ref[i][j]}_pr, "
else:
kij0_c += f"{self.kij_0[i][j]}_pr]\n"
kijinf_c += f"{self.kij_inf[i][j]}_pr]\n"
lij_c += f"{self.lij[i][j]}_pr]\n"
t_ref_c += f"{self.t_ref[i][j]}_pr]\n"
fcode += kij0_c + "\n"
fcode += kijinf_c + "\n"
fcode += lij_c + "\n"
fcode += t_ref_c + "\n"
fcode += "mixrule = QMRTD(k=kij_inf, k0=kij_0, Tref=t_ref, l=lij)\n\n"
return fcode
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
fcode = (
"type(QMRTD) :: mixrule"
"\n"
"real(pr) :: kij_inf(nc, nc), kij_0(nc, nc), lij(nc, nc)\n"
"real(pr) :: t_ref(nc, nc)\n\n"
)
return fcode
[docs]
class MHV(CubicMixRule):
"""Modified Huron-Vidal mixing rule.
Parameters
----------
ge : GeModel
Excess Gibbs energy model
q : float
q parameter. Use:
q = -0.594 for Soave-Redlich-Kwong
q = -0.53 for Peng-Robinson
q = -0.85 for Van der Waals
lij : array_like, optional
lij binary interaction parameters matrix, by default None
Attributes
----------
ge : GeModel
Excess Gibbs energy model
q : float
q parameter
lij : array_like
lij binary interaction parameters matrix
Example
-------
.. code-block:: python
from yaeos import MHV, SoaveRedlichKwong, NRTL
tc = [647.14, 513.92] # critical temperature [K]
pc = [220.64, 61.48] # critical pressure [bar]
w = [0.344, 0.649] # acentric factor
a = [[0, 3.458], [-0.801, 0]] # NRTL aij parameters
b = [[0, -586.1], [246.2, 0]] # NRTL bij parameters
c = [[0, 0.3], [0.3, 0]] # NRTL cij parameters
ge_model = NRTL(a, b, c)
mixrule = MHV(ge_model, q=-0.53)
model_mhv = PengRobinson76(tc, pc, w, mixrule)
"""
name = "MHV"
def __init__(self, ge: GeModel, q: float, lij=None) -> None:
nc = ge.size()
self.ge = ge
self.q = q
if lij is None:
self.lij = np.zeros((nc, nc), order="F")
else:
self.lij = lij
[docs]
def set_mixrule(self, ar_model_id: int) -> None:
"""Set modified Huron-Vidal mix rule method.
Parameters
----------
ar_model_id : int
ID of the cubic EoS model
"""
yaeos_c.set_mhv(
ar_id=ar_model_id, ge_id=self.ge.id, q=self.q, lij=self.lij
)
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
fcode = ""
fcode += self.ge._model_params_as_str()
lij_c = ""
for i in range(self.ge.size()):
lij_c += f"lij({i + 1}, :) = ["
for j in range(self.ge.size()):
if j < self.ge.size() - 1:
lij_c += f"{self.lij[i][j]}_pr, "
else:
lij_c += f"{self.lij[i][j]}_pr]\n"
fcode += lij_c + "\n"
fcode += (
f"mixrule = MHV(ge=ge_model, q={self.q}_pr, b=ar_model%b, lij=lij)\n\n"
)
return fcode
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
fcode = ""
fcode += self.ge._model_params_declaration_as_str()
# Just in case all possible replaces
fcode = fcode.replace(
f"integer, parameter :: nc={self.ge.size()}\n\n", ""
)
fcode = fcode.replace(
f"integer, parameter :: nc={self.ge.size()}\n", ""
)
fcode = fcode.replace(f"integer, parameter :: nc={self.ge.size()}", "")
# Mixrule setup
fcode += "type(MHV) :: mixrule" "\n"
fcode += "real(pr) :: lij(nc, nc)\n\n"
return fcode
[docs]
class HV(CubicMixRule):
"""Huron-Vidal mixing rule.
Parameters
----------
ge : GeModel
Excess Gibbs energy model
Attributes
----------
ge : GeModel
Excess Gibbs energy model
Example
-------
.. code-block:: python
from yaeos import HV, SoaveRedlichKwong, NRTL
tc = [647.14, 513.92] # critical temperature [K]
pc = [220.64, 61.48] # critical pressure [bar]
w = [0.344, 0.649] # acentric factor
a = [[0, 3.458], [-0.801, 0]] # NRTL aij parameters
b = [[0, -586.1], [246.2, 0]] # NRTL bij parameters
c = [[0, 0.3], [0.3, 0]] # NRTL cij parameters
ge_model = NRTL(a, b, c)
mixrule = HV(ge_model)
model_hv = PengRobinson76(tc, pc, w, mixrule)
"""
name = "HV"
def __init__(self, ge: GeModel) -> None:
self.ge = ge
[docs]
def set_mixrule(self, ar_model_id: int) -> None:
"""Set modified Huron-Vidal mix rule method.
Parameters
----------
ar_model_id : int
ID of the cubic EoS model
"""
yaeos_c.set_hv(ar_model_id, self.ge.id)
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
fcode = ""
fcode += self.ge._model_params_as_str()
fcode += "mixrule%ge = ge_model \n"
fcode += "mixrule%bi = ar_model%b \n"
fcode += "mixrule%del1 = ar_model%del1 \n"
return fcode
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
fcode = ""
fcode += self.ge._model_params_declaration_as_str()
# Just in case all possible replaces
fcode = fcode.replace(
f"integer, parameter :: nc={self.ge.size()}\n\n", ""
)
fcode = fcode.replace(
f"integer, parameter :: nc={self.ge.size()}\n", ""
)
fcode = fcode.replace(f"integer, parameter :: nc={self.ge.size()}", "")
fcode += "type(HV) :: mixrule" "\n"
return fcode
[docs]
class HVNRTL(CubicMixRule):
"""Huron-Vidal mixing rule using an NRTL model.
Huron Vidal mixing rule coupled with the excess Gibbs energy model defined
by Huron-Vidal. This model can use the classic VdW1f mixing rules kij
parameters when desired.
Parameters
----------
alpha : matrix_like
NRTL alpha parameters matrix
gji : matrix_like
NRTL gij parameters matrix
gjiT: matrix_like
NRTL gij temperature coefficient parameters matrix
use_kij : matrix_like
Boolean matrix indicating whether to use kij parameters
kij : matrix_like
kij binary interaction parameters matrix
Attributes
----------
alpha : matrix_like
NRTL alpha parameters matrix
gji : matrix_like
NRTL gij parameters matrix
gjiT: matrix_like
NRTL gij temperature coefficient parameters matrix
use_kij : matrix_like
Boolean matrix indicating whether to use kij parameters
kij : matrix_like
kij binary interaction parameters matrix
"""
name = "HVNRTL"
def __init__(self, alpha, gji, gjiT=None, use_kij=None, kij=None) -> None:
nc = np.shape(gji)[0]
self.alpha = alpha
self.gji = gji
if gjiT is None:
self.gjiT = np.zeros((nc, nc), order="F")
else:
self.gjiT = np.array(gjiT, order="F")
if use_kij is None and kij is None:
self.use_kij = np.zeros((nc, nc), dtype=bool, order="F")
self.kij = np.zeros((nc, nc), order="F")
else:
self.use_kij = np.array(use_kij, order="F")
self.kij = np.array(kij, order="F")
[docs]
def set_mixrule(self, ar_model_id: int) -> None:
"""Set modified Huron-Vidal mix rule method.
Parameters
----------
ar_model_id : int
ID of the cubic EoS model
"""
yaeos_c.set_hvnrtl(
ar_id=ar_model_id,
alpha=self.alpha,
gji0=self.gji,
gjit=self.gjiT,
use_kij=self.use_kij,
kij=self.kij,
)
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
nc = len(self.gji)
fcode = ""
gji0_c = ""
gjit_c = ""
use_kij_c = ""
kij_c = ""
alpha_c = ""
for i in range(nc):
gji0_c += f"gji0({i + 1}, :) = ["
gjit_c += f"gjiT({i + 1}, :) = ["
use_kij_c += f"use_kij({i + 1}, :) = ["
kij_c += f"kij({i + 1}, :) = ["
alpha_c += f"alphaij({i + 1}, :) = ["
for j in range(nc):
sep = ", " if j < nc - 1 else "]\n"
gji0_c += f"{self.gji[i][j]}_pr{sep}"
gjit_c += f"{self.gjiT[i][j]}_pr{sep}"
kij_c += f"{self.kij[i][j]}_pr{sep}"
alpha_c += f"{self.alpha[i][j]}_pr{sep}"
use_kij_c += (
".true." if self.use_kij[i][j] else ".false."
) + sep
fcode += alpha_c + "\n"
fcode += kij_c + "\n"
fcode += gji0_c + "\n"
fcode += gjit_c + "\n"
fcode += use_kij_c + "\n"
fcode += (
"mixrule = init_hvnrtl("
"b=ar_model%b, del1=ar_model%del1, "
"alpha=alphaij, gji0=gji0, gjiT=gjiT, "
"use_kij=use_kij, kij=kij)\n\n"
)
return fcode
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
fcode = (
"type(HV_NRTL) :: mixrule\n"
"\n"
"real(pr) :: gji0(nc,nc), gjiT(nc,nc), kij(nc,nc), alphaij(nc,nc)"
"\n"
"logical :: use_kij(nc,nc)\n\n"
)
return fcode
[docs]
@classmethod
def setup_from_kij(cls, model, temperatures, kij):
r"""
Generalizes the calculation of gji parameters and returns linear
correlation matrices gji0 and gjiT for any number of components.
.. warning::
Internally this function depends on the :math:`\delta_1` of the
model. For the case of the RKPR EoS where there are different
:math:`\delta_1` parameters for each component the method will
only use a mean value. Which will make the predictions further
away of the predictions being made with :math:`k_{ij}` compared
to using classic Cubic EoS. This method is intended only as a way
to ease the initialization of parameters on optimization methods.
Parameters
----------
model: yaeos.CubicEoS
Thermodynamic model intended to be used as reference model.
temperatures: array_like
Array of temperature values on which make the correlation
kij: float
kij value to replicate.
Returns
-------
yaeos.HVNRTL: Instance of a HVNRTL mixing rule with parameters
set to be close the parameters that would replicate the classic
QMR mixing rule.
"""
from yaeos.lib import yaeos_c
def lambd(d1=1+np.sqrt(2)):
f = d1 + 1
g = (d1 + 1)*d1 + d1 - 1
h = np.log((d1+1)**2 / 2)
L = f/g * h
return L
# 1. Setup basic parameters
nc = model.size()
ac, b, del1, del2, k = yaeos_c.get_ac_b_del1_del2_k(model.id, nc)
nt = len(temperatures)
ais = model.get_attractive_parameter(temperatures)["a"] # Shape: (nt, nc)
bis = model.get_repulsive_parameter() # Shape: (nc,)
lambd = lambd(np.mean(del1))
# Storage for all gji matrices across all temperatures
# Shape: (Temperature, Component_j, Component_i)
gjiss = np.zeros((nt, nc, nc))
# 2. Calculate gji at each Temperature
for t_idx, T in enumerate(temperatures):
ais_i = ais[t_idx]
# Calculate gii (pure component interaction) for this T
# gii = - (a_i / b_i) * lambda
gii = - (ais_i / bis) * lambd
gji_temp = np.zeros((nc, nc))
for i in range(nc):
for j in range(nc):
# Cross-interaction term logic
term1 = -2 * np.sqrt(bis[i] * bis[j]) / (bis[i] + bis[j])
term2 = np.sqrt(gii[i] * gii[j]) * (1 - kij)
# gji = (Interaction Term) - gii[i]
gji_temp[j, i] = (term1 * term2) - gii[i]
gjiss[t_idx] = gji_temp
# 3. Linear Correlation: gji = gji0 + gjiT * T
gji0 = np.zeros((nc, nc))
gjiT = np.zeros((nc, nc))
for i in range(nc):
for j in range(nc):
if i != j:
# Extract the vector of values for parameter g[j,i] over time
y = gjiss[:, j, i]
# Polyfit returns [slope, intercept]
slope, intercept = np.polyfit(temperatures, y, deg=1)
gjiT[j, i] = slope
gji0[j, i] = intercept
return cls(alpha=np.zeros((nc,nc)), gji=gji0, gjiT=gjiT)
[docs]
class sDDLC(CubicMixRule):
r"""segmented Density Dependent Local Composition mixing rule.
..math::
k_{ij}(T) = k_{ij}^{\infty} + k_{ij}^0 \exp{\left(-T/T^{ref}\right)}
Parameters
----------
q : array_like
Segment size for each component
kij_0 : matrix_like
kij_0 binary interaction parameters matrix
kij_inf : matrix_like
kij_inf binary interaction parameters matrix
t_ref: matrix_like
Reference temperature
lij : matrix_like
lij binary interaction parameters matrix
Attributes
----------
q : array_like
Segment size for each component
kij_0, optional : matrix_like
kij_0 binary interaction parameters matrix. Defaults to zero.
kij_inf, optional : matrix_like
kij_inf binary interaction parameters matrix. Defaults to zero.
t_ref, optional: matrix_like
Reference temperature. Defaults to zero.
lij, optional : matrix_like
lij binary interaction parameters matrix. Defaults to zero.
Example
-------
.. code-block:: python
from yaeos import sDDLC, SoaveRedlichKwong
qs = [1.16, 6.0]
kij_0 = [[0.0, 0.1], [0.1, 0.0]]
kij_inf = [[0.0, 0.1], [0.1, 0.0]]
Tref = [[0.0, 390], [390, 0.0]]
lij = [[0.0, 0.02], [0.02, 0.0]]
# sDDLC
mixrule = sDDLC(qs, kij_0, kij_inf, Tref, lij)
tc = [305.32, 469.7] # critical temperature [K]
pc = [48.72, 33.7] # critical pressure [bar]
w = [0.0995, 0.152] # acentric factor
model = SoaveRedlichKwong(tc, pc, w, mixrule)
"""
name = "sDDLC"
def __init__(
self, qs, kij_0=None, kij_inf=None, t_ref=None, lij=None
) -> None:
nc = len(qs)
self.qs = np.array(qs, order="F")
# kij constant term
if kij_inf is None:
self.kij_inf = np.zeros((nc, nc))
else:
self.kij_inf = np.array(kij_inf, order="F")
if t_ref is not None and kij_0 is None:
raise ValueError("kij_0 should be provided with t_ref")
if kij_0 is not None and t_ref is None:
raise ValueError("t_ref should be provided with kij_0")
# kij temperature dependance term
if kij_0 is None:
self.kij_0 = np.zeros((nc, nc))
self.t_ref = np.zeros((nc, nc))
else:
self.kij_0 = np.array(kij_0, order="F")
self.t_ref = np.array(t_ref, order="F")
if lij:
self.lij = np.array(lij, order="F")
else:
self.lij = np.zeros((nc, nc))
[docs]
def set_mixrule(self, ar_model_id: int) -> None:
"""Set mix rule method."""
yaeos_c.set_sddlc(
ar_model_id,
qs=self.qs,
kij_0=self.kij_0,
kij_inf=self.kij_inf,
t_star=self.t_ref,
lij=self.lij,
)
def _model_params_as_str(self) -> str:
"""Return the model parameters as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters. This string should be valid
Fortran code that assigns the model variables.
"""
fcode = ""
kij0_c = ""
kijinf_c = ""
lij_c = ""
t_ref_c = ""
qs_c = "qs = ["
for q in self.qs:
qs_c += f"{q:.5f}_pr,"
qs_c += "]"
for i in range(len(self.kij_0)):
kij0_c += f"kij_0({i + 1}, :) = ["
kijinf_c += f"kij_inf({i + 1}, :) = ["
lij_c += f"lij({i + 1}, :) = ["
t_ref_c += f"t_ref({i + 1}, :) = ["
for j in range(len(self.kij_0)):
if j < len(self.kij_0) - 1:
kij0_c += f"{self.kij_0[i][j]}_pr, "
kijinf_c += f"{self.kij_inf[i][j]}_pr, "
lij_c += f"{self.lij[i][j]}_pr, "
t_ref_c += f"{self.t_ref[i][j]}_pr, "
else:
kij0_c += f"{self.kij_0[i][j]}_pr]\n"
kijinf_c += f"{self.kij_inf[i][j]}_pr]\n"
lij_c += f"{self.lij[i][j]}_pr]\n"
t_ref_c += f"{self.t_ref[i][j]}_pr]\n"
fcode += qs_c + "\n"
fcode += kij0_c + "\n"
fcode += kijinf_c + "\n"
fcode += lij_c + "\n"
fcode += t_ref_c + "\n"
fcode += """
mixrule = sDDLC(q=qs, k=kij_inf, k0=kij_0, l=lij, tref=t_ref)
"""
return fcode
def _model_params_declaration_as_str(self) -> str:
"""Return the model parameters declaration as a string.
This method should be implemented by subclasses to return a string
representation of the model parameters declaration. This string should
be valid Fortran code that declares the model variables.
"""
fcode = (
"type(sDDLC) :: mixrule"
"\n"
"real(pr) :: kij_inf(nc, nc), kij_0(nc, nc), lij(nc, nc)\n"
"real(pr) :: t_ref(nc, nc), qs(nc)\n\n"
)
return fcode