The Mathematics of Ethylene Oligomerisation and Polymerisation

Abstract

Linear α-olefins or LAOs are produced by the catalytic oligomerisation of ethylene on a multimillion ton scale annually. A range of LAOs is typically obtained with varying chain lengths which follow a distribution. Depending on the catalyst, various types of distributions have been identified, such as Schulz–Flory, Poisson, alternating and selective oligomerisations such as ethylene trimerisation to 1-hexene and tetramerisation to 1-octene. A comprehensive mathematical analysis for all oligomer distributions is presented, showing the relations between the various distributions and with ethylene polymerisation, as well as providing mechanistic insight into the underlying chemical processes.

1 Introduction

Linear alpha olefins (LAOs, or 1-alkenes) find extensive use as co-monomers in olefin polymerisation and as intermediates to detergents, lubricants, and plasticisers [1,2,3]. Global consumption of LAOs grew at an average annual rate of 5.6% from 2012 to 2016, and global capacity for LAO production has been estimated at 6.2 million tonnes per year (mtpy) for 2018 [4, 5]. LAOs are generally prepared via the oligomerisation of ethylene, as is shown below in Scheme 1 [6].

Scheme 1

Generic scheme for LAO synthesis via the oligomerisation of ethylene

LAO processes may be grouped into two types: wide-range and on-purpose [7]. Wide-range processes give a distribution of LAO products from C₄ to C₂₀₊, whilst on-purpose processes predominantly give C₄, C₆ or C₈ LAOs [7]. Examples of catalytic systems used in industry are given in Figs. 1 and 2 [6,7,8,9].

Fig. 1

Wide-range industrial processes giving LAO distributions from C₄ to C₂₀₊ (X = appropriate anion). Note: INEOS’s process gives a Poisson distribution of LAOs, whilst the SHOP and Gulfene processes give Schulz–Flory distributions of LAOs

Fig. 2

On-purpose industrial processes predominantly giving C₄, C₆ or C₈ LAOs (X = appropriate anion)

Three of the most common types of distributions from ethylene oligomerisation and polymerisation processes—namely Schulz–Flory, selective-LAO, and Poisson—are observed with the systems in Figs. 1 and 2. In addition, a number of alternating LAO distributions with chromium-based systems have been reported in the literature [11,12,13,14,15,16,17]. These four types of ethylene product distribution are depicted in Fig. 3.

Fig. 3

Clockwise from top-left: illustrations of Schulz–Flory, alternating, selective-LAO, and Poisson distributions

Following on from the ground-breaking work 80 years ago by Schulz and Flory [18, 19], novel mathematical approaches for analysing product distributions from ethylene oligomerisation and polymerisation processes are presented in this paper. These approaches provide straightforward, applicable methods for characterising key mechanistic features of underlying catalytic processes, and demonstrate how Schulz–Flory, alternating, selective-LAO, and Poisson distributions may be considered as parts of a connected product distribution landscape.

2 Modelling Steady-State Oligomerisation Processes

2.1 First-Order Oligomerisation Processes

In 1935, Schulz characterised the mathematics behind first-order oligomerisation processes [18]. He noted the direct proportionality between mol(n)% (mol percentage for an oligomer of n ethylene units) and the product of a sequence featuring a constant propagation probability α:

$$mol (n)\% = c_{1} \alpha^{n}$$

(1a)

Following a continuous function approximation, this implies (see SI.1.1 for derivation):

$$mol\left( n \right)\% = - \left( {\ln \alpha } \right)\alpha^{n}$$

(1b)

The relationship between wt(n)%, mol(n)% and n is:

$$wt\left( n \right)\% \propto mol\left( n \right)\% \left( {RFM {\text{C}}_{2} {\text{H}}_{4} } \right)n$$

(2a)

where (RFM C₂H₄) is the Relative Formula Mass of ethylene.

Using the expression for mol(n)% from Eq. 1a and condensing:

$$wt\left( n \right)\% = c_{2} n\alpha^{n}$$

(2b)

Again making use of a continuous function approximation (see SI.1.2 for derivation), this implies:

$$wt\left( n \right)\% = \left( {ln\alpha } \right)^{2} n\alpha^{n}$$

(2c)

The following year, in 1936, Flory presented precise derivations without making use of continuous function approximations [19]. Flory noted the dependence of mol(n), the amount of moles for each oligomer of n ethylene units, on the product of (n − 1) propagation steps and 1 elimination step along with R, the total mol of oligomers produced during the reaction time:

$$mol\left( n \right) = (1 - \alpha )\alpha^{n - 1} R$$

(3a)

Consequently:

$$mol\left( n \right)\% = \frac{{(1 - \alpha )\alpha^{n - 1} R}}{R} = (1 - \alpha )\alpha^{n - 1}$$

(3b)

It follows that:

$$wt\left( n \right) = n\left( {\text{RFM}\, \text{C}_{2} \text{H}_{4} } \right)mol\left( n \right) = n\left( {\text{RFM}\, \text{C}_{2} \text{H}_{4} } \right)R(1 - \alpha )\alpha^{n - 1}$$

(4a)

which gives (see SI.1.3 for derivation):

$$wt\left( n \right)\% = n\left( {1 - \alpha } \right)^{2} \alpha^{n - 1}$$

(4b)

Comparing Eq. 1b with Eq. 3b, it is evident that both Schulz’s equation and Flory’s equation express:

$$mol\left( n \right)\% = c\alpha^{n}$$

(5a)

The two Eqs. 1b and 3b are subtly different, as shown in Fig. 4, because $- \left( {ln\alpha } \right)$ is the proportionality constant for Schulz and $(1 - \alpha)/\alpha$ is the proportionality constant for Flory. However, in practice, this makes little difference because $\alpha$ is typically determined by taking natural logarithms to each side of Eq. 5a, plotting ln(mol(n) %) vs. n, performing a linear regression analysis to find the gradient lnα and then exponentiating the result as a power of e [20]. The value for $\alpha$ determined in this way will be the same in both cases.

Fig. 4

Oligomer distribution in mol(n) % vs. n using Eqs. 1b (Schulz) and 3b (Flory) ($\alpha$ = 0.70)

By comparing Eq. 2c with Eq. 4b, it may be seen that:

$$wt\left( n \right)\% = c^{\prime}n\alpha^{n}$$

(5b)

with $\left( {ln\alpha } \right)^{2}$ as the proportionality constant for Schulz and $\left( {1 - \alpha } \right)^{2} /\alpha$ as the proportionality constant for Flory (Table 1).

Table 1 Comparison of Schulz’s equations and Flory’s equations

In the next section, a new derivation for Eq. 3a will be presented which utilises recurrence relations, followed by the presentation of a new back-solving method to determine α and R. The analysis lays the groundwork for an analogous consideration of second-order oligomerisation processes.

2.1.1 Derivation of the Schulz–Flory Equation via Recurrence Relations

Catalytic oligomerisation systems which give Schulz–Flory LAO distributions propagate either via a Cossee chain growth mechanism or via a metallacyclic mechanism (Schemes 2 and 3). Both mechanisms may be generically represented by Scheme 4. Deuterium labelling studies have shown that many chromium systems operate via the metallacyclic pathway [21,22,23].

Scheme 2

Cossee-type chain growth mechanism for Schulz–Flory LAO distributions. Active species and monomer shown in black, eliminated products shown in blue

Scheme 3

Metallacyclic mechanism for Schulz–Flory LAO distributions. Active species and monomer shown in black, eliminated products shown in blue

Scheme 4

Generic scheme for Schulz–Flory LAO distributions. Active species shown in black, rates and eliminated products shown in blue, probabilities of propagation and elimination pathways shown in green

In Scheme 4:

• R is the total amount in mol of oligomers produced (from T₁ onwards) during the reaction time i.e. the rate at which intermediate species arrive at , in units of mol t⁻¹.

• α is the propagation probability from to etc.

  • So (1 − α) is the elimination probability from etc.

• $T_{n}$ is the total mol of a particular LAO with chain length (C₂)_n produced during the reaction time (i.e. the rate at which that particular LAO is eliminated, in units of mol t⁻¹), with T₁ as the starting point in the sequence.

Expressions for each LAO fraction may be found by inspection of Scheme 4, or by using a recurrence relation methodology. These two techniques are demonstrated in the centre and right columns of Table 2 respectively.

Table 2 Expressions for ${\text{T}}_{\text{n}}$ via inspection (centre) or via recurrence relation (right)

The aim is to express T_n in terms of α and R. From inspection of the central column of Table 2, it may be inferred that the expression is:

$$T_{n + 1} = \left( {1 - \alpha } \right)\alpha^{n} R \leftrightarrow T_{n} = \left( {1 - \alpha } \right)\alpha^{n - 1} R \leftrightarrow mol\left( n \right) = \left( {1 - \alpha } \right)\alpha^{n - 1} R$$

(6)

Equation 6 is exactly the unnormalised Schulz–Flory equation, as per Eq. 3a.

Via the recurrence relation method, the general expression is:

$$T_{n + 1} = \alpha T_{n} \leftrightarrow T_{n + 1} - \alpha T_{n} = 0$$

(7)

Following a derivation (see SI.1.4), the general solution for this expression is:

$$T_{n} = mol\left( n \right) = \frac{{\left( {1 - \alpha } \right)}}{\alpha }R\alpha^{n} = \left( {1 - \alpha } \right)\alpha^{n - 1} R$$

(8)

i.e. the unnormalised Schulz–Flory equation, as per Eqs. 3a and 6.

With constant values for α and R, it is implicitly assumed that:

• The sequence begins at a suitable T₁ (e.g. a longer 1-alkene).

• The oligomerisation proceeds under steady state conditions (i.e. negligible change in ethylene pressure, reaction temperature and quantity of catalytic species).

  • Under steady state conditions, the standing quantities of each catalytic species are constant. The sum of rates to any particular active species “box” in Scheme 4 are equal to the sum of rates from that “box”. This is equivalent to saying that the sum of probabilities for any particular active species “box” is 1 e.g. for $R = \alpha R + \left( {1 - \alpha } \right)R \leftrightarrow 1 = \alpha + \left( {1 - \alpha } \right)$, for $\alpha R = \alpha^{2} R + \alpha \left( {1 - \alpha } \right)R \leftrightarrow 1 = \alpha + \left( {1 - \alpha } \right)$, etc.

2.1.2 Back-Solving to Determine $\alpha$ and R, and Checking for the Best Fit

In practice, the T_n terms are known (i.e. the amount in mol of each oligomer produced during the reaction time), and the aim is to find α and R. See SI.1.5 for a discussion on the method used to determine α and R, as well as a checking method to determine that the best fit has been found.

2.1.3 Application of the Fitting Procedure to Experiment

Bis(benzimidazole)pyridine chromium(III) chloride, in combination with MAO co-catalyst, gives a Schulz–Flory distribution of LAOs via a metallacyclic mechanism (Fig. 5) [24]. Fitting α either via a ln(mol(n)%) vs. n plot or via the recurrence relation method gives a resulting value of α = 0.84.

Fig. 5

Schulz − Flory LAO distribution as mol(n) % vs n (units of ethylene) obtained with bis(benzimidazole)pyridine chromium(III) chloride (1 μmol) [16, 24]. Inset: $ln\left( {mol\left( n \right)\% } \right)$ vs. $n$ plot. Conditions: 4 bar of ethylene pressure, 30 °C, MAO (8 mmol), toluene (220 mL), 70 min. Activity = 5314 g mmol⁻¹ h⁻¹ bar⁻¹

It is worth noting that the determination of α via the recurrence relation method offers a subtle advantage over ln(mol(n)%) vs. n plots, in that accurate fits can only be obtained for distributions with reasonably constant α values as per Scheme 4. Linear regression analysis via the ln(mol(n)%) vs. n method may permit adequate fits for data sets displaying non-constant α behaviour, leading to Schulz–Flory equations which deviate from experimental plots.

2.1.4 Derivations of Further Experimentally-Relevant Metrics

For polymer distributions, M_n, M_w, and PDI are defined as follows:

$$M_{n} = \frac{{\mathop \sum \nolimits_{n = 1}^{n = \infty } wt\left( n \right)}}{{\mathop \sum \nolimits_{n = 1}^{n = \infty } mol\left( n \right)}}$$

(9a)

$$M_{w} = \mathop \sum \limits_{n = 1}^{n = \infty } wt\left( n \right)\% n\left( {\text{RFM}\,\text{C}_{2} \text{H}_{4} } \right)$$

(9b)

$$PDI = \frac{{M_{w} }}{{M_{n} }}$$

(9c)

These definitions can be applied to Schulz–Flory distributions. Using Eq. 3a, 4a and 4b the following expressions are obtained (see SI.1.6 for derivations):

$$M_{n} = \frac{{\left( {\text{RFM}\, \text{C}_{2} \text{H}_{4} } \right)}}{{\left( {1 - \alpha } \right)}}$$

(9d)

$$M_{w} = \frac{{\left( {\text{RFM} \text{C}_{2} \text{H}_{4} } \right)\left( {1 + \alpha } \right)}}{{\left( {1 - \alpha } \right)}}$$

(9e)

$$PDI = \left( {1 + \alpha } \right)$$

(9f)

Note that α values are between 0 and 1. Consequently, PDI values are between 1 and 2 for Schulz–Flory distributions.

2.2 Second-Order Oligomerisation Processes

In 2006, Tomov et al. discovered a highly active chromium-based oligomerisation system, with activities in excess of 100,000 g mmol⁻¹ h⁻¹ bar⁻¹, which gave an alternating distribution of 1-alkenes via a metallacyclic mechanism (Fig. 6) [15].

Fig. 6

Distribution of 1-alkene oligomers as mol(n)% vs n (units of ethylene) produced via N,N-bis((1H-benzimidazol-2-yl)methyl)-N-methylamine chromium(III) chloride/MAO [16, 17]. Conditions: 32 nmol Cr, 4 bar continuous ethylene pressure, 50 °C, MAO (7 mmol), toluene (200 mL), 1 h. Activity = 102,300 g mmol⁻¹ h⁻¹ bar⁻¹

The shape of this distribution cannot be described by first-order linear homogeneous recurrence relations with constant coefficients e.g. the Schulz–Flory equation. Variations on the first-order theme with orderly coefficients are likewise unable to accurately describe the experimentally-observed results. For instance, in a hypothetical alternating coefficients case, the quantities of the series of oligomers 1-hexene, 1-decene, 1-tetradecene etc. always decrease with n. Yet the distribution in Fig. 6 evidently displays some sort of periodicity.

All oligomerisations may be mathematically considered in terms of recurrence relations, so from an analytic point of view it is logical to proceed from a first-order treatment to a second-order treatment. Chemically speaking, this implies a mechanism involving both single and double coordination and insertion of ethylene. The double coordination and insertion of olefins has been proposed previously. Ystenes proposed a “trigger” mechanism whereby a coordinated olefin only inserts once a second olefin has coordinated to the metal [25]. Lemos and co-workers proposed a related mechanism based on kinetic studies of Ziegler–Natta catalysts [26], and a similar mechanism has been implicated in other oligomerisation processes [27, 28].

Building on the analysis from Sect. 2.1, a new “alternating equation” is derived here, which may be viewed as a direct analogue to the Schulz–Flory equation for second-order oligomerisation processes and which accurately describes alternating distributions as per Fig. 6. This is followed by the presentation of a new back-solving method for the determination of key mechanistic parameters.

2.2.1 Derivation of the Alternating Equation via Recurrence Relations

In Scheme 5:

• R₁ is the rate at which intermediate species arrive at , in units of mol t⁻¹.

• R₂ is the rate at which intermediate species originating from before arrive at , in units of mol t⁻¹.

  • So the total amount in mol of oligomers, R, produced (from T₁ onwards) during the reaction time is expressed by $R = R_{1} + R_{2}$

• α is the propagation probability from to etc.

• β is the propagation probability from to etc.

  • So (1 − α − β) is the probability of elimination from etc.

• T_n is the total mol of a particular 1-alkene with chain length (C₂)_n produced during the reaction time (i.e. the rate at which that particular 1-alkene is eliminated, in units of mol t⁻¹), with T₁ as the starting point in the sequence.

Scheme 5

Catalytic scheme for the production of alternating 1-alkene oligomer distributions. Active species shown in black, rates and eliminated products shown in blue, probabilities of propagation and elimination pathways shown in green

With reference to Scheme 5, Table 3 expresses each LAO via inspection (central column) or via recurrence relation (right column). The aim is to express T_n in terms of α, β, R₁ and R₂. In contrast to Sect. 2.1, obtaining a general expression via inspection is no longer feasible. The method of solving via a recurrence relation is much more straightforward, as will be demonstrated below.

Table 3 Expressions for ${\text{T}}_{\text{n}}$ via inspection (centre) or via recurrence relation (right)

The recurrence relation is:

$$T_{n + 2} = {\alpha }T_{n + 1} + {\beta }T_{n} \leftrightarrow T_{n + 2} - {\alpha}T_{n + 1} - {\beta }T_{n} = 0$$

(10)

Following a derivation (see SI.2.1), the general solution to this expression is: [16]

$$T_{n} = mol\left( n \right) = c_{1} r_{1}^{n} + c_{2} r_{2}^{n}$$

with:

$$r_{1} = \frac{{{\alpha } + \sqrt {{\alpha }^{2} + 4{\beta }} }}{2}$$

$$r_{2} = \frac{{{\alpha } - \sqrt {{\alpha}^{2} + 4{\beta }} }}{2}$$

$$c_{1} = \frac{{\left( {1 - {\alpha} - {\beta }} \right)\left( {R_{1} r_{1} + R_{2} } \right)}}{{r_{1} \left( {r_{1} - r_{2} } \right)}}$$

$$c_{2} = \frac{{\left( {1 - {\alpha} - {\beta }} \right)\left( {R_{1} r_{2} + R_{2} } \right)}}{{r_{2} \left( {r_{2} - r_{1} } \right)}}$$

(11)

i.e. the unnormalised alternating equation.

With constant values for α, β, R₁ and R₂, it is implicitly assumed that:

• The sequence begins at a suitable T₁ (e.g. a longer 1-alkene)

• The oligomerisation proceeds under steady state conditions (i.e. negligible change in ethylene pressure, reaction temperature and quantity of catalytic species)

• Under steady state conditions, the standing quantities of each catalytic species are constant. The sum of rates to any particular active species “box” are equal the sum of rates from that “box”. This is equivalent to saying that the sum of probabilities from any particular active species “box” is 1 e.g. for $R_{1} = \alpha R_{1} + \beta R_{1} + \left( {1 - \alpha - \beta } \right)R_{1} \leftrightarrow 1 = \alpha + \beta + \left( {1 - \alpha - \beta } \right)$, for $\alpha R_{1} + R_{2} = \alpha \left( {\alpha R_{1} + R_{2} } \right) + \beta \left( {\alpha R_{1} + R_{2} } \right) + \left( {1 - \alpha - \beta } \right)\left( {\alpha R_{1} + R_{2} } \right) \leftrightarrow \alpha R_{1} = \alpha^{2} R_{1} + \beta \alpha R_{1} + \left( {1 - \alpha - \beta } \right)\alpha R_{1} \leftrightarrow 1 = \alpha + \beta + \left( {1 - \alpha - \beta } \right)$, etc.

To convert the unnormalised alternating equation, Eq. 11, into the normalised form:

$$mol\left( n \right)\% = \frac{mol\left( n \right)}{{\mathop \sum \nolimits_{n = 1}^{n = \infty } mol\left( n \right)}} = \frac{{c_{1} r_{1}^{n} + c_{2} r_{2}^{n} }}{{R_{1} + R_{2} }}$$

(12)

or equivalently (see SI.2.2):

$$mol\left( n \right)\% = \frac{{c_{1} r_{1}^{n} + c_{2} r_{2}^{n} }}{{\frac{{c_{1} r_{1} }}{{1 - r_{1} }} + \frac{{c_{2} r_{2} }}{{1 - r_{2} }}}}$$

(13)

Converting Eq. 11 into an unnormalised wt(n) equivalent:

$$wt\left( n \right) = n\left( {RFM C_{2} H_{4} } \right)mol\left( n \right) = n\left( {RFM C_{2} H_{4} } \right)(c_{1} r_{1}^{n} + c_{2} r_{2}^{n} )$$

(14)

which gives (see SI.2.3):

$$wt\left( n \right)\% = \frac{{n(c_{1} r_{1}^{n} + c_{2} r_{2}^{n} )}}{{\frac{{c_{1} r_{1} }}{{\left( {1 - r_{1} } \right)^{2} }} + \frac{{c_{2} r_{2} }}{{\left( {1 - r_{2} } \right)^{2} }}}}$$

(15)

2.2.2 Back-Solving to Determine α, β, R₁, and R₂, and Checking for the Best Fit

In practice, the T_n terms are known (i.e. the amount in mol of each oligomer produced during the reaction time), and the aim is to find α, β, R₁, and R₂. See SI.2.4 for a discussion on the method used to determine α and R, as well as a checking method to determine that the best fit has been found.

2.2.3 Application of the Fitting Procedure to Experiment

See Sect. 2.4 for examples of fittings using the described recurrence relation method with N,N-bis((1H-benzimidazol-2-yl)methyl)-N-methylamine chromium(III) chloride/MMAO-12 systems.

2.2.4 Derivations of Further Experimentally-Relevant Metrics

As per Schulz–Flory distributions, the definitions of M_n, M_w and PDI in Eqs. 9a, 9b, 9c can be applied to alternating distributions. Using Eqs. 11, 14 and 15, and Eqs. 12ii and 14iii in SI, the following expressions are obtained (see SI.2.5 for derivations):

$$M_{n} = \frac{{\left( {RFM C_{2} H_{4} } \right)}}{{\frac{{c_{1} r_{1} }}{{1 - r_{1} }} + \frac{{c_{2} r_{2} }}{{1 - r_{2} }}}}\left[ {\frac{{c_{1} r_{1} }}{{\left( {1 - r_{1} } \right)^{2} }} + \frac{{c_{2} r_{2} }}{{\left( {1 - r_{2} } \right)^{2} }}} \right]$$

(16a)

$$M_{w} = \frac{{\left( {RFM C_{2} H_{4} } \right)}}{{\frac{{c_{1} r_{1} }}{{\left( {1 - r_{1} } \right)^{2} }} + \frac{{c_{2} r_{2} }}{{\left( {1 - r_{2} } \right)^{2} }}}}\left[ {\frac{{c_{1} r_{1} \left( {1 + r_{1} } \right)}}{{\left( {1 - r_{1} } \right)^{3} }} + \frac{{c_{2} r_{2} \left( {1 + r_{2} } \right)}}{{\left( {1 - r_{2} } \right)^{3} }}} \right]$$

(16b)

$$PDI = \frac{py + qx}{{\left[ {py^{2} + qx^{2} } \right]^{2} }}\left[ {py^{3} \left( {2 - x} \right) + qx^{3} \left( {2 - y} \right)} \right]$$

with:

$$p = c_{1} r_{1}$$

$$q = c_{2} r_{2}$$

$$x = 1 - r_{1}$$

$$y = 1 - r_{2}$$

(16c)

Note that with $q = 0, y = 1, \text{and}\,r_{1} = \alpha$, PDI matches the expression previously obtained with Schulz–Flory distributions (see Eq. 9f):

$$PDI = \frac{p}{{p^{2} }}p\left( {2 - x} \right) = 2 - x = 1 + r_{1} = 1 + \alpha$$

As with Schulz–Flory distributions, PDI values for alternating distributions lie between 1 and 2. This may be demonstrated by examining four limiting cases (see SI.2.6).

2.3 Connecting First-Order and Second-Order Oligomerisation Processes

By inspection of Eq. 11, two cases may be identified where the second-order general solution collapses to a first-order solution as per Eq. 10. In the first case, if β = 0, then from Eq. 10iv $r_{1} = \alpha$ and from Eq. 10v $r_{2} = 0$. With $r_{2} = 0$, Eq. 10xiii gives $R_{2} = 0$. Then, using Eqs. 10x and 10xi (in SI.2.1 and SI.2.4):

$$c_{1} r_{1}^{n} = \frac{{\left( {1 - \alpha } \right)R_{1} r_{1} }}{{r_{1} \left( {r_{1} - r_{2} } \right)}}r_{1}^{n} = \frac{{\left( {1 - \alpha } \right)}}{{\left( {r_{1} - r_{2} } \right)}}r_{1}^{n} R_{1} = \frac{{\left( {1 - \alpha } \right)}}{{\left( {\alpha - 0} \right)}}\alpha^{n} R_{1} = \left( {1 - \alpha } \right)\alpha^{n - 1} R_{1}$$

(17a)

$$c_{2} r_{2}^{n} = \frac{{\left( {1 - \alpha } \right)R_{1} r_{2} }}{{r_{2} \left( {r_{2} - r_{1} } \right)}}r_{2}^{n} = \frac{{\left( {1 - \alpha } \right)}}{{\left( {r_{2} - r_{1} } \right)}}r_{2}^{n} R_{1} = \frac{{\left( {1 - \alpha } \right)}}{{\left( {0 - \alpha } \right)}}0^{n} R_{1} = 0$$

(17b)

Since $R_{2} = 0$, then $R_{1} = R$. Using Eq. 11, the general solution is therefore:

$$T_{n} = \left( {1 - \alpha } \right)\alpha^{n - 1} R$$

(18)

i.e. the unnormalised Schulz–Flory equation.

In the second, more general case, if $c_{2} = 0$ but $\beta$ is not necessarily 0, then with Eqs. 10xii and 10xiii in Supplementary Information (SI.2.4):

$$R_{1} = \frac{{r_{1} c_{1} }}{{\left( {1 - \alpha - \beta } \right)}}$$

$$R_{2} = \frac{{ - r_{1} r_{2} c_{1} }}{{\left( {1 - \alpha - \beta } \right)}}$$

$$\to R_{2} = - r_{2} R_{1} \leftrightarrow \frac{{R_{2} }}{{R_{1} }} = - r_{2}$$

(19)

Using Eq. 10v (Supplementary Information SI.2.1), Eq. 19 becomes:

$$\frac{{R_{2} }}{{R_{1} }} = \frac{{ - \alpha + \sqrt {\alpha^{2} + 4\beta } }}{2}$$

(20)

Using Eq. 19 with Eqs. 10x and 10xi (in SI.2.1):

$$c_{1} r_{1}^{n} = \frac{{\left( {1 - \alpha - \beta } \right)\left( {R_{1} r_{1} - R_{1} r_{2} } \right)}}{{r_{1} \left( {r_{1} - r_{2} } \right)}}r_{1}^{n} = \frac{{\left( {1 - \alpha - \beta } \right)\left( {r_{1} - r_{2} } \right)}}{{\left( {r_{1} - r_{2} } \right)}}r_{1}^{n - 1} R_{1} = \left( {1 - \alpha - \beta } \right)r_{1}^{n - 1} R_{1}$$

(21a)

$$c_{2} r_{2}^{n} = \frac{{\left( {1 - \alpha - \beta } \right)\left( {R_{1} r_{2} - R_{1} r_{2} } \right)}}{{r_{2} \left( {r_{2} - r_{1} } \right)}}r_{2}^{n} = \frac{{\left( {1 - \alpha - \beta } \right)\left( {r_{2} - r_{2} } \right)}}{{\left( {r_{2} - r_{1} } \right)}}r_{2}^{n - 1} R_{1} = \frac{{\left( {1 - \alpha - \beta } \right)0}}{{\left( {r_{2} - r_{1} } \right)}}r_{2}^{n - 1} R_{1} = 0$$

(21b)

From Eq. 11, the general solution is therefore:

$$T_{n} = \left( {1 - \alpha - \beta } \right)r_{1}^{n - 1} R_{1}$$

(22)

which gives the unnormalised Schulz–Flory equation if β = 0. The second case may be summarised by stating that given particular values for α and β, R₂/R₁ may be tuned to generate first-order behaviour.

In this sub-section, the mathematics of the second case is illustrated, followed by a method to describe the degree of similarity between second-order recurrence relations and first-order ones. This analysis has important implications for the fitting of more “gently” alternating oligomer distributions.

2.3.1 The First-Order Surface in R₂/R₁, α and β Space

Equation 20 features three variables, R₂/R₁, $\alpha$ and β, which are physically interpretable as per Scheme 5. Conveniently, then, it is possible to plot in 3D space the combinations of $R_{2} /R_{1}$, $\alpha$ and β values which give rise to first-order behaviour. Furthermore, since r₂ does not feature in Eq. 22 and since β is not necessarily 0, Eq. 10iv (Supplementary Information SI.2.1) now operates as a parametric equation where different combinations of α and β may give the same result for r₁. A visualisation of the first-order surface in R₂/R₁, $\alpha$ and $\beta$ space, constructed by lines along which r₁ is constant, is shown in Fig. 7 and 8. Note that first-order distributions following Eq. 18 are restricted to the alpha axis with β = 0 and R₂/R₁ = 0.

Fig. 7

Visualisation of the first-order surface in $R_{2} /R_{1}$, $\alpha$ and $\beta$ space along the $\alpha$ axis (left) and along the $\beta$ axis (right). Colours correspond to particular values of $r_{1}$ (key centre)

Fig. 8

Visualisation of the first-order surface in $R_{2} /R_{1}$, $\alpha$ and $\beta$ space along the $R_{2} /R_{1}$ axis. Optically identical first-order distributions with $r_{1} = 0.50$ illustrated

2.3.2 Degree of Alternation

By sight, it is evident that second-order distributions can possess varying degrees of alternating character. To arrive at a definition for the degree of alternation in a second-order distribution, we need to consider the pictorial representation of Eq. 11.

Figure 9 illustrates how an alternating distribution is mathematically constructed by adding a constantly decreasing first-order recurrence relation to an alternating first-order recurrence relation (recall that with 0 ≤ α < 1, 0 ≤ β < 1, and 0 < (1 − α − β) < 1, then 0 < r₁ < 1, $- 1 < r_{2} < 0, {\text{and also}}\left| {r_{2} } \right| \le \left| {r_{1} } \right|$). The degree of alternation is defined as:

$$Deg.\,of\,Alt. = \frac{{\mathop \sum \nolimits_{n = 1}^{n = \infty } \left| {c_{2} } \right|\left| {r_{2} } \right|^{n} }}{{\mathop \sum \nolimits_{n = 1}^{n = \infty } c_{1} r_{1}^{n} }} = \frac{{\left| {c_{2} } \right|}}{{c_{1} }}\frac{{\mathop \sum \nolimits_{n = 1}^{n = \infty } \left( { - r_{2} } \right)^{n} }}{{\mathop \sum \nolimits_{n = 1}^{n = \infty } r_{1}^{n} }}$$

(23a)

Fig. 9

Pictorial representation of Eq. 11

Absolute values are taken in the numerator to avoid cancellation of terms on summation (refer to Fig. 9). 0 < c₁ (see Eq. 10viii in SuppIementary Information SI.2.1) so $0 < \mathop \sum \nolimits_{n = 1}^{n = \infty } c_{1} r_{1}^{n}$. Given that $- 1 < r_{2} < 0$, then $\left| {r_{2} } \right|^{n} = \left( { - r_{2} } \right)^{n}$.

Evaluation of Eq. 23a (in SI.3.1) gives the following result:

$$Deg.\;of\;Alt. = \frac{{ + \sqrt {\left( {r_{2} + \frac{{R_{2} }}{{R_{1} }}} \right)^{2} } \left( {1 - r_{1} } \right)}}{{\left( {r_{1} + \frac{{R_{2} }}{{R_{1} }}} \right)\left( {1 + r_{2} } \right)}}$$

(23b)

Since $0 < r_{1} < 1$, $- 1 < r_{2} < 0, {\text{and }}\left| {r_{2} } \right| \le \left| {r_{1} } \right|$, then $0 \le Deg.\;of\;Alt. < 1$ (i.e. $Deg.\;of\;Alt.$ may be considered as a percentage).

As $r_{1}$ and $r_{2}$ may be expressed in terms of α and β, Eq. 23b gives the degree of alternation as a function of R₂/R₁,$\alpha$ and β. In other words, every point in R₂/R₁,$\alpha$ and β space has an associated degree of alternation. A visualisation of the degree of alternation map in R₂/R₁,$\alpha$ and β space is shown in Fig. 10. Note that although $0 \le R_{2} /R_{1}$, for convenience Fig. 10 runs only from $0 \le \frac{{R_{2} }}{{R_{1} }} \le 1$.

Fig. 10

Front (left) and back (right) view of the degree of alternation map in $R_{2} /R_{1}$, $\alpha$ and $\beta$ space. Colours correspond to particular values for degree of alternation (key centre)

The red areas in Fig. 10 represent those areas with the largest degree of alternation, i.e. where $R_{2} /R_{1}$ or $\beta$ are highest. The light blue layer corresponds to the first-order surface shown in Figs. 7 and 8, where distributions display no degree of alternation. This layer acts as a crossing-point, above which distributions display “normal” alternating behaviour, and below which distributions display “inverse” alternation. Figure 11 illustrates these results. Note that in Fig. 11, ${\alpha}$ and ${\beta}$ remain constant while $R_{2} /R_{1}$ varies. While experimentally so far only “normal” alternating distributions have been found, the discovery of “inverse” alternating behaviour is theoretically possible.

Fig. 11

Examples of distributions obtained above the first-order surface (top) showing “normal” alternating behaviour, on the first-order surface (middle), and below the first-order surface (bottom) showing “inverse” alternating behaviour

2.3.3 Implications for Second-Order Fitting with “Gently” Alternating Distributions

“Gently” alternating distributions display significantly less $c_{2} r_{2}^{n}$ character than $c_{1} r_{1}^{n}$ character. Such distributions lie close to the first-order surface in $R_{2} /R_{1}$, $\alpha$ and $\beta$ space as illustrated in Figs. 7 and 8. In back-solving to determine $\alpha$, $\beta$, $R_{1}$, and $R_{2}$, “gently” alternating distributions present a unique difficulty in the fitting process when experimental error is taken into account.

The first-order surface in Figs. 7 and 8 is constructed by lines along which $r_{1}$ is constant and along which distributions appear identical. Second-order distributions which are similar in appearance, but which lie close to the first-order surface, may accurately be described using quite different $R_{2} /R_{1}$, $\alpha$ and $\beta$ values.

By way of illustration, consider the hypothetical second-order distribution shown in Fig. 12, as may be detected via GC from 1-hexene onward. If, due to experimental error, a quantity of 1-hexene evaporates prior to GC measurement, then the back-solving method produces a different fit for the data set slightly perturbed at C₆.

Fig. 12

Hypothetical second-order “gently” alternating distribution

Comparing the fitted parameters in Fig. 12 with those in Fig. 13, there are significant differences between values for $\alpha$, $\beta$, $R_{2} /R_{1}$, $c_{2}$ and $r_{2}$. The only reliable metric for “gently” alternating distributions is the first-order $c_{1} r_{1}^{n}$ term. In borderline cases where a second-order distribution lies too close to the first-order surface to obtain reliable values for $\alpha$, $\beta$, $R_{1}$, and $R_{2}$, pressure studies (see Sect. 2.4) and/or repeat experiments will be necessary.

Fig. 13

Hypothetical second-order “gently” alternating distribution with some 1-hexene missing

2.4 The Effect of Pressure in Second-Order Steady-State Oligomerisation

In Sect. 2.2, a catalytic scheme (Scheme 5) for the production of alternating 1-alkene oligomer distributions was presented. The terms $M\left( {C_{2} } \right)_{n}$, $M\left( {C_{2} } \right)_{n + 1}$ etc. group together all intermediates involving $n$ ethylene units. In order to investigate the effect of changing ethylene pressure on $\alpha$ and $\beta$, it is necessary to propose a more detailed scheme with the individual elementary reaction steps that occur within each active species “box” (i.e. between each type of $M\left( {C_{2} } \right)_{n}$ species). This so-called candidate scheme is kept as general as possible, for application to any processes proceeding via a metallacyclic mechanism (Scheme 6):

Scheme 6

Candidate scheme proposed for individual reaction steps occurring between each type of metallacyclic $M\left( {C_{2} } \right)_{n}$ species, with selected rate constants labelled

Following a derivation (see SI.4.1), it can be shown that the ratio between the probabilities for double ethylene insertion β and single insertion α is defined as:

$$\frac{{p\left( {2 Et inserted} \right)}}{{p\left( {1 Et inserted} \right)}} = \frac{{\frac{{k_{2} k_{7} P_{et} }}{{k_{H} }}}}{{k_{5} \left( {k_{3} + k_{6} + k_{7} + k_{10} } \right) + \frac{{k_{2} k_{6} P_{et} }}{{k_{H} }}}}$$

(24a)

So long as $k_{2}, k_{7} \ne 0$, the relationship between the probabilities $\alpha$ and $\beta$ and ethylene pressure $P_{et}$ can be expressed as:

$$\frac{\beta }{\alpha } = \frac{{P_{et} }}{{g_{1} + g_{2} P_{et} }}$$

(24b)

where:

$$g_{1} = \frac{{k_{5} \left( {k_{3} + k_{6} + k_{7} + k_{10} } \right)k_{H} }}{{k_{2} k_{7} }} = {\text{a constant across}} P_{et}$$

$$g_{2} = \frac{{k_{6} }}{{k_{7} }} = {\text{a constant across}} P_{et}$$

2.4.1 Comparison of Experimentally Determined Results vs. Fitted Data

A series of oligomerisation experiments have been performed at different ethylene pressures, using N,N-bis((1H-benzimidazol-2-yl)methyl)-N-methylamine chromium(III) chloride/MMAO-12 as the catalyst (experimental details as per ref. 17). The $\alpha$ and $\beta$ values for each run were obtained using the back-solving method described in Sect. 2.2, and the values of ${\beta }/{\alpha }$ were plotted vs. ethylene pressure. A curve following Eq. 24b was fitted to these calculated values, using a minimisation of the sum of squared differences to determine the most suitable values for $g_{1}$ and $g_{2}$. The results are given in Table 4 and Fig. 14.

Table 4 Key parameters in the fitting of second-order oligomerisation parameters vs. ethylene pressure

Fig. 14

Experimentally determined $\beta /\alpha$ values vs. pressure, together with fitted curve according to Eq. 24b

The fitting of the data in Table 4 with Eq. 24b gives values for $g_{1} = 0.177$ and $g_{2} = 0.176$. Note that $g_{2}$ is equivalent to the ratio of rate constants for single insertion vs. double insertion from a double-ethylene associated intermediate. There is a good degree of correlation between experimentally determined and the fitted values of ${\beta}/{\alpha }$ in Fig. 14 (goodness of fit R² = 0.974).

A final note should be added, that for ethylene oligomerisations at lower pressures with these catalysts, it is possible that the most abundant oligomeric product, 1-butene in this case, is incorporated during the oligomerisation reaction. We have observed this with some of the most active catalysts, resulting in a minor distribution of 2-ethyl-1-alkenes with the opposite alternating pattern, alongside the major distribution of 1-alkenes with the normal alternating behaviour [17].

2.5 Steady-State Oligomerisation with Nonconstant Propagation Probabilities

In assuming constant values for the propagation probabilities $\alpha$ and $\beta$ in Schemes 4 and 5 (Sects. 2.1 and 2.2 respectively), recurrence relations with constant coefficients are generated. The solution of these may be expressed in equation format (see Eqs. 8 and 11). If, however, the propagation probabilities adopt nonconstant values in the catalytic cycle, then recurrence relations with variable coefficients are generated. The aim in such cases is to determine a unique value for each unknown parameter, and/or to elucidate the key mechanistic features of the system.

In this sub-section, a method for handling first-order recurrence relations with variable coefficients is presented. The method finds application in oligomerisation systems where $\alpha$ changes at some point in the catalytic cycle. Second-order systems with nonconstant propagation probabilities will be discussed later, with a particular focus on selective-tetramerisation.

2.5.1 First-Order Recurrence Relations with Nonconstant Propagation Probabilities

Experimental LAO distributions which deviate from Schulz–Flory have been observed, for example chromium(III) bis(carbene)pyridine systems reported by McGuinness and co-workers (Fig. 15) [29]. The observation has been rationalised by less favourable elimination at the early stages of metallacyclic growth [29]. In Fig. 15, $\alpha$ is variable before C₈, then adopts a constant value from C₈ onwards.

Fig. 15

Distribution of 1-alkene oligomers as relative mol(n) vs n (units of ethylene) produced with complex 2,6-bis(1-(2,6-diisopropylphenyl)imidazol-2-ylidene)pyridine chromium(III) chloride. Conditions: 10 μmol Cr, 5 bar continuous ethylene pressure, 25 °C, MAO (5 mmol), toluene (50 mL), 30 min. Activity (C₄–C₂₀ LAOs) = 159 g mmol⁻¹ h⁻¹ bar⁻¹. Variable $\alpha$ and constant $\alpha$ regions indicated

Data from C₈ (n = 4) onward may be treated with the back-solving method described in Sect. 2.1 to find the constant ${\alpha}$ value and R, the total mol of C₈₊ oligomers produced during the reaction time (i.e. the rate at which intermediate species arrive at . Adapting Scheme 4 to show rates, the first nonconstant term $T_{1}^{\prime}$ (i.e. C₆) may be appended, along with $R^{\prime}$, the rate at which intermediate species arrive at (Scheme 7).

To find $\alpha^{\prime}$ (i.e. the $\alpha$ value at C₆) and $R^{\prime}$, two distinct equations involving the two unknown variables may be identified (note that $T_{1}^{\prime}$ is the experimental data point for 1-hexene):

$$\alpha^{\prime}R^{\prime} = R$$

(25a)

$$\left( {1 - \alpha^{\prime}} \right)R^{\prime} = T_{1}^{\prime}$$

(25b)

These equations sum to give:

$$R^{\prime} = T_{1}^{\prime} + R$$

(25c)

Given that there are two distinct equations for two unknowns ($\alpha^{\prime}$ and $R^{\prime}$), a unique solution does exist. $R^{\prime}$ may be found with Eq. 25c, then substituted into Eqs. 25a or 25b to find $\alpha^{\prime}$. The second non-constant term $T_{2}^{\prime}$ (i.e. C₄) may now be appended to Scheme 7, and the $\alpha$ value at C₄ found using the same methodology. The $\alpha$ values calculated for Fig. 15 via this methodology are 0.86, 0.65 and 0.54 from C₄, C₆ and C₈₊ respectively. In other words, for this particular chromium catalyst, the termination probability is less favourable at the early stages and increases as the metallacycle grows from a metallacyclopentane to—heptane and—nonane. Once the metallacycle has reached a certain size, the termination probability becomes constant.

Scheme 7

Catalytic scheme for Fig. 15 at the boundary of the variable $\alpha$ and constant $\alpha$ regions. Active species shown in black, rates and eliminated products shown in blue

2.5.2 Second-Order Recurrence Relations with Nonconstant Propagation Probabilities

McGuinness, Britovsek and co-workers produced a detailed experimental and theoretical study of single- and double-coordination mechanisms in ethylene tri- and tetramerisation with Cr/PNP catalysts [30]. With PNP = Ph₂PN(ⁱPr)PPh₂, high selectivity for 1-octene is observed (~ 70 mol%), together with approximately 16 mol% 1-hexene, 3 mol% C₁₀₊ α-olefins, and various other coproducts (Fig. 16).

Fig. 16

Selective α-olefin distribution as mol(n)% vs n (units of ethylene) obtained with Cr/Ph₂PN(ⁱPr)PPh₂/MAO [16, 30]. Conditions: 10 μmol Cr, 30 bar, 30 °C, MAO (3 mmol), toluene (100 mL), 30 min. Inset: enlargement of distribution from C₁₆₊. Inset: Expanded mol(n)% vs. n for the C₁₆–C₃₂ fraction

Computational studies employing the M06L functional on simplified model compounds have indicated that a single- and double-coordination mechanism operates (see Fig. 17) [30,31,32,33]. Elimination from a chromacyclopentane intermediate to give 1-butene is calculated to be energetically unfavourable [30]. Competing processes to 1-hexene and 1-octene formation diverge at a singly-coordinated chromacyclopentane intermediate 7, with 1-hexene formed predominately via single insertion from intermediate 7, and 1-octene formed predominantly via a bis(ethylene) metallacyclic intermediate 17, with additional involvement of a mono(ethylene) route through intermediate 18 [30].

Fig. 17

Condensed energy profiles for 1-hexene and 1-octene formation from the singly-coordinated chromacyclopentane intermediate 7 [30]. Gaussian09, M06L/quadruple-ζ valence def2-QZVP basis set and SDD ECP on Cr, 6-311 + G(2d,p) basis set on all other atoms

Given the varying differences in energy between the early intermediates in Fig. 17, the selective-tetramerisation process in Fig. 16 may be mathematically considered as a second-order recurrence relation with variable coefficients. Data from C₁₆ onwards corresponds to the tail end of an alternating distribution, where $\alpha$ and $\beta$ values are constant by virtue of negligible geometry and energy differences between larger chromacyclic species. For the purpose of the present discussion, setting aside the difficulties in fitting gently alternating distributions as outlined in Sect. 2.3, the region from C₁₆ onwards may conceivably be treated with the back-solving method described in Sect. 2.2 to find the constant values for $\alpha$ and $\beta$ as well as $R_{1}$, the rate at which intermediate species arrive at , and $R_{2}$, the rate at which intermediate species originating from before arrive at . Adapting Scheme 5 to show rates, the first non-constant term $T_{1}^{\prime}$ (i.e. C₁₄) may be appended, along with $R_{1}^{\prime}$, the rate at which intermediate species arrive at , and $R_{2}^{\prime}$, the rate at which intermediate species originating from before arrive at (Scheme 8).

Scheme 8

Catalytic scheme for Fig. 16 at the boundary of the variable $\alpha, \beta$ and constant $\alpha , \beta$ regions

In an effort to find $\alpha^{\prime}$ and $\beta^{\prime}$ (i.e. the $\alpha$ and $\beta$ values at C₁₄) and $R_{1}^{\prime}$ and $R_{2}^{\prime}$, three distinct equations involving the four unknown variables may be identified (note that $T_{1}^{\prime}$ is the experimental data point for 1-tetradecene):

$$\beta^{\prime}R_{1}^{\prime} = R_{2}$$

(26a)

$$\alpha^{\prime} R_{1}^{\prime} + R_{2}^{\prime} = R_{1}$$

(26b)

$$\left( {1 - \alpha^{\prime} - \beta^{\prime}} \right)R_{1}^{\prime} = T_{1}^{\prime}$$

(26c)

These equations sum to give:

$$R_{1}^{\prime} + R_{2}^{\prime} = R_{1} + R_{2} + T_{1}^{\prime}$$

(26d)

Given that there are three distinct equations for four unknowns ($\alpha^{\prime}, \beta^{\prime}, R_{1}^{\prime} and\,R_{2}^{\prime}$), a unique solution does not exist in this case. If no further information is available/applied from e.g. theoretical studies, an assumption (such as $\frac{{\alpha^{\prime}}}{{\beta^{\prime}}} = \frac{\alpha }{\beta }$) is required as a fourth equation to determine unique solutions for the four unknown variables.

3 Modelling Non-Steady-State Polymerisation Processes

3.1 Connecting Schulz–Flory Distributions and Poisson Distributions

The mathematics behind Schulz–Flory ethylene oligomer distributions was first accurately described in 1936 [19]. In 1940, Flory went on to give a kinetic derivation of Poisson’s equation for polymer distributions obtained from a monofunctional monomer (ethylene oxide) [34]. Yet the mathematical connection between Schulz–Flory distributions and Poisson distributions has not been fully elucidated. In this section, the relationship between the two is made explicit from first principles.

It should be noted that the treatment presented below is described with ethylene as the monomer, but is generally applicable to oligomerisation and polymerisation processes involving the sequential incorporation of monomers in forming linear product chains.

3.1.1 First Principles

The mechanics of a first-order ethylene oligomerisation or polymerisation process whereby propagation proceeds via metal alkyl or metallacyclic species is summarised in Scheme 9. In this model and the following mathematical treatment, we assume that the rate constants k_p and k_e are constant over time and independent of chain length. Furthermore, all active species must start the polymerisation instantaneously and there is no catalyst deactivation.

Scheme 9

Catalytic scheme for first-order propagation

$M\left( n \right)$ represents catalyst intermediates involving $n$ ethylene units (e.g. $M\left( 0 \right)$ = metal hydride species, $M\left( 1 \right)$ = metal ethyl or metallacyclopropane species, $M\left( 2 \right)$ = metal butyl or metallacyclopentane species etc.), $k_{p}$ is the rate constant of propagation, $k_{e}$ is the rate constant of elimination, and $L$ is the total precatalyst loading.

A general expression for $M\left( n \right)$, as a function of time, may be derived (see SI.5.1):

(27)

The label “non-steady-state” refers to time dependent terms, whilst the label “steady-state” refers to time independent terms. It is worth noting again that Scheme 9 assumes $k_{p}$ and $k_{e}$ are constant over time and independent of $M\left( n \right)$, i.e. no distinction in propagation or elimination kinetics is made between the various intermediate groups in the catalytic scheme. In addition, by combining $M\left( n \right)$ species in forming the kinetics equations, it is implicitly assumed that a single species within each $M\left( n \right)$ group dominates. For instance, the $M\left( 0 \right)$ group includes both metal hydrides with a vacant site and ethylene-associated metal hydrides. This assumption is no longer required for the Schulz–Flory case, given that steady-state kinetics operate.

3.1.2 Derivation of the Schulz–Flory Equation for First-Order Steady-State Oligomerisation

Referring to Scheme 9, the scenario where steady state is achieved from $t_{0}$ (i.e. instantaneously from the beginning of reaction) with $k_{e} \ne 0$ corresponds to the limiting case of first-order steady-state oligomerisation. In such a case, by ignoring non-steady-state terms, Eq. 27simplifies to:

$$M\left( n \right) = \frac{{Lk_{p}^{n} k_{e} }}{{(k_{p} + k_{e} )^{{\left( {n + 1} \right)}} }}$$

(28)

Given that elimination occurs at a constant rate, the interest lies only in the sum of linear alpha olefins of chain length $n$ ethylene units continuously eliminated from $t_{0}$ to $t_{end}$:

$$mol(n)_{{t_{end}}} = k_{e} \mathop \int \limits_{0}^{{t_{end} }} M\left( n \right)dt = \frac{{Lk_{p}^{n} k_{e}^{2} t_{end} }}{{(k_{p} + k_{e} )^{{\left( {n + 1} \right)}} }} = \frac{{Lk_{e}^{2} t_{end} }}{{k_{p} }}\alpha^{{\left( {n + 1} \right)}} \text{where } \alpha = \frac{{k_{p} }}{{\left( {k_{p} + k_{e} } \right)}}$$

(29a)

$$mol(n)\% = \frac{{\frac{{Lk_{e}^{2} t_{end} }}{{k_{p} }}\alpha^{{\left( {n + 1} \right)}} }}{{\frac{{Lk_{e}^{2} t_{end} }}{{k_{p} }}\left[ {\alpha^{2} + \alpha^{3} + \alpha^{4} + \ldots} \right]}} = \frac{{\alpha^{{\left( {n + 1} \right)}} }}{{\alpha^{2} + \alpha^{3} + \alpha^{4} + \ldots}}$$

$$\alpha^{2} + \alpha^{3} + \alpha^{4} + \ldots = \frac{{\alpha^{2} }}{{\left( {1 - \alpha } \right)}} \text{when} \left| \alpha \right| < 1, \text{so}:$$

$$mol\left( n \right)\% = \frac{{\left( {1 - \alpha } \right)}}{{\alpha^{2} }}\alpha^{{\left( {n + 1} \right)}} = \left( {1 - \alpha } \right)\alpha^{{\left( {n - 1} \right)}}$$

(29b)

Equation 29a is the unnormalised Schulz–Flory equation, and Eq. 29b is the normalised Schulz–Flory equation. Note in Eq. 29a that the total mol of oligomers produced during the reaction time $R = L\alpha k_{e} t_{end}$. Given that both $\frac{{Lk_{e}^{2} t_{end} }}{{k_{p} }}$ (where $L$ and $t_{end}$ are known) and $\alpha$ may be fitted, $k_{p}$ and $k_{e}$ may be inferred.

3.1.3 Derivation of Poisson’s Equation for First-Order Living Polymerisation

Referring again to Scheme 9, setting $k_{e} = 0$, i.e. there is no chain termination process, corresponds to the limiting case of first-order non-steady-state living polymerisation. In such a case, Eq. 27 simplifies to:

$$M\left( n \right)_{t} = \frac{L}{{n!k_{p} }}\left[ {e^{{ - k_{p} t}} \left( {k_{p}^{{\left( {n + 1} \right)}} t^{n} } \right)} \right] \leftrightarrow M\left( n \right)_{t} = \frac{{L(k_{p} t)^{n} }}{n!}e^{{ - k_{p} t}}$$

(30)

Given that no termination occurs throughout, the interest lies only in the distribution of intermediate species at $t_{end}$ (i.e. the end of reaction), which would give n-alkene products on thermal termination (as per the INEOS (Ethyl) process [8]) or n-alkane products upon hydrolysis. At this point:

$$M\left( n \right)_{{t_{end} }} = mol\left( n \right)_{{t_{end} }} = \frac{{L(k_{p} t_{end} )^{n} }}{n!}e^{{ - k_{p} t_{end} }} = \frac{{L\lambda^{n} }}{n!}e^{ - \lambda }$$

(31a)

$$mol\left( n \right)_{{t_{end} }} \% = \frac{{\frac{{L\lambda^{n} }}{n!}e^{ - \lambda } }}{{\mathop \sum \nolimits_{n = 1}^{n = \infty } \frac{{L\lambda^{n} }}{n!}e^{ - \lambda } }} = \frac{{L\lambda^{n} e^{ - \lambda } }}{{n!Le^{ - \lambda } \mathop \sum \nolimits_{n = 1}^{n = \infty } \frac{{\lambda^{n} }}{n!}}} = \frac{{\lambda^{n} e^{ - \lambda } }}{{n!e^{ - \lambda } \left( {e^{\lambda } - 1} \right)}} = \frac{{\lambda^{n} e^{ - \lambda } }}{{n!\left( {1 - e^{ - \lambda } } \right)}}$$

(31b)

Equation 31a is Poisson’s equation scaled by the precatalyst loading $L$, and Eq. 31b is Poisson’s equation scaled by the factor $\left( {1 - e^{ - \lambda } } \right)^{ - 1}$ (which quickly approaches 1 as $\lambda$ increases). Note that $\lambda$ is the mean number of ethylene units per chain. Given that $\lambda$ may be fitted (where $t_{end}$ is known), $k_{p}$ may be inferred.

3.2 Non-Steady-State Living Polymerisation with Catalysed Chain Growth

In 1952, Ziegler et al. reported the Aufbaureaktion involving the formation of long chain trialkylaluminium species via stepwise insertion of ethylene monomers into Al–C bonds (starting with triethylaluminium, for instance) [35]. The effect of transition metals, such as colloidal Ni, in suppressing chain growth in the Aufbaureaktion to selectively form 1-butene was subsequently discovered [36]. This prompted a screening of the d-block transitional metals for other metal effects and, in the case of titanium, this led to the first reported transition metal catalysed ethylene polymerisation [37, 38].

In their early work, Ziegler et al. described the process as a transition metal catalysed Aufbaureaktion, with the alkylaluminium functioning as the polymerisation catalyst and the transition metal compound functioning as the co-catalyst [37]. Later it was found that the polymer chain grows on the transition metal centre, rather than on aluminium. More recently however, several processes have been reported where the polymer chain is transferred during the course of the reaction from a transition metal to a main group metal such as aluminium, magnesium or zinc [39,40,41,42,43,44]. If this chain transfer is very fast and reversible, then it appears as if the polymer chains grow on the main group metal centres instead. Such a reaction is generally described as a catalysed chain growth (CCG) process or a coordinative chain-transfer polymerisation (CCT), or, to employ the original terminology, a transition metal catalysed Aufbaureaktion [44,45,46].

Tremendous advances have been made in CCG during the past 2 decades, but the mathematics of generating product distributions via CCG has not yet been reported. Such a treatment from first principles is presented here and then applied next.

3.2.1 First Principles

Scheme 10 summarises the mechanics of a first-order ethylene oligomerisation or polymerisation process whereby propagation proceeds via metal alkyl or metallacyclic species, and chain transfer to and amongst a population of non-propagating agents is possible. A similar caveat as in Sect. 3.1 applies in that in this model and the following mathematical treatment, we assume that the rate constants k_p, k_e and k_ct are constant over time and independent of chain length. Furthermore, all active species must start the polymerisation instantaneously and there is no catalyst deactivation.

Scheme 10

Catalytic scheme for first-order propagation with chain transfer

$M\left( n \right)$ represents transition metal intermediates involving $n$ ethylene units, $CT\left( n \right)$ represents the population of chain transfer species bearing an alkyl chain of length $n$ ethylene units, $k\left( {M\left( n \right)} \right)_{p}$ is the rate constant of propagation for $M\left( n \right)$, $k\left( {M\left( n \right)} \right)_{e}$ is the rate constant of elimination for $M\left( n \right)$, and $k_{ct}$’s are rate constants of chain transfer between species. Scheme 10 is not labelled with individual $k_{ct}$’s given the general nature of the present treatment.

It is possible to show (see SI.5.2) that:

$$\frac{d}{dt}M\left( n \right)_{t} = k(M\left( {n - 1} \right))_{p} M\left( {n - 1} \right)_{t} - \left[ {k(M\left( n \right))_{p} + k(M\left( n \right))_{e} } \right]M\left( n \right)_{t}$$

(32)

This expression for change does not show dependence on chain transfer terms. In other words, the underlying $M\left( n \right)$ system drives towards the distribution that it would adopt in the absence of chain transfer agent. However, the nature of instantaneous steady state for CT species slows down the overall rate at which the $M\left( n \right)$ system evolves vs. the case with no chain transfer.

3.2.2 Numerical Simulations

Two illustrative examples are provided to demonstrate the preceding mathematics. Numerical simulations of evolving oligomerisation or polymerisation systems which involve chain transfer, as per Scheme 10, are possible with Microsoft Excel (see SI.5.3 for details).

3.2.3 Illustrative Example 1: Living Polymerisation with Catalysed Chain Growth

See Scheme 11, Figs. 18, 19, 20.

Scheme 11

Scheme for a living polymerisation system with CCG

Fig. 18

System evolution and product distribution with CT agent loading = 0 μmol dm⁻³ (0 eqv. vs. transition metal catalyst). Note: M(n) and CT(n) up to n = 5 illustrated left and centre, all n captured right

Fig. 19

System evolution and product distribution with CT agent loading = 10 μmol dm⁻³ (1 eqv. vs. transition metal catalyst). Note: M(n) and CT(n) up to n = 5 illustrated left and centre, all n captured right

Fig. 20

System evolution and product distribution with CT agent loading = 100 μmol dm⁻³ (10 eqv. vs. transition metal catalyst). Note: M(n) and CT(n) up to n = 5 illustrated left and centre, all n captured right

3.2.4 Illustrative Example 2: First-Order Oligomerisation with Catalysed Chain Growth

See Scheme 12, Figs. 21, 22, 23.

Scheme 12

Scheme for a first-order oligomerisation system with CCG

Fig. 21

System evolution and product distribution with CT agent loading = 0 μmol dm⁻³ (0 eqv. vs. transition metal catalyst). Note: M(n) and CT(n) up to n = 5 illustrated left and centre, all n captured right

Fig. 22

System evolution and product distribution with CT agent loading = 10 μmol dm⁻³ (1 eqv. vs. transition metal catalyst). Note: M(n) and CT(n) up to n = 5 illustrated left and centre, all n captured right

Fig. 23

System evolution and product distribution with CT agent loading = 100 μmol dm⁻³ (10 eqv. vs. transition metal catalyst). Note: M(n) and CT(n) up to n = 5 illustrated left and centre, all n captured right

3.2.5 The Fitting Process

The preceding discussion demonstrated how Poisson distributions or Schulz–Flory distributions are still generated for first-order systems involving CCG. Section 2.1 provided the details for modelling first-order steady state oligomerisation systems, which by extension may be applied to first-order oligomerisation systems involving CCG. Section 3.1 gave a derivation to the Poisson distribution for living polymerisation systems, ending at Eqs. 31a and 31b.

When Eqs. 31a and 31b are applied to living polymerisation systems with CCG, in order to capture all intermediate alkyl chains, $L$ is now defined as:

$$L = Transition~metal~catalyst~loading + \left( {CT loading} \right)\left( {Alkyl~chains~per~CT~species } \right)$$

(33)

Furthermore, $k_{p}$ may now be interpreted as the effective propagation rate constant for the entire system. The process for modelling living polymerisation systems with CCG is described next, along with derivations of the experimentally-relevant metrics $M_{n}$, $M_{w}$ and $PDI$. For illustration, a comparison is also made between an experimentally determined result vs. the fitted data.

3.2.6 Back-Solving to Determine L and λ, and Checking for the Best Fit

In practice, the $mol\left( n \right)_{{t_{end} }}$ values in Eq. 31a are known (i.e. mol of each oligomer produced during the reaction time), and the aim is to find $L$ and $\lambda$. See SI.5.4 for a discussion on the method used to determine $\alpha$ and $R$, as well as a checking method to determine that the best fit has been found.

3.2.7 Derivations of Further Experimentally-Relevant Metrics

For a distribution of n-alkenes, the following expressions can be derived (see SI.5.5 and SI.5.6):

$$wt(n)_{{t_{end} }} \% = \frac{{n\lambda^{n - 1} e^{-{\lambda}} }}{n!}$$

(34)

$$M_{{n}{t_{end}}} = \left( {RFM C_{2} H_{4} } \right)\lambda(1-e^{-\lambda})^{-1}$$

(35)

$$M_{{w}{{t_{end}} }} = \left( {RFM C_{2} H_{4} } \right)\left( {\lambda + 1} \right)$$

(36)

$$PDI = \frac{{M_{{w}{{t_{end}} }} }}{{M_{{n}{{t_{end}}}} }} = \frac{{\left( {\lambda + 1} \right)\left( {1 - e^{ - \lambda } } \right) }}{\lambda }$$

(37)

As $\lambda$ increases, the factor $\left( {1 - e^{ - \lambda } } \right)$ approaches 1. Under such circumstances, Eqs. 35 and 37 simplify.

$$M_{{n}{{t_{end}} }} \to \left( {RFM C_{2} H_{4} } \right)\lambda$$

(38)

$$PDI \to \frac{{\left( {\lambda + 1} \right) }}{\lambda } = 1 + \frac{1}{\lambda } \to 1$$

(39)

With high values of $\lambda$, $PDI$ tends towards unity, as per the literature [47].

3.2.8 Comparison of an Experimentally Determined Result vs. Fitted Data

CCG polymerisations have been performed, using ZnEt₂ as the CT agent with 2,6-diacetylpyridinebis(2,6-diisopropylanil)FeCl₂/MAO (experimental details, see SI.5.7). Termination of the polymerisation reaction is carried out by hydrolysis, resulting in a mixture of long-chain alkanes. Given a highly symmetric Poisson distribution of the resulting n-alkane oligomers, quantitative ¹³C{¹H} NMR analysis may be used to estimate the mean chain length by comparing integrals for CH₃ with those for CH₂. For the given experiment, quantitative ¹³C{¹H} NMR spectroscopy gives a CH₃:CH₂ ratio of 2:45.38, implying a mean chain length between 46 and 48 carbons (see Fig. 24).

Fig. 24

Quantitative ¹³C{¹H} NMR spectrum of n-alkane oligomer product (solvent: d₂-1,1,2,2-tetrachloroethane/1,2,4-trichlorobenzene 50/50 v/v), 100 MHz, 100 °C

Experimental characterisation of the alkane distribution in such mixtures is non-trivial. The average molecular weight is typically too low for accurate GPC analysis, and the longer alkanes are often insufficiently volatile for conventional quantitative GC-FID analysis. Only the first part of the n-alkane oligomer distribution is usually quantifiable via GC. These signals have been used here as inputs for the fitting model. Following the procedure previously described in this sub-section, the values in Table 5 are found.

Fig. 25

Modelled outputs (see Table 5) vs. experimental GC data

3.2.9 Model Outputs

Table 5 Modelled outputs using experimental GC data

There is a high degree of correlation between calculated and experimental values in Table 5 (goodness of fit R² = 0.996), and the calculated mean chain length of 46.6 carbons is within the range previously specified by quantitative ¹³C{¹H} NMR spectroscopy. This method therefore allows the entire product distribution to be determined from a subset of experimental data obtained by GC analysis.

4 Conclusions and Summary

A variety of different distributions—including Schulz–Flory, alternating, selective-LAO, and Poisson distributions—of linear hydrocarbon products can be obtained from ethylene oligomerisation and polymerisation processes. Oligomerisation processes are steady-state in nature, and can be mathematically described using recurrence relations. First order oligomerisation processes propagate via sequential single-insertions of ethylene, whereas second order oligomerisation processes propagate via sequential single- and double-insertions of ethylene. From an analytical point of view, an important consideration with regards to second order processes is the degree of alternation exhibited by a given distribution, which is generally more pronounced at higher ethylene pressure. The degree of alternation for a particular distribution may be defined as the relative ratio of alternating vs. Schulz–Flory character. An overview of the key analytical results from first order and second order processes with constant probabilities of propagation is given in Table 6.

Table 6 Key analytical results from 1st order and 2nd order processes with constant probabilities of propagation

Examples of first order and second order oligomerisation processes which involve varying probabilities of propagation have been identified from the literature. Such systems include selective oligomerisation catalysts such as the well-known chromium-based tetramerisation catalysts. An analytical method is available for handling first order processes with varying probabilities of propagation. However, there is a theoretical limitation to a similar approach for analogous second order processes.

The connection between distributions from first order steady-state oligomerisations (Schulz–Flory) and first order non-steady-state polymerisations (Poisson) is made through a consideration of kinetic differential equations. The addition of a chain transfer agent (e.g. ZnEt₂) to such systems results in a slowing down of the overall propagation rate. Accurate modelling of a CCG polymerisation system featuring ZnEt₂ as the CT agent has been achieved, allowing for the derivation of several key mechanistic parameters.

The mathematics of ethylene oligomerisation presented here can be easily applied to other systems where chain growth proceeds through an oligomerisation process. Examples of such systems would be the oligomerisation of higher α-olefins to poly(α-olefin) oils, the dehydropolymerisation of aminoboranes, chain growth reactions in the Fischer–Tropsch process or the formation of oligonucleotides. By fitting the experimental results to the formulas described herein for the various oligomerisation scenarios, valuable information may be extracted regarding the underlying reaction mechanisms of the oligomerisation process.