polykin.copolymerization

PenultimateModel

Penultimate binary copolymerization model.

According to this model, the reactivity of a macroradical depends on the nature of the penultimate and terminal repeating units. A binary system is, thus, described by eight propagation reactions:

\[\begin{matrix} P^{\bullet}_{11} + M_1 \overset{k_{111}}{\rightarrow} P^{\bullet}_{11} \\ P^{\bullet}_{11} + M_2 \overset{k_{112}}{\rightarrow} P^{\bullet}_{12} \\ P^{\bullet}_{12} + M_1 \overset{k_{121}}{\rightarrow} P^{\bullet}_{21} \\ P^{\bullet}_{12} + M_2 \overset{k_{122}}{\rightarrow} P^{\bullet}_{22} \\ P^{\bullet}_{21} + M_1 \overset{k_{211}}{\rightarrow} P^{\bullet}_{11} \\ P^{\bullet}_{21} + M_2 \overset{k_{212}}{\rightarrow} P^{\bullet}_{12} \\ P^{\bullet}_{22} + M_1 \overset{k_{221}}{\rightarrow} P^{\bullet}_{21} \\ P^{\bullet}_{22} + M_2 \overset{k_{222}}{\rightarrow} P^{\bullet}_{22} \\ \end{matrix}\]

where \(k_{iii}\) are the homo-propagation rate coefficients and \(k_{ijk}\) are the cross-propagation coefficients. The six cross-propagation coefficients are specified via an equal number of reactivity ratios, which are divided in two categories. There are four monomer reactivity ratios, defined as \(r_{11}=k_{111}/k_{112}\), \(r_{12}=k_{122}/k_{121}\), \(r_{21}=k_{211}/k_{212}\) and \(r_{22}=k_{222}/k_{221}\). Additionally, there are two radical reactivity ratios defined as \(s_1=k_{211}/k_{111}\) and \(s_2=k_{122}/k_{222}\). The latter influence the average propagation rate coefficient, but have no effect on the copolymer composition.

PARAMETER	DESCRIPTION
r11	Monomer reactivity ratio, \(r_{11}=k_{111}/k_{112}\). TYPE: float
r12	Monomer reactivity ratio, \(r_{12}=k_{122}/k_{121}\). TYPE: float
r21	Monomer reactivity ratio, \(r_{21}=k_{211}/k_{212}\). TYPE: float
r22	Monomer reactivity ratio, \(r_{22}=k_{222}/k_{221}\). TYPE: float
s1	Radical reactivity ratio, \(s_1=k_{211}/k_{111}\). TYPE: float
s2	Radical reactivity ratio, \(s_2=k_{122}/k_{222}\). TYPE: float
k1	Homopropagation rate coefficient of M1, \(k_1 \equiv k_{111}\). TYPE: Arrhenius | None DEFAULT: None
k2	Homopropagation rate coefficient of M2, \(k_2 \equiv k_{222}\). TYPE: Arrhenius | None DEFAULT: None
M1	Name of M1. TYPE: str DEFAULT: 'M1'
M2	Name of M2. TYPE: str DEFAULT: 'M2'
name	Name. TYPE: str DEFAULT: ''

Source code in src/polykin/copolymerization/penultimate.py

class PenultimateModel(CopoModel):
    r"""Penultimate binary copolymerization model.

    According to this model, the reactivity of a macroradical depends on the
    nature of the _penultimate_ and _terminal_ repeating units. A binary system
    is, thus, described by eight propagation reactions:

    \begin{matrix}
    P^{\bullet}_{11} + M_1 \overset{k_{111}}{\rightarrow} P^{\bullet}_{11} \\
    P^{\bullet}_{11} + M_2 \overset{k_{112}}{\rightarrow} P^{\bullet}_{12} \\
    P^{\bullet}_{12} + M_1 \overset{k_{121}}{\rightarrow} P^{\bullet}_{21} \\
    P^{\bullet}_{12} + M_2 \overset{k_{122}}{\rightarrow} P^{\bullet}_{22} \\
    P^{\bullet}_{21} + M_1 \overset{k_{211}}{\rightarrow} P^{\bullet}_{11} \\
    P^{\bullet}_{21} + M_2 \overset{k_{212}}{\rightarrow} P^{\bullet}_{12} \\
    P^{\bullet}_{22} + M_1 \overset{k_{221}}{\rightarrow} P^{\bullet}_{21} \\
    P^{\bullet}_{22} + M_2 \overset{k_{222}}{\rightarrow} P^{\bullet}_{22} \\
    \end{matrix}

    where $k_{iii}$ are the homo-propagation rate coefficients and $k_{ijk}$
    are the cross-propagation coefficients. The six cross-propagation
    coefficients are specified via an equal number of reactivity ratios, which
    are divided in two categories. There are four monomer reactivity ratios,
    defined as $r_{11}=k_{111}/k_{112}$, $r_{12}=k_{122}/k_{121}$,
    $r_{21}=k_{211}/k_{212}$ and $r_{22}=k_{222}/k_{221}$. Additionally, there
    are two radical reactivity ratios defined as $s_1=k_{211}/k_{111}$ and
    $s_2=k_{122}/k_{222}$. The latter influence the average propagation rate
    coefficient, but have no effect on the copolymer composition.

    Parameters
    ----------
    r11 : float
        Monomer reactivity ratio, $r_{11}=k_{111}/k_{112}$.
    r12 : float
        Monomer reactivity ratio, $r_{12}=k_{122}/k_{121}$.
    r21 : float
        Monomer reactivity ratio, $r_{21}=k_{211}/k_{212}$.
    r22 : float
        Monomer reactivity ratio, $r_{22}=k_{222}/k_{221}$.
    s1 : float
        Radical reactivity ratio, $s_1=k_{211}/k_{111}$.
    s2 : float
        Radical reactivity ratio, $s_2=k_{122}/k_{222}$.
    k1 : Arrhenius | None
        Homopropagation rate coefficient of M1, $k_1 \equiv k_{111}$.
    k2 : Arrhenius | None
        Homopropagation rate coefficient of M2, $k_2 \equiv k_{222}$.
    M1 : str
        Name of M1.
    M2 : str
        Name of M2.
    name : str
        Name.
    """

    r11: float
    r12: float
    r21: float
    r22: float
    s1: float
    s2: float

    _pnames = ('r11', 'r12', 'r21', 'r22', 's1', 's2')

    def __init__(self,
                 r11: float,
                 r12: float,
                 r21: float,
                 r22: float,
                 s1: float,
                 s2: float,
                 k1: Optional[Arrhenius] = None,
                 k2: Optional[Arrhenius] = None,
                 M1: str = 'M1',
                 M2: str = 'M2',
                 name: str = ''
                 ) -> None:

        check_bounds(r11, 0., np.inf, 'r11')
        check_bounds(r12, 0., np.inf, 'r12')
        check_bounds(r21, 0., np.inf, 'r21')
        check_bounds(r22, 0., np.inf, 'r22')
        check_bounds(s1, 0., np.inf, 's1')
        check_bounds(s2, 0., np.inf, 's2')

        self.r11 = r11
        self.r12 = r12
        self.r21 = r21
        self.r22 = r22
        self.s1 = s1
        self.s2 = s2
        super().__init__(k1, k2, M1, M2, name)

    def ri(self,
            f1: Union[float, FloatArray]
           ) -> tuple[Union[float, FloatArray], Union[float, FloatArray]]:
        r"""Pseudoreactivity ratios.

        In the penultimate model, the pseudoreactivity ratios depend on the
        instantaneous comonomer composition according to:

        \begin{aligned}
           \bar{r}_1 &= r_{21}\frac{f_1 r_{11} + f_2}{f_1 r_{21} + f_2} \\
           \bar{r}_2 &= r_{12}\frac{f_2 r_{22} + f_1}{f_2 r_{12} + f_1}
        \end{aligned}

        where $r_{ij}$ are the monomer reactivity ratios.

        Parameters
        ----------
        f1 : float | FloatArray
            Molar fraction of M1.

        Returns
        -------
        tuple[FloatOrArray, FloatOrArray]
            Pseudoreactivity ratios, ($\bar{r}_1$, $\bar{r}_2$).
        """
        f2 = 1. - f1
        r11 = self.r11
        r12 = self.r12
        r21 = self.r21
        r22 = self.r22
        r1 = r21*(f1*r11 + f2)/(f1*r21 + f2)
        r2 = r12*(f2*r22 + f1)/(f2*r12 + f1)
        return (r1, r2)

    def kii(self,
            f1: Union[float, FloatArray],
            T: float,
            Tunit: Literal['C', 'K'] = 'K',
            ) -> tuple[Union[float, FloatArray], Union[float, FloatArray]]:
        r"""Pseudohomopropagation rate coefficients.

        In the penultimate model, the pseudohomopropagation rate coefficients
        depend on the instantaneous comonomer composition according to:

        \begin{aligned}
        \bar{k}_{11} &= k_{111}\frac{f_1 r_{11} + f_2}{f_1 r_{11} + f_2/s_1} \\
        \bar{k}_{22} &= k_{222}\frac{f_2 r_{22} + f_1}{f_2 r_{22} + f_1/s_2}
        \end{aligned}

        where $r_{ij}$ are the monomer reactivity ratios, $s_i$ are the radical
        reactivity ratios, and $k_{iii}$ are the homopropagation rate
        coefficients.

        Parameters
        ----------
        f1 : float | FloatArray
            Molar fraction of M1.
        T : float
            Temperature. Unit = `Tunit`.
        Tunit : Literal['C', 'K']
            Temperature unit.

        Returns
        -------
        tuple[float | FloatArray, float | FloatArray]
            Pseudohomopropagation rate coefficients,
            ($\bar{k}_{11}$, $\bar{k}_{22}$).
        """
        if self.k1 is None or self.k2 is None:
            raise ValueError(
                "To use this feature, `k1` and `k2` cannot be `None`.")
        f2 = 1. - f1
        r11 = self.r11
        r22 = self.r22
        s1 = self.s1
        s2 = self.s2
        k1 = self.k1(T, Tunit)
        k2 = self.k2(T, Tunit)
        k11 = k1*(f1*r11 + f2)/(f1*r11 + f2/s1)
        k22 = k2*(f2*r22 + f1)/(f2*r22 + f1/s2)
        return (k11, k22)

    def transitions(self,
                    f1: Union[float, FloatArrayLike]
                    ) -> dict[str, Union[float, FloatArray]]:
        r"""Calculate the instantaneous transition probabilities.

        For a binary system, the transition probabilities are given by:

        \begin{aligned}
            P_{iii} &= \frac{r_{ii} f_i}{r_{ii} f_i + (1 - f_i)} \\
            P_{jii} &= \frac{r_{ji} f_i}{r_{ji} f_i + (1 - f_i)} \\
            P_{iij} &= 1 -  P_{iii} \\
            P_{jij} &= 1 -  P_{jii}
        \end{aligned}

        where $i,j=1,2, i \neq j$.

        **References**

        * NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization
        process modeling, Wiley, 1996, p. 181.

        Parameters
        ----------
        f1 : float | FloatArrayLike
            Molar fraction of M1.

        Returns
        -------
        dict[str, float | FloatArray]
            Transition probabilities,
            {'111': $P_{111}$, '211': $P_{211}$, '121': $P_{121}$, ... }.
        """

        if isinstance(f1, (list, tuple)):
            f1 = np.array(f1, dtype=np.float64)
        check_bounds(f1, 0., 1., 'f1')

        f2 = 1. - f1
        r11 = self.r11
        r12 = self.r12
        r21 = self.r21
        r22 = self.r22

        P111 = r11*f1/(r11*f1 + f2)
        P211 = r21*f1/(r21*f1 + f2)
        P112 = 1. - P111
        P212 = 1. - P211

        P222 = r22*f2/(r22*f2 + f1)
        P122 = r12*f2/(r12*f2 + f1)
        P221 = 1. - P222
        P121 = 1. - P122

        result = {'111': P111,
                  '211': P211,
                  '112': P112,
                  '212': P212,
                  '222': P222,
                  '122': P122,
                  '221': P221,
                  '121': P121}

        return result

    def triads(self,
               f1: Union[float, FloatArrayLike],
               ) -> dict[str, Union[float, FloatArray]]:
        r"""Calculate the instantaneous triad fractions.

        For a binary system, the triad fractions are given by:

        \begin{aligned}
            F_{iii} &\propto  P_{jii} \frac{P_{iii}}{1 - P_{iii}} \\
            F_{iij} &\propto 2 P_{jii} \\
            F_{jij} &\propto 1 - P_{jii}
        \end{aligned}

        where $P_{ijk}$ is the transition probability
        $i \rightarrow j \rightarrow k$, which is a function of the monomer
        composition and the reactivity ratios.

        **References**

        * NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization
        process modeling, Wiley, 1996, p. 181.

        Parameters
        ----------
        f1 : float | FloatArrayLike
            Molar fraction of M1.

        Returns
        -------
        dict[str, float | FloatArray]
            Triad fractions,
            {'111': $F_{111}$, '112': $F_{112}$, '212': $F_{212}$, ... }.
        """

        P111, P211, _, _, P222, P122, _, _ = self.transitions(f1).values()

        F111 = P211*P111/(1. - P111 + eps)
        F112 = 2*P211
        F212 = 1. - P211
        Fsum = F111 + F112 + F212
        F111 /= Fsum
        F112 /= Fsum
        F212 /= Fsum

        F222 = P122*P222/(1. - P222 + eps)
        F221 = 2*P122
        F121 = 1. - P122
        Fsum = F222 + F221 + F121
        F222 /= Fsum
        F221 /= Fsum
        F121 /= Fsum

        result = {'111': F111,
                  '112': F112,
                  '212': F212,
                  '222': F222,
                  '221': F221,
                  '121': F121}

        return result

    def sequence(self,
                 f1: Union[float, FloatArrayLike],
                 k: Optional[Union[int, IntArrayLike]] = None,
                 ) -> dict[str, Union[float, FloatArray]]:
        r"""Calculate the instantaneous sequence length probability or the
        number-average sequence length.

        For a binary system, the probability of finding $k$ consecutive units
        of monomer $i$ in a chain is:

        $$ S_{i,k} = \begin{cases}
            1 - P_{jii} & \text{if } k = 1 \\
            P_{jii}(1 - P_{iii})P_{iii}^{k-2}& \text{if } k \ge 2
        \end{cases} $$

        and the corresponding number-average sequence length is:

        $$ \bar{S}_i = \sum_k k S_{i,k} = 1 + \frac{P_{jii}}{1 - P_{iii}} $$

        where $P_{ijk}$ is the transition probability
        $i \rightarrow j \rightarrow k$, which is a function of the monomer
        composition and the reactivity ratios.

        **References**

        * NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization
        process modeling, Wiley, 1996, p. 180.

        Parameters
        ----------
        f1 : float | FloatArrayLike
            Molar fraction of M1.
        k : int | IntArrayLike | None
            Sequence length, i.e., number of consecutive units in a chain.
            If `None`, the number-average sequence length will be computed.

        Returns
        -------
        dict[str, float | FloatArray]
            If `k is None`, the number-average sequence lengths,
            {'1': $\bar{S}_1$, '2': $\bar{S}_2$}. Otherwise, the
            sequence probabilities, {'1': $S_{1,k}$, '2': $S_{2,k}$}.
        """

        P111, P211, _, _, P222, P122, _, _ = self.transitions(f1).values()

        if k is None:
            S1 = 1. + P211/(1. - P111 + eps)
            S2 = 1. + P122/(1. - P222 + eps)
        else:
            if isinstance(k, (list, tuple)):
                k = np.array(k, dtype=np.int32)
            S1 = np.where(k == 1, 1. - P211, P211*(1. - P111)*P111**(k - 2))
            S2 = np.where(k == 1, 1. - P122, P122*(1. - P222)*P222**(k - 2))

        return {'1': S1, '2': S2}

F1

F1(
    f1: Union[float, FloatArrayLike]
) -> Union[float, FloatArray]

Calculate the instantaneous copolymer composition, \(F_1\).

The calculation is handled by inst_copolymer_binary.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArrayLike

RETURNS	DESCRIPTION
float | FloatArray	Instantaneous copolymer composition, \(F_1\).

Source code in src/polykin/copolymerization/terminal.py

def F1(self,
       f1: Union[float, FloatArrayLike]
       ) -> Union[float, FloatArray]:
    r"""Calculate the instantaneous copolymer composition, $F_1$.

    The calculation is handled by
    [`inst_copolymer_binary`](inst_copolymer_binary.md).

    Parameters
    ----------
    f1 : float | FloatArrayLike
        Molar fraction of M1.

    Returns
    -------
    float | FloatArray
        Instantaneous copolymer composition, $F_1$.
    """

    if isinstance(f1, (list, tuple)):
        f1 = np.array(f1, dtype=np.float64)
    check_bounds(f1, 0., 1., 'f1')
    return inst_copolymer_binary(f1, *self.ri(f1))

add_data

add_data(
    data: Union[CopoDataset, list[CopoDataset]]
) -> None

Add a copolymerization dataset for subsequent analysis.

PARAMETER	DESCRIPTION
data	Experimental dataset(s). TYPE: Union[CopoDataset, list[CopoDataset]]

Source code in src/polykin/copolymerization/terminal.py

def add_data(self,
             data: Union[CopoDataset, list[CopoDataset]]
             ) -> None:
    r"""Add a copolymerization dataset for subsequent analysis.

    Parameters
    ----------
    data : Union[CopoDataset, list[CopoDataset]]
        Experimental dataset(s).
    """
    if not isinstance(data, list):
        data = [data]

    valid_monomers = {self.M1, self.M2}
    for ds in data:
        ds_monomers = {ds.M1, ds.M2}
        if ds_monomers != valid_monomers:
            raise ValueError(
                f"Monomers defined in dataset `{ds.name}` are invalid: "
                f"{ds_monomers}!={valid_monomers}")
        if ds not in set(self.data):
            self.data.append(ds)
        else:
            print(f"Warning: duplicate dataset '{ds.name}' was skipped.")

azeotrope property

azeotrope: Optional[float]

Calculate the azeotrope composition.

RETURNS	DESCRIPTION
float | None	If an azeotrope exists, it returns its composition in terms of \(f_1\).

drift

drift(
    f10: Union[float, FloatVectorLike],
    x: Union[float, FloatVectorLike],
) -> FloatArray

Calculate drift of comonomer composition in a closed system for a given total monomer conversion.

In a closed binary system, the drift in monomer composition is given by the solution of the following differential equation:

\[ \frac{\textup{d} f_1}{\textup{d}x} = \frac{f_1 - F_1}{1 - x} \]

with initial condition \(f_1(0)=f_{1,0}\), where \(f_1\) and \(F_1\) are, respectively, the instantaneous comonomer and copolymer composition of M1, and \(x\) is the total molar monomer conversion.

PARAMETER	DESCRIPTION
f10	Initial molar fraction of M1, \(f_{1,0}=f_1(0)\). TYPE: float | FloatVectorLike
x	Value(s) of total monomer conversion values where the drift is to be evaluated. TYPE: float | FloatVectorLike

RETURNS	DESCRIPTION
FloatArray	Monomer fraction of M1 at a given conversion, \(f_1(x)\).

Source code in src/polykin/copolymerization/terminal.py

def drift(self,
          f10: Union[float, FloatVectorLike],
          x: Union[float, FloatVectorLike]
          ) -> FloatArray:
    r"""Calculate drift of comonomer composition in a closed system for a
    given total monomer conversion.

    In a closed binary system, the drift in monomer composition is given by
    the solution of the following differential equation:

    $$ \frac{\textup{d} f_1}{\textup{d}x} = \frac{f_1 - F_1}{1 - x} $$

    with initial condition $f_1(0)=f_{1,0}$, where $f_1$ and $F_1$ are,
    respectively, the instantaneous comonomer and copolymer composition of
    M1, and $x$ is the total molar monomer conversion.

    Parameters
    ----------
    f10 : float | FloatVectorLike
        Initial molar fraction of M1, $f_{1,0}=f_1(0)$.
    x : float | FloatVectorLike
        Value(s) of total monomer conversion values where the drift is to
        be evaluated.

    Returns
    -------
    FloatArray
        Monomer fraction of M1 at a given conversion, $f_1(x)$.
    """
    f10, x = convert_FloatOrVectorLike_to_FloatVector([f10, x], False)

    def df1dx(x_, f1_):
        F1_ = inst_copolymer_binary(f1_, *self.ri(f1_))
        return (f1_ - F1_)/(1. - x_ + eps)

    sol = solve_ivp(df1dx,
                    (0., max(x)),
                    f10,
                    t_eval=x,
                    method='LSODA',
                    vectorized=True,
                    atol=1e-4,
                    rtol=1e-4)
    if sol.success:
        result = sol.y
        result = np.maximum(0., result)
        if result.shape[0] == 1:
            result = result[0]
    else:
        result = np.empty_like(x)
        result[:] = np.nan
        print(sol.message)

    return result

kii

kii(
    f1: Union[float, FloatArray],
    T: float,
    Tunit: Literal["C", "K"] = "K",
) -> tuple[
    Union[float, FloatArray], Union[float, FloatArray]
]

Pseudohomopropagation rate coefficients.

In the penultimate model, the pseudohomopropagation rate coefficients depend on the instantaneous comonomer composition according to:

\[\begin{aligned} \bar{k}_{11} &= k_{111}\frac{f_1 r_{11} + f_2}{f_1 r_{11} + f_2/s_1} \\ \bar{k}_{22} &= k_{222}\frac{f_2 r_{22} + f_1}{f_2 r_{22} + f_1/s_2} \end{aligned}\]

where \(r_{ij}\) are the monomer reactivity ratios, \(s_i\) are the radical reactivity ratios, and \(k_{iii}\) are the homopropagation rate coefficients.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArray
T	Temperature. Unit = Tunit. TYPE: float
Tunit	Temperature unit. TYPE: Literal['C', 'K'] DEFAULT: 'K'

RETURNS	DESCRIPTION
tuple[float | FloatArray, float | FloatArray]	Pseudohomopropagation rate coefficients, (\(\bar{k}_{11}\), \(\bar{k}_{22}\)).

Source code in src/polykin/copolymerization/penultimate.py

def kii(self,
        f1: Union[float, FloatArray],
        T: float,
        Tunit: Literal['C', 'K'] = 'K',
        ) -> tuple[Union[float, FloatArray], Union[float, FloatArray]]:
    r"""Pseudohomopropagation rate coefficients.

    In the penultimate model, the pseudohomopropagation rate coefficients
    depend on the instantaneous comonomer composition according to:

    \begin{aligned}
    \bar{k}_{11} &= k_{111}\frac{f_1 r_{11} + f_2}{f_1 r_{11} + f_2/s_1} \\
    \bar{k}_{22} &= k_{222}\frac{f_2 r_{22} + f_1}{f_2 r_{22} + f_1/s_2}
    \end{aligned}

    where $r_{ij}$ are the monomer reactivity ratios, $s_i$ are the radical
    reactivity ratios, and $k_{iii}$ are the homopropagation rate
    coefficients.

    Parameters
    ----------
    f1 : float | FloatArray
        Molar fraction of M1.
    T : float
        Temperature. Unit = `Tunit`.
    Tunit : Literal['C', 'K']
        Temperature unit.

    Returns
    -------
    tuple[float | FloatArray, float | FloatArray]
        Pseudohomopropagation rate coefficients,
        ($\bar{k}_{11}$, $\bar{k}_{22}$).
    """
    if self.k1 is None or self.k2 is None:
        raise ValueError(
            "To use this feature, `k1` and `k2` cannot be `None`.")
    f2 = 1. - f1
    r11 = self.r11
    r22 = self.r22
    s1 = self.s1
    s2 = self.s2
    k1 = self.k1(T, Tunit)
    k2 = self.k2(T, Tunit)
    k11 = k1*(f1*r11 + f2)/(f1*r11 + f2/s1)
    k22 = k2*(f2*r22 + f1)/(f2*r22 + f1/s2)
    return (k11, k22)

kp

kp(
    f1: Union[float, FloatArrayLike],
    T: float,
    Tunit: Literal["C", "K"] = "K",
) -> Union[float, FloatArray]

Calculate the average propagation rate coefficient, \(\bar{k}_p\).

The calculation is handled by kp_average_binary.

Note

This feature requires the attributes k11 and k22 to be defined.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArrayLike
T	Temperature. Unit = Tunit. TYPE: float
Tunit	Temperature unit. TYPE: Literal['C', 'K'] DEFAULT: 'K'

RETURNS	DESCRIPTION
float | FloatArray	Average propagation rate coefficient. Unit = L/(mol·s)

Source code in src/polykin/copolymerization/terminal.py

def kp(self,
       f1: Union[float, FloatArrayLike],
       T: float,
       Tunit: Literal['C', 'K'] = 'K'
       ) -> Union[float, FloatArray]:
    r"""Calculate the average propagation rate coefficient, $\bar{k}_p$.

    The calculation is handled by
    [`kp_average_binary`](kp_average_binary.md).

    Note
    ----
    This feature requires the attributes `k11` and `k22` to be defined.

    Parameters
    ----------
    f1 : float | FloatArrayLike
        Molar fraction of M1.
    T : float
        Temperature. Unit = `Tunit`.
    Tunit : Literal['C', 'K']
        Temperature unit.

    Returns
    -------
    float | FloatArray
        Average propagation rate coefficient. Unit = L/(mol·s)
    """

    if isinstance(f1, (list, tuple)):
        f1 = np.array(f1, dtype=np.float64)
    check_bounds(f1, 0., 1., 'f1')
    return kp_average_binary(f1, *self.ri(f1), *self.kii(f1, T, Tunit))

plot

plot(
    kind: Literal["drift", "kp", "Mayo", "triads"],
    show: Literal["auto", "all", "data", "model"] = "auto",
    M: Literal[1, 2] = 1,
    f0: Optional[Union[float, FloatVectorLike]] = None,
    T: Optional[float] = None,
    Tunit: Literal["C", "K"] = "K",
    title: Optional[str] = None,
    axes: Optional[Axes] = None,
    return_objects: bool = False,
) -> Optional[tuple[Optional[Figure], Axes]]

Generate a plot of instantaneous copolymer composition, monomer composition drift, or average propagation rate coefficient.

PARAMETER	DESCRIPTION
kind	Kind of plot to be generated. TYPE: Literal['drift', 'kp', 'Mayo', 'triads']
show	What informatation is to be plotted. TYPE: Literal['auto', 'all', 'data', 'model'] DEFAULT: 'auto'
M	Index of the monomer to be used in input argument f0 and in output results. Specifically, if M=i, then f0 stands for \(f_i(0)\) and plots will be generated in terms of \(f_i\) and \(F_i\). TYPE: Literal[1, 2] DEFAULT: 1
f0	Initial monomer composition, \(f_i(0)\), as required for a monomer composition drift plot. TYPE: float | FloatVectorLike | None DEFAULT: None
T	Temperature. Unit = Tunit. TYPE: float | None DEFAULT: None
Tunit	Temperature unit. TYPE: Literal['C', 'K'] DEFAULT: 'K'
title	Title of plot. If None, a default title with the monomer names will be used. TYPE: str | None DEFAULT: None
axes	Matplotlib Axes object. TYPE: Axes | None DEFAULT: None
return_objects	If True, the Figure and Axes objects are returned (for saving or further manipulations). TYPE: bool DEFAULT: False

RETURNS	DESCRIPTION
tuple[Figure | None, Axes] | None	Figure and Axes objects if return_objects is True.

Source code in src/polykin/copolymerization/terminal.py

def plot(self,
         kind: Literal['drift', 'kp', 'Mayo', 'triads'],
         show: Literal['auto', 'all', 'data', 'model'] = 'auto',
         M: Literal[1, 2] = 1,
         f0: Optional[Union[float, FloatVectorLike]] = None,
         T: Optional[float] = None,
         Tunit: Literal['C', 'K'] = 'K',
         title: Optional[str] = None,
         axes: Optional[Axes] = None,
         return_objects: bool = False
         ) -> Optional[tuple[Optional[Figure], Axes]]:
    r"""Generate a plot of instantaneous copolymer composition, monomer
    composition drift, or average propagation rate coefficient.

    Parameters
    ----------
    kind : Literal['drift', 'kp', 'Mayo', 'triads']
        Kind of plot to be generated.
    show : Literal['auto', 'all', 'data', 'model']
        What informatation is to be plotted.
    M : Literal[1, 2]
        Index of the monomer to be used in input argument `f0` and in
        output results. Specifically, if `M=i`, then `f0` stands for
        $f_i(0)$ and plots will be generated in terms of $f_i$ and $F_i$.
    f0 : float | FloatVectorLike | None
        Initial monomer composition, $f_i(0)$, as required for a monomer
        composition drift plot.
    T : float | None
        Temperature. Unit = `Tunit`.
    Tunit : Literal['C', 'K']
        Temperature unit.
    title : str | None
        Title of plot. If `None`, a default title with the monomer names
        will be used.
    axes : Axes | None
        Matplotlib Axes object.
    return_objects : bool
        If `True`, the Figure and Axes objects are returned (for saving or
        further manipulations).

    Returns
    -------
    tuple[Figure | None, Axes] | None
        Figure and Axes objects if return_objects is `True`.
    """
    check_in_set(M, {1, 2}, 'M')
    check_in_set(kind, {'Mayo', 'kp', 'drift', 'triads'}, 'kind')
    check_in_set(show, {'auto', 'all', 'data', 'model'}, 'show')

    label_model = None
    if axes is None:
        fig, ax = plt.subplots()
        if title is None:
            titles = {'Mayo': "Mayo-Lewis diagram",
                      'drift': "Monomer composition drift",
                      'kp': "Average propagation coefficient",
                      'triads': "Triad fractions"}
            title = titles[kind] + f" {self.M1}(1)-{self.M2}(2)"
        if title:
            fig.suptitle(title)
    else:
        ax = axes
        fig = None
        if self.name:
            label_model = self.name

    unit_range = (0., 1.)
    npoints = 1000
    Mname = self.M1 if M == 1 else self.M2
    ndataseries = 0

    if show == 'auto':
        if not self.data:
            show = 'model'
        else:
            show = 'all'

    if kind == 'Mayo':

        ax.set_xlabel(fr"$f_{M}$")
        ax.set_ylabel(fr"$F_{M}$")
        ax.set_xlim(*unit_range)
        ax.set_ylim(*unit_range)

        ax.plot(unit_range, unit_range, color='black', linewidth=0.5)

        if show in {'model', 'all'}:
            x = np.linspace(*unit_range, npoints, dtype=np.float64)
            y = self.F1(x)
            if M == 2:
                x[:] = 1. - x
                y[:] = 1. - y  # type: ignore
            ax.plot(x, y, label=label_model)

        if show in {'data', 'all'}:
            for ds in self.data:
                if not isinstance(ds, MayoDataset):
                    continue
                x = ds.getvar('f', Mname)
                y = ds.getvar('F', Mname)
                ndataseries += 1
                ax.scatter(x, y,
                           label=ds.name if ds.name else None)

    elif kind == 'drift':

        ax.set_xlabel(r"Total molar monomer conversion, $x$")
        ax.set_ylabel(fr"$f_{M}$")
        ax.set_xlim(*unit_range)
        ax.set_ylim(*unit_range)

        if show == 'model':

            if f0 is None:
                raise ValueError("`f0` is required for a `drift` plot.")
            else:
                if isinstance(f0, (int, float)):
                    f0 = [f0]
                f0 = np.array(f0, dtype=np.float64)
                check_bounds(f0, *unit_range, 'f0')

            x = np.linspace(*unit_range, 1000)
            if M == 2:
                f0[:] = 1. - f0
            y = self.drift(f0, x)
            if M == 2:
                y[:] = 1. - y
            if y.ndim == 1:
                y = y[np.newaxis, :]
            for i in range(y.shape[0]):
                ax.plot(x, y[i, :], label=label_model)

        if show in {'data', 'all'}:

            for ds in self.data:
                if not isinstance(ds, DriftDataset):
                    continue
                x = ds.getvar('x')
                y = ds.getvar('f', Mname)
                ndataseries += 1
                ax.scatter(x, y,
                           label=ds.name if ds.name else None)

                if show == 'all':
                    x = np.linspace(*unit_range, npoints)
                    y = self.drift(ds.getvar('f', self.M1)[0], x)
                    if M == 2:
                        y[:] = 1. - y
                    ax.plot(x, y, label=label_model)

    elif kind == 'kp':

        ax.set_xlabel(fr"$f_{M}$")
        ax.set_ylabel(r"$\bar{k}_p$")
        ax.set_xlim(*unit_range)

        if show == 'model':

            if T is None:
                raise ValueError("`T` is required for a `kp` plot.")

            x = np.linspace(*unit_range, npoints)
            y = self.kp(x, T, Tunit)
            if M == 2:
                x[:] = 1. - x
            ax.plot(x, y, label=label_model)

        if show in {'data', 'all'}:

            for ds in self.data:
                if not isinstance(ds, kpDataset):
                    continue
                x = ds.getvar('f', Mname)
                y = ds.getvar('kp')
                ndataseries += 1
                ax.scatter(x, y,
                           label=ds.name if ds.name else None)
                if show == 'all':
                    x = np.linspace(*unit_range, npoints)
                    y = self.kp(x, ds.T, ds.Tunit)
                    if M == 2:
                        x[:] = 1. - x
                    ax.plot(x, y, label=label_model)

    ax.grid(True)

    if axes is not None or ndataseries:
        ax.legend(bbox_to_anchor=(1.05, 1.), loc="upper left")

    if return_objects:
        return (fig, ax)

ri

ri(
    f1: Union[float, FloatArray]
) -> tuple[
    Union[float, FloatArray], Union[float, FloatArray]
]

Pseudoreactivity ratios.

In the penultimate model, the pseudoreactivity ratios depend on the instantaneous comonomer composition according to:

\[\begin{aligned} \bar{r}_1 &= r_{21}\frac{f_1 r_{11} + f_2}{f_1 r_{21} + f_2} \\ \bar{r}_2 &= r_{12}\frac{f_2 r_{22} + f_1}{f_2 r_{12} + f_1} \end{aligned}\]

where \(r_{ij}\) are the monomer reactivity ratios.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArray

RETURNS	DESCRIPTION
tuple[FloatOrArray, FloatOrArray]	Pseudoreactivity ratios, (\(\bar{r}_1\), \(\bar{r}_2\)).

Source code in src/polykin/copolymerization/penultimate.py

def ri(self,
        f1: Union[float, FloatArray]
       ) -> tuple[Union[float, FloatArray], Union[float, FloatArray]]:
    r"""Pseudoreactivity ratios.

    In the penultimate model, the pseudoreactivity ratios depend on the
    instantaneous comonomer composition according to:

    \begin{aligned}
       \bar{r}_1 &= r_{21}\frac{f_1 r_{11} + f_2}{f_1 r_{21} + f_2} \\
       \bar{r}_2 &= r_{12}\frac{f_2 r_{22} + f_1}{f_2 r_{12} + f_1}
    \end{aligned}

    where $r_{ij}$ are the monomer reactivity ratios.

    Parameters
    ----------
    f1 : float | FloatArray
        Molar fraction of M1.

    Returns
    -------
    tuple[FloatOrArray, FloatOrArray]
        Pseudoreactivity ratios, ($\bar{r}_1$, $\bar{r}_2$).
    """
    f2 = 1. - f1
    r11 = self.r11
    r12 = self.r12
    r21 = self.r21
    r22 = self.r22
    r1 = r21*(f1*r11 + f2)/(f1*r21 + f2)
    r2 = r12*(f2*r22 + f1)/(f2*r12 + f1)
    return (r1, r2)

sequence

sequence(
    f1: Union[float, FloatArrayLike],
    k: Optional[Union[int, IntArrayLike]] = None,
) -> dict[str, Union[float, FloatArray]]

Calculate the instantaneous sequence length probability or the number-average sequence length.

For a binary system, the probability of finding \(k\) consecutive units of monomer \(i\) in a chain is:

\[ S_{i,k} = \begin{cases} 1 - P_{jii} & \text{if } k = 1 \\ P_{jii}(1 - P_{iii})P_{iii}^{k-2}& \text{if } k \ge 2 \end{cases} \]

and the corresponding number-average sequence length is:

\[ \bar{S}_i = \sum_k k S_{i,k} = 1 + \frac{P_{jii}}{1 - P_{iii}} \]

where \(P_{ijk}\) is the transition probability \(i \rightarrow j \rightarrow k\), which is a function of the monomer composition and the reactivity ratios.

References

• NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization process modeling, Wiley, 1996, p. 180.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArrayLike
k	Sequence length, i.e., number of consecutive units in a chain. If None, the number-average sequence length will be computed. TYPE: int | IntArrayLike | None DEFAULT: None

RETURNS	DESCRIPTION
dict[str, float | FloatArray]	If k is None, the number-average sequence lengths, {'1': \(\bar{S}_1\), '2': \(\bar{S}_2\)}. Otherwise, the sequence probabilities, {'1': \(S_{1,k}\), '2': \(S_{2,k}\)}.

Source code in src/polykin/copolymerization/penultimate.py

def sequence(self,
             f1: Union[float, FloatArrayLike],
             k: Optional[Union[int, IntArrayLike]] = None,
             ) -> dict[str, Union[float, FloatArray]]:
    r"""Calculate the instantaneous sequence length probability or the
    number-average sequence length.

    For a binary system, the probability of finding $k$ consecutive units
    of monomer $i$ in a chain is:

    $$ S_{i,k} = \begin{cases}
        1 - P_{jii} & \text{if } k = 1 \\
        P_{jii}(1 - P_{iii})P_{iii}^{k-2}& \text{if } k \ge 2
    \end{cases} $$

    and the corresponding number-average sequence length is:

    $$ \bar{S}_i = \sum_k k S_{i,k} = 1 + \frac{P_{jii}}{1 - P_{iii}} $$

    where $P_{ijk}$ is the transition probability
    $i \rightarrow j \rightarrow k$, which is a function of the monomer
    composition and the reactivity ratios.

    **References**

    * NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization
    process modeling, Wiley, 1996, p. 180.

    Parameters
    ----------
    f1 : float | FloatArrayLike
        Molar fraction of M1.
    k : int | IntArrayLike | None
        Sequence length, i.e., number of consecutive units in a chain.
        If `None`, the number-average sequence length will be computed.

    Returns
    -------
    dict[str, float | FloatArray]
        If `k is None`, the number-average sequence lengths,
        {'1': $\bar{S}_1$, '2': $\bar{S}_2$}. Otherwise, the
        sequence probabilities, {'1': $S_{1,k}$, '2': $S_{2,k}$}.
    """

    P111, P211, _, _, P222, P122, _, _ = self.transitions(f1).values()

    if k is None:
        S1 = 1. + P211/(1. - P111 + eps)
        S2 = 1. + P122/(1. - P222 + eps)
    else:
        if isinstance(k, (list, tuple)):
            k = np.array(k, dtype=np.int32)
        S1 = np.where(k == 1, 1. - P211, P211*(1. - P111)*P111**(k - 2))
        S2 = np.where(k == 1, 1. - P122, P122*(1. - P222)*P222**(k - 2))

    return {'1': S1, '2': S2}

transitions

transitions(
    f1: Union[float, FloatArrayLike]
) -> dict[str, Union[float, FloatArray]]

Calculate the instantaneous transition probabilities.

For a binary system, the transition probabilities are given by:

\[\begin{aligned} P_{iii} &= \frac{r_{ii} f_i}{r_{ii} f_i + (1 - f_i)} \\ P_{jii} &= \frac{r_{ji} f_i}{r_{ji} f_i + (1 - f_i)} \\ P_{iij} &= 1 - P_{iii} \\ P_{jij} &= 1 - P_{jii} \end{aligned}\]

where \(i,j=1,2, i \neq j\).

References

• NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization process modeling, Wiley, 1996, p. 181.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArrayLike

RETURNS	DESCRIPTION
dict[str, float | FloatArray]	Transition probabilities, {'111': \(P_{111}\), '211': \(P_{211}\), '121': \(P_{121}\), ... }.

Source code in src/polykin/copolymerization/penultimate.py

def transitions(self,
                f1: Union[float, FloatArrayLike]
                ) -> dict[str, Union[float, FloatArray]]:
    r"""Calculate the instantaneous transition probabilities.

    For a binary system, the transition probabilities are given by:

    \begin{aligned}
        P_{iii} &= \frac{r_{ii} f_i}{r_{ii} f_i + (1 - f_i)} \\
        P_{jii} &= \frac{r_{ji} f_i}{r_{ji} f_i + (1 - f_i)} \\
        P_{iij} &= 1 -  P_{iii} \\
        P_{jij} &= 1 -  P_{jii}
    \end{aligned}

    where $i,j=1,2, i \neq j$.

    **References**

    * NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization
    process modeling, Wiley, 1996, p. 181.

    Parameters
    ----------
    f1 : float | FloatArrayLike
        Molar fraction of M1.

    Returns
    -------
    dict[str, float | FloatArray]
        Transition probabilities,
        {'111': $P_{111}$, '211': $P_{211}$, '121': $P_{121}$, ... }.
    """

    if isinstance(f1, (list, tuple)):
        f1 = np.array(f1, dtype=np.float64)
    check_bounds(f1, 0., 1., 'f1')

    f2 = 1. - f1
    r11 = self.r11
    r12 = self.r12
    r21 = self.r21
    r22 = self.r22

    P111 = r11*f1/(r11*f1 + f2)
    P211 = r21*f1/(r21*f1 + f2)
    P112 = 1. - P111
    P212 = 1. - P211

    P222 = r22*f2/(r22*f2 + f1)
    P122 = r12*f2/(r12*f2 + f1)
    P221 = 1. - P222
    P121 = 1. - P122

    result = {'111': P111,
              '211': P211,
              '112': P112,
              '212': P212,
              '222': P222,
              '122': P122,
              '221': P221,
              '121': P121}

    return result

triads

triads(
    f1: Union[float, FloatArrayLike]
) -> dict[str, Union[float, FloatArray]]

Calculate the instantaneous triad fractions.

For a binary system, the triad fractions are given by:

\[\begin{aligned} F_{iii} &\propto P_{jii} \frac{P_{iii}}{1 - P_{iii}} \\ F_{iij} &\propto 2 P_{jii} \\ F_{jij} &\propto 1 - P_{jii} \end{aligned}\]

where \(P_{ijk}\) is the transition probability \(i \rightarrow j \rightarrow k\), which is a function of the monomer composition and the reactivity ratios.

References

• NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization process modeling, Wiley, 1996, p. 181.

PARAMETER	DESCRIPTION
f1	Molar fraction of M1. TYPE: float | FloatArrayLike

RETURNS	DESCRIPTION
dict[str, float | FloatArray]	Triad fractions, {'111': \(F_{111}\), '112': \(F_{112}\), '212': \(F_{212}\), ... }.

Source code in src/polykin/copolymerization/penultimate.py

def triads(self,
           f1: Union[float, FloatArrayLike],
           ) -> dict[str, Union[float, FloatArray]]:
    r"""Calculate the instantaneous triad fractions.

    For a binary system, the triad fractions are given by:

    \begin{aligned}
        F_{iii} &\propto  P_{jii} \frac{P_{iii}}{1 - P_{iii}} \\
        F_{iij} &\propto 2 P_{jii} \\
        F_{jij} &\propto 1 - P_{jii}
    \end{aligned}

    where $P_{ijk}$ is the transition probability
    $i \rightarrow j \rightarrow k$, which is a function of the monomer
    composition and the reactivity ratios.

    **References**

    * NA Dotson, R Galván, RL Laurence, and M Tirrel. Polymerization
    process modeling, Wiley, 1996, p. 181.

    Parameters
    ----------
    f1 : float | FloatArrayLike
        Molar fraction of M1.

    Returns
    -------
    dict[str, float | FloatArray]
        Triad fractions,
        {'111': $F_{111}$, '112': $F_{112}$, '212': $F_{212}$, ... }.
    """

    P111, P211, _, _, P222, P122, _, _ = self.transitions(f1).values()

    F111 = P211*P111/(1. - P111 + eps)
    F112 = 2*P211
    F212 = 1. - P211
    Fsum = F111 + F112 + F212
    F111 /= Fsum
    F112 /= Fsum
    F212 /= Fsum

    F222 = P122*P222/(1. - P222 + eps)
    F221 = 2*P122
    F121 = 1. - P122
    Fsum = F222 + F221 + F121
    F222 /= Fsum
    F221 /= Fsum
    F121 /= Fsum

    result = {'111': F111,
              '112': F112,
              '212': F212,
              '222': F222,
              '221': F221,
              '121': F121}

    return result