General Applicable Residence Time Distribution Model to Estimate Reaction Rates in a Rotor–Stator Spinning Disc Reactor

Abstract

Many designs of industrial reactors stem from designs from the 1960s–1970s. For a wide range of reactions, these designs lead to suboptimal reaction configurations due to limitations in heat- or mass-transfer. Process intensification has come up with a different approach, resulting in micro- and mini-reactors being safer and more cost-effective on a full industrial scale. However, based on the experience in the suboptimal reactor designs, the reaction rates of these reactions seem too low for full-scale reactions in a mini reactor. We suggest a test for the reaction rate based on a generalized model in combination with a specific type of mini-reactor: the rotor–stator spinning disc reactor. The generalized model is based on existing models on residence time distribution in such reactors. It does not need to be tailor-fitted for a specific rotor–stator spinning disc reactor that is used for the test, as is the case with current models. In this article, we show that our simplifications induce a difference in outcome in reaction rate of less than 10% with the existing models. Experiments with the well-studied chemical reaction of the hydrolysis of acetic anhydride show that the reaction rates calculated based on this scan fall within the range of reported data from the literature.

1. Introduction

Many designs of industrial plants that are currently in use stem from designs from the 1960s–1970s, when large-scale reactors were the only way to realize large-scale production. For a wide range of reactions, large-scale reactors lead to suboptimal reaction configurations due to limitations in heat- or mass-transfer, causing significant energy usage during the production of the required product. An example is highly exothermic reactions, in which the reaction speed is kept low by diluting the reactants to prevent run-away reactions from occurring, causing significant energy consumption for separating the required product from the dilutant. Another example is the occurrence of unwanted side reactions due to inefficiencies in mixing. Such mixing inefficiencies are practically impossible to avoid in large-scale reactors. These side products also need to be separated from the product stream. Since separation processes account for 40–80% of all energy costs in the process industry, it is relevant to prevent or reduce these separation processes if possible.

Since the 1990s, a different approach to realizing large-scale production has emerged—process intensification. Instead of focusing on increasing the scale of the production, the focus is on increasing the speed and control over the reactions occurring. This can be explained in terms of improved mixing on a micro-scale. If mixing occurs at high shear, reactants (either liquid/liquid or liquid/gas) show on a micro-scale, good or at least improved mixing. This means that the actual concentration of molecules close enough to each other to react will be equal or close to the concentration expected based on the concentration in the bulk, which increases the overall reaction rate.

Striving for better micro-mixing process intensification has resulted in a whole myriad of new reactors that allow for industrial-scale production, using reactors on a micro- or mini-scale. In case that the small scale alone does not result in a sufficient increase in micro-mixing, HiGee-reactors can offer a solution. HiGee reactors are mini reactors that use high gravitational forces, usually by spinning or moving parts, thereby creating high levels of shear. This high shear increases both mass- and heat-transfer on a micro level. In this way, process intensification offers a high potential to improve the efficiency of several industrial processes and, thus, decrease the environmental footprint of the chemical industry [1,2].

Currently, there is no list of which processes could benefit from which process-intensive reactor types and what the impact of this benefit would be on the total performance of a factory. A major barrier to creating this insight is that for many industrial processes, the conversion rate is determined in situations where the speed of reaction is severely limited by mass- and/or heat-transfer and, therefore, gives the false impression that a reaction is too slow to realize full-scale industrial production in the small size process intensification reactors. A specific HiGee reactor that is suited to overcome most mass- and heat-transfer limitations is a spinning-disc reactor.

There are models available to determine the reaction rates based on experiments using a spinning-disc reactor. However, these models either require specific knowledge and equipment or need to be tailor-fitted to the specific reactor they apply to. Therefore, the use of these models requires serious resources. While in the (pre-)feasibility phase of a reactor replacement, such resources are generally not available yet. Resources are only freed after a brief list of potential solutions is made. This means that, in general, process intensification options do not make it to the shortlist.

To tackle this problem at the source, we decided to develop a test to determine reaction rates in a rotor–stator spinning disc reactor using a generalized model that can be applied to all rotor–stator spinning disc reactors without additional fitting. This test allows for the determination of kinetic data for an industrial chemical reaction at a minimum limitation in mass- and heat-transfer, under conditions that can be reproduced at an industrial scale. The rotor–stator spinning disc reaction is a HiGee reactor that has excellent mass- and heat-transfer and is accessible to both educational institutions and industrial parties. In this way, we want to make it more feasible to show the potential of micro and mini reactors in full-scale industrial production processes in a low-budget way in order to raise the interest required to realize full-scale investigations for the optimum reactor settings.

2. Current Understanding of Rotor–Stator Spinning Disc Reactors

In this chapter, the literature on heat- and mass-transfer in closed rotor–stator spinning disc reactors is discussed, including the phenomena causing these characteristics, the influences on reaction speed and selectivity, and the way to model these phenomena.

2.1. Mechanical Description of Rotor–Stator Spinning Disc Reactors

A rotor–stator spinning disc reactor is a reactor that spins one or more discs in a space between stators. The gap between the stator and rotor is small, typically in the order of 0.5–2 mm. A schematic picture of a multi-stage rotor–stator spinning disc reactor is shown in Figure 1. Reactants (gasses, liquids, and/or suspensions) are pumped in at the top and flow around the rotors, exiting as one stream at the bottom while the rotors spin at high speed. This is different from some other types of spinning disc reactors described in the literature, where the liquid flows as a film on the rotor without a stator and/or exits at the side of the disc [3]. The rotor–stator spinning disc reactor has higher shear forces (so better mixing) and better heat transfer than other spinning disc reactors because of the presence of the stator and can be applied as a multi-rotor device (typically three rotors) because the sides are closed. In this article, we mean a multi-stage rotor–stator spinning disc reactor when we mention spinning disc reactors.

2.2. Hydrodynamics of Single-Phase Reactions in a Spinning Disc Reactor

As mentioned in the Introduction, mass and heat-transfer need to be very good when determining the reaction rate. Normally, this requires a certain degree of turbulence. Turbulence is defined as a function of the rotor–stator configuration.

The flow of fluids in a rotor–stator configuration is called a Couette flow. Daily and Nece (1960) [6] described the different flow regimes. They differentiated four different flow regimes characterized by the degree of turbulence in terms of the Reynolds number, Re, and the dimensions of the rotor–stator configuration in terms of the Couette number, G.

The Reynolds number is defined as:

$$Re={\displaystyle \frac{\omega ·R^2}{\nu }}$$

where ω = angular momentum (rad/s], R is the rotor radius (m), and ν is the kinematic viscosity (m²/s).

The Couette number is defined as:

$$G={\displaystyle \frac{h}{R}}$$

where R is the rotor radius (m), and h is the distance between the rotor blade and the stator above or below the rotor blade (m).

The third type of flow regime defined by Daily and Nece (1960) [6] represents the turbulent flow that is specific for the flow between tightly spaced rotors and stators, as is the case for a spinning disc reactor. It is referred to as Regime III: turbulent torsional Couette flow.

A visual representation of the four flow regimes together with an illustration of the different dimensions of a rotor–stator are provided by a large number of more recent publications, for example Haddadi and Poncet (2008) [5], De Beer et al. (2014) [7] and Hop et al. (2023) [8].

In the years after the publication by Daily and Nece (1960) [6], this picture is further detailed with the Batchelor and Stewartson flow and their respective effects on Von Kármán and Bödewadt boundary layers. The Stewartson flow is observed in the case of superimposed flow C_w [7]:

$$C_W={\displaystyle \frac{{\mathsf{\Phi }}_V}{\mathsf{\upsilon }·R}}$$

where Φ_V is the flow (m³/s) that is pumped through the spinning disc reactor, R is the rotor radius (m), and υ is the kinematic viscosity (m²/s).

Several authors present CFD calculations that provide insight into when these phenomena occur and how they vary with G, Re, and C_w. For example, Haddadi and Poncet (2008) [5] presents the results of CFD calculations characterizing the flow of a closed single rotor–stator configuration under circumstances that are typical for a spinning disc reactor rotating at high angular speed: Couette number G ~0.02–0.03, a Reynolds number Re > 2 × 10⁵. They calculated that turbulence is almost homogeneous along the axial direction and vanishes towards the disk. They showed that homogeneity of the turbulence in the axial direction between the rotor and stator is a function of position on the rotor, Reynolds number, and C_w.

Other authors have tried to make outcomes more quantitative by combining experimental and numerical work like Kleiner et al. (2018) [9] and Hop et al. (2023) [8]. Kleiner et al. (2018) [9] reported that average Nusselt numbers (as a measure of heat transfer) converged with experimentally determined Nusselt with increasing Reynolds numbers to agree with Reynolds numbers Re = 3 × 10⁵. The measurements were conducted on an enclosed single disc spinning disc (at a C_w = 400). Comparison with data from a multi-stage spinning disc [7] gave qualitative agreement. The difference was attributed to the difference in disc numbers.

Hop et al. (2023) [8] modeled a complete spinning disc reactor and validated the results by calculating residence time distributions for a wide range of Reynolds numbers at a G = 0.03 and a C_w = 150. He achieved good agreement with experimentally determined residence time distributions.

Based on the above, it is clear that CFD provides a deep understanding of the hydrodynamics of the spinning disc reactor. It is even possible to obtain the residence time distribution from such a generally applicable model. However, such models require very specialized equipment and knowledge to operate them and are therefore not suitable for our low-budget application.

Furthermore, they do not provide insight into whether sufficient micro-mixing occurs, is that directly upon entering the turbulent flow regime or do the Batchelor and Stewartson flow with their respective effects on Von Kármán and Bödewadt boundary layers also have a considerable influence? Therefore, in the following paragraphs, we will investigate, respectively, micro-mixing and a macro-scale approach to the determination of the residence time distribution.

2.3. Degree of Micro-Mixing in a Spinning Disc Reactor

As mentioned in the Introduction, the spinning disc has proven to be very efficient in terms of allowing reactants from so-called “immiscible” fluids. This is explained in terms of increased mixing on the microscale of the mixture, decreasing the distance molecules must diffuse to the interface between the two fluids, also referred to as increased micro-mixing.

Manzano Martinez et al. (2017 and 2022) [10,11] measured the increase of micro-mixing by performing the Villermaux–Dushman parallel competitive reaction scheme. He used the relation between mixing time and local energy dissipation rate that is proportional to the frictional torque as proposed by Daily and Nece (1960) [6]. The derivation of this equation is written out in Appendix A of an article by de Beer et al. (2016) [12]. Manzano Martinez et al. (2017 and 2022) [10,11] showed that the mixing time is directly proportional to the local energy dissipation rate divided by the kinetic viscosity:

$$t_{mix}={\displaystyle \frac{{\epsilon }_{loc}}{\mathsf{\upsilon }}}$$

where t_mix is the micro-mixing time (s), ε_loc is the local energy dissipation rate (m²/s) and υ the kinematic viscosity (m²/s). The advantage is that the relation is the same, independent of C_w, G, and Re.

Furthermore, micro-mixing is complete when ε_loc/υ > 1 × 10⁸/s². Which appears to occur, i.e., in this study at Re > 1 × 10⁵, G = 0.007 and C_w~1 [11], which is in line with the Reynolds and Couette number for the turbulent flow indicated at the transition from Regime II to Regime III of Daily and Nece (1960) [6].

Chaudhuri et al. (2021) [13], Kleiner and Hinrichsen (2019) [14], Van Kouwen et al. (2021) [15], Wietelmann et al. (2022) [16] and Magosso et al. (2021) [17] studied the reaction of several types of industrially relevant chemical reactions that suffer from mass transfer limitations due to a limitation in the interfacial surface between the reactants. They showed that the measured reaction speeds and selectivity benefitted from the increased shear and explained that this was caused by an increase in the interfacial surface. The experiments also showed excellent heat control due to the high turbulent flow between the rotor and the stator and the high flow of coolant around the stator in combination with the high surface area compared to the reactor volume. The reactants are never more than 1 mm away from a wall that can be cooled/heated. Due to the limited reaction volume, excellent mixing, and high temperature control spinning disc reactors are inherently safer [13,14,15,16].

Summarizing the work by Manzano Martinez et al. [10,11] shows that for a single-phase reaction, the optimal reactor conditions are reached at an ε_loc/v > 1 × 10⁸/s², which appears to occur when the flow is turbulent (in general at Re > 1 × 10⁵).

2.4. Residence Time Distribution Modeling in a Spinning Disc Reactor

The residence time is determined by the mixing pattern of a reactor on a macro-scale. The residence time is an important variable since the residence time, together with the reaction rate, influences the concentration of the reactants and, thus, determines the conversion of a reaction. Therefore, residence time distributions are used to characterize the mixing pattern on the macro-scale of a reactor. A traditional engineering approach to model the residence time distribution is by combining the models of the residence time distribution of idealized tank models with a known concentration profile: a plug flow reactor (PFR) and one or several continuous stirred tank reactors (CSTR). Modeling residence time distribution with several CSTRs in series is referred to as the tanks-in-series model, where a PFR is modeled as an infinite amount of CSTRs in series. For an introduction to this subject, we refer to the textbook Chemical Reaction Engineering [18].

Table 1. Reaction equations with reaction rates for all modeled reactions. A, B, and C are reactants with $C_{A0}\le C_{B0}$, P is the desired product, and S is the unwanted side product. Conversion is calculated with $X={\displaystyle \frac{C_P}{C_{A0}}}$ (with A:P = 1:1) for residence times of $\tau ={\displaystyle \frac{C_0-C_1}{-r_A}}$ in (each) CSTR and $\tau ={{\displaystyle \int }}_{C_1}^{C_0}{\displaystyle \frac{dC}{-r_A}}$ in PFR [18].

1. First-order reaction	$A\to P$	$-r_A=r_P=k·C_A$
2. Equimolar second-order reaction	$A+A\to 2P$	$-r_A=r_P=k·C_A^2$
3. General second-order reaction	$A+B\to P$	$-r_A=-r_B=r_p=k·C_A·C_B$
4. Equilibrium reaction	$A+B\rightleftarrows P+S$	$-r_A=-r_B=r_P=k_1·C_A·C_B-k_2·C_P·C_S$
5. Reaction in series	$A+B\to P$ $P+C\to S$	$-r_A=-r_B=k_1·C_A·C_B$ $r_P=k_1·C_A·C_B-k_2·C_P·C_C$ $-r_C=k_2·C_P·C_C$
6. Reaction in parallel	$A+B\to P$ $A+C\to S$	$-r_A=k_1·C_A·C_B+k_2·C_A·C_C$ $-r_B=r_P=k_1·C_A·C_B$ $-r_C=k_2·C_A·C_C$

Symbols: r_i = reaction rate of i and C_i = concentration of i, in which i can be A, B, C, S, or P, k = reaction rate constant, C₀ = ingoing concentration of A, and C₁ is the outgoing concentration of A.

of Section 3.2. shows equations used for this type of calculation.

The experimental determination of residence time distribution in a rotor–stator spinning disc reactor and modeling this as multiple CSTRs and/or PFR has been reported in the literature. Visscher et al. (2013) [19] modeled a one-disc rotor–stator spinning disc reactor with a variation on the tanks-in-series model: a measured volume of PFR (determined by the first appearance of the output signal) was combined with a variable number of CSTRs.

Visscher et al. [19] modeled residence time as PFRs and CSTRs in a spinning disc reactor, in this case, a spinning disc reactor with only one disc and three different gap sizes: 1 mm, 2 mm, 3 mm (G = 0.008, 0.015, and 0.0230). Results of Visscher et al. [19] show a range of 7–68% for PFR volume combined with 1.9–2.4 CSTRs. This model was also used by De Beer et al. (2014) [4] for a three-disc rotor–stator spinning disc reactor with both a variable PFR volume and a variable number of CSTRs, where good fits were obtained by combining a fixed number of three CSTRs with a variable volume of PFRs. The reason for a combination of PFR and CSTR (instead of only CSTRs in the standard tanks-in-series model) is that the flow in a spinning disc reactor near the axis of the rotor was determined to be similar to plug flow, while more turbulent mixing is present near the outer rim of the disc (see Figure 1) [4,5].

Haseidl et al. (2016) [20] also modeled spinning disc reactors with residence time distribution using the same method as Visscher et al. [19] but added a convolution model, which is slightly more accurate. Hop et al. (2023) [8] used a combination of a non-ideal PFRs with one CSTR, which also gives a slightly better fit. Both models are harder to use because of the higher (mathematical) complexity than the simple and commonly used tanks-in-series model used by Visscher et al. [19] and De Beer et al. [4].

Figure 2 shows the combinations of a PFR and three CSTRs used to model the residence time distribution of the multiple rotor–stator spinning disc reactor by De Beer et al. [4] for different spinning rates, gap size, and flow rates. The y-axes of this figure explicitly show that the percentage of PFRs used in the model varies strongly with these variables. Since the relationship between these variables and the percentage of PFR is not fully understood, there also may be an influence of other dimensions of the used spinning disc reactor. This means that the insights from these graphs cannot be extrapolated to other spinning disc reactors or to other spinning rates or flow rates. For each reactor, the residence time distribution has to be determined experimentally to use these experimental data to fit a residence time distribution model.

To summarize, the general “engineering approach” to determining the residence time distribution of a spinning disc reactor is to do residence time experiments and fit the data to a combination of a certain percentage of PFR and a number of CSTRs in series. The presented models are not generally applicable.

2.5. Modeling of Gas–Liquid Reactions in a Spinning Disc Reactor

In the previous paragraphs, the literature results presented all applied to single-phase reactions. Some of them described liquid reactions between liquids that are not likely to mix on a micro scale but no gas–liquid reactions. In this paragraph, we discuss experimental results presented in the literature on gas–liquid flow and reactions in a spinning disc reactor.

Meeuwse et al. (2009) [21] presented a paper describing the flow in a closed single rotor–stator spinning disc configuration. The paper showed that gas hold-up increased and bubble size strongly decreased with increasing rotor speed. As a result, the volumetric mass transfer (expressed in cubic meter liquid per cubic meter reactor per second) initially declined a bit (due to decreasing gas surface compared to gas volume due to increasing gas hold-up), to strongly increase with further increase in rotor speed due to a strong decrease in bubble size.

The maximum rotor speed of the equipment used was 179 rad/s. It seems that at even higher rotor speeds, the volumetric gas transfer would have increased even further. However, the experiments described in this paper were conducted at two gas flows, a “low” gas flow of 1.8 × 10⁻⁶ m³/s and one high gas flow of 7.3 × 10⁻⁶ m³/s. The volumetric mass transfer corrected for gas hold-up (cubic meter liquid per cubic meter gas per second) for the “low” gas flow seems to become independent of rotor speed at approximately 80 rad/s, while the mass transfer at the high flow seems to increase until 150 rad/s (2009). This is attributed to the equilibrium between gas and liquid.

Meeuwse et al. (2012) [22] presented scale-up experiments showing similar flow patterns scaling from an enclosed rotor–stator to a multi-stage spinning disc reactor, comparing volumetric mass transfer results for one, two, and three rotor blades. This similarity was expected since he argued that the same combination of film and dispersed flow, as observed with a single enclosed rotor blade, would occur after each passage at the rotor due to the coalescence of the gas bubbles near the center of the reactor flowing out the reactor in a single stage configuration and into the next stage in case of a multi-stage configuration.

De Beer (2016) [12] conducted a similar test but at significantly higher gas volumetric flow rates; the gas volumetric flow rate was 40 times higher, and the ratio of gas to liquid was 100 times higher than reported by Meeuwse et al. (2009 and 2012) [21,22]. The tests were conducted in in a single-rotor spinning disc reactor; reactants were added from the top and removed from the bottom. At the top, he also found a bubble formation region, and below the rotor, a dispersed bubble region. The size of the bubbles in the dispersed region depends heavily on the gas flow, the Reynolds number, and, to a lesser extent, the liquid flow. In addition to the film and dispersion regime identified by Meeuwse et al. (2012) [22], a transition from a heterogenous to a homogenous dispersion regime was reported to occur at a gas low of 30 × 10⁻⁶ m³/s and liquid flow of 6 × 10⁻⁶ m³/s and Re = 2.4 × 10⁶. In the homogeneous regime, all bubbles are smaller than the height between the rotor and the stator, changing the behavior of the fluid mixture from a heterogeneous mixture to a pseudo fluid and is therefore identified as the optimal regime to operate a spinning disc reactor for multiphase reactions. The transition from heterogeneous to homogeneous is experimentally determined in two ways: visually and by means of pressure drop measurements (bubbles with a diameter > h cause pressure waves). Both measures are in reasonable agreement.

Furthermore, the decrease in bubble diameter is shown to be linearly proportional to the local energy dissipation rate [m²/s] ε_loc and a power function of the Reynolds number. The authors suggest this direct relation with the local energy dissipation rate suggests that the local energy dissipation closed to the rim of the reactor, i.e., where the film flow regime transits to the dispersion regime by shearing of gas bubbles, determines the maximum bubble diameter in the dispersion phase. Therefore, the local energy dissipation rate is expected to allow for scaling the current findings to scaled-up reactors and different reaction fluids.

2.6. Conclusions of the Literature Studies on Spinning Disc Reactors

Based on the literature findings, it is shown that the current way to determine residence time distribution in a spinning disc can be realized by CFD modeling [8] or by tailor-fitting experimental residence time distribution measurements to a combination of a PFR and CSTRs in series [4,5,19,20]. Both methods are too complicated for our aim.

To have optimal micro-mixing single-phase reactions in a spinning-disc reactor should be carried in the turbulent regime (ε_loc/υ > 1 × 10⁸ or Re > 1 × 10⁵).

Multiphase reaction in a spinning disc reactor should be operated in the homogenous regime, i.e., the shear should be so high that all gas bubble diameters are smaller than the height below the rotor. This can be determined experimentally since bubbles with a larger diameter cause pressure waves in the reactor.

3. Development of a General Model to Estimate Residence Time Distribution

The next step is to study the options at hand to simplify this approach and determine whether the inaccuracies that are inevitably caused by such simplifications are acceptable compared to the level of accuracy that normally is realized during feasibility studies.

Therefore, the question is how large is the effect of the inaccuracies that are induced when you do extrapolate these results to other spinning disc reactors, and how does that compare to the inaccuracies of other commonly used methods to determine the kinetic data for the use in feasibility studies?

3.1. Generalization of the Tanks-in-Series Model for a Spinning Disc Reactor

Section 2.4 establishes that the model by De Beer et al. (2014) [4] combines simplicity with good predictions, resulting in residence time distribution described by the combination of a variable percentage of PFR combined with three CSTRs in series. Therefore, we use this model as the basis for our approach.

We explored which variables prevent the current model from being applicable for the determination of kinetic parameters. We found two variables:

• The exact percentages of CSTRs and PFR;
• The sequence of CSTRs and a PFR in the model.

The exact percentages of CSTRs and PFR

Variables like spinning rate and flow rate have a considerable influence on the percentage of PFR, as modeled by De Beer et al. (2014) [4]. Their results show that the percentage of PFR varies from 0% to 40% combined with 3 CSTRs, depending on operation variables such as gap size, spinning rate, and flow rate. The question is how relevant the difference in percentage PFR is for the conversion of a reaction occurring in the rotor–stator spinning disc reactor. This then determines how large the inaccuracy is if we use a generalized model independent of variables such as gap size, spinning rate, and flow rate.

The sequence of CSTRs and PFR in the model

All residence time distribution-based models use a combination of PFR and CSTRs models [4,8,19,20]. According to the theory of reaction kinetics in reactors, conversion after a series of different reactors depends on the sequence in which different reactors are connected to each other [18]. This sequence does not influence residence time distribution and was therefore not included in the literature models. We, therefore, want to check how large the differences are when changing the sequence of the different reactors in the reaction kinetics and if, in our general model, a certain sequence must be used.

To determine the error induced by generalizing these variables, we adopted the following four-step approach:

3.

For several common reaction types, we calculated conversion in a PFR system and in three CSTRs in series, as is used in the model of De Beer et al. (2014) [4]. We then compare the conversion of both systems and draw conclusions on the differences.

4.

We calculated conversion for different sequences of the PFR and the CSTRs in the modeling of the conversion and the differences that are found in this way.

5.

We calculated conversion based on the percentage of the volume of the spinning disc reactor modeled as a PFR and the percentage of three CSTRs as reported by De Beer et al. (2014) [4], and conversion based on the combination of a PFR and 2–2.5× CSTRs as reported by Visscher et al. (2013) [19] to see how large differences in conversion are for the different variable values in these data sets.

6.

To determine whether the variations found in the three steps above are relevant compared to the variations that are normally caused by the experimental determination of kinetic values, we experimentally determined the conversion rate of acetic anhydride after hydrolysis at different temperatures in a spinning disc reactor and validated the model by comparison of our own dataset to conversion rates in literature data on this reaction. In addition, we compared the accuracy of our model to the accuracy of our own experimental data and the kinetic data from the literature.

We describe our approach to these four steps in the four following paragraphs. The outcomes of these steps are presented in Section 4.

3.2. Step 1: Modeling Reactions in a PFR and Three CSTRs in Series

The reaction rates and conversion were calculated in a PFR and in three CSTRs in series for different reaction types, according to the general theory on the reaction kinetics and reactors [18]. We chose the following reaction types to use in our model:

7.

A first-order reaction where the concentration of only one reactant is important (in a first-order manner).

8.

A second-order reaction occurs when one reactant reacts in a second-order manner, or two reactants with the same initial concentration react with each other.

9.

A general second-order reaction is where two reactants react with each other, and both concentrations are important.

10.

A second-order equilibrium reaction that reacts both ways until an equilibrium is reached (if the residence time is long enough).

11.

A second-order reaction, followed by a second-order reaction in which the desired product is involved as a reactant and therefore reacts away: a reaction in series.

12.

A second-order reaction where the limiting reactant also reacts to an unwanted side product: reactions in parallel.

The used reaction equations and reaction rates as a function of reactant concentration, are listed for the different types of reactions in Table 1. Reactions 3–6 are all described as second-order reactions but will also function as (pseudo) first-order reactions when the difference between the initial concentrations of the reactants is large. Conversions in PFR and in three CSTRs in series were calculated by choosing a wide range of values using a random generator for the reactant concentrations, the reaction rate constants, and the residence time.

3.3. Step 2: Effect of the Sequence of CSTRs and a PFR on the Conversion Rate

The effect of the sequence of CSTRs and a PFR was investigated by calculating conversion for a second-order reaction with a maximum difference between three CSTRs and a PFR; variable values are explained with the results (Section 4.1).

As described in the Introduction, a spinning disc reactor can be calculated as a combination of different percentages of a PFR and three CSTRs in series. Since the effect of different sequences and percentages of a PFR and three CSTRs were not clear, we assessed several percentages and sequences.

The sequences used were:

13.

First three CSTRs, then a PFR;

14.

First a PFR, then three CSTRs;

15.

One-disc model: a PFR, three CSTRs, a PFR;

16.

Three-disc model: a PFR, one CSTR, a PFR, one CSTR, a PFR, and one CSTR, a PFR.

The one-disc model and the three-disc model are based on the sequence in which PFR and CSTR are assumed to occur in a spinning disc with one or three discs according to Figure 1. In addition to the three different sequences the percentages of PFR and CSTR were also varied. Conversion was calculated for the different sequences for equal volumes of a PFR and three CSTRs (50%/50%), for 75% PFR/25% three CSTRs, and for 25% PFR/75% three CSTRs on a volume basis. A visual representation of the different sequences and varying volumes of a PFR and CSTRs is shown with the results (Section 4.2).

3.4. Step 3: Calculation of the Conversion from the Residence Time Distribution Data

De Beer et al. (2014) [4] and Visscher et al. (2013) [19] modeled residence time distribution in a spinning disc reactor with a combination of a PFR and CSTRs. We used these data to calculate reaction conversion based on the percentages of a PFR and the percentage and number of CSTRs. Calculations were done for a second-order reaction with variable values that gave the highest difference between a PFR and three CSTRs, as shown in the results (Section 4.1).

De Beer et al. (2014) [4] modeled a three-disc reactor with a variable percentage of a PFR combined with three CSTRs in series; the results of this are shown in Figure 2. For calculating conversion, we used the sequence of the three-disc model as described above (see Section 3.3). Visscher et al. (2013) [19] used a one-disc reactor and modeled this with a percentage PFR combined with a variable number of CSTRs. In their results, this number of CSTRs varied between 1.9 and 2.4. Because it is not possible to model conversion in part of a CSTR, conversion was calculated for one, two, and three CSTRs in series, and a weighted average of these values was used to calculate conversion for the reported amount of CSTRs in between 1.9 and 2.4. The one-disc sequence was used for the calculations for the combination of a PFR and CSTRs (see Section 3.3).

3.5. Step 4: Experimental Validation of the Model

To validate the model, a well-studied reaction was chosen to see if the conversions found in a spinning disc reactor were indeed between a PFR and three CSTRs in series. The reaction chosen is the hydrolysis of acetic anhydride. This reaction is chosen since this reaction is well documented, and the degree of conversion can be determined using aniline titration. Furthermore, no large deviations in reported kinetic data are expected between previously established reaction rates and the reaction rate determined in the spinning disc reactor, since no limitations in mass and heat transfer are expected when using low concentrations of acetic anhydride.

3.6. Literature Data on the Hydrolysis of Acetic Anhydride

The hydrolysis of acetic anhydride is a second-order reaction between acetic anhydride and water, having acetic acid as a result, as shown in Figure 3. Because water is also the solvent, it can also be described as a pseudo-first-order reaction, as the concentration of water hardly changes.

The hydrolysis reaction has been studied for years, and the kinetics are measured in separate ways and reported in the literature [23,24,25,26,27,28,29]. To compare the experimental data with all the literature found describing this reaction according to second-order kinetics is used. The Arrhenius parameters describing the kinetics at different temperatures, which are the activation energy (Ea) and the natural logarithm of the Arrhenius constant (A), are shown in

Table 2. Arrhenius parameters from the literature for the hydrolysis of acetic anhydride, calculated with second-order kinetics. The reaction rate constant can be calculated with $k=A·e^{-{\displaystyle \frac{E_a}{R·T}}}$ and this equation can be linearized with $\ln \left(k\right)=\ln \left(A\right)-{\displaystyle \frac{E_a}{R}}·{\displaystyle \frac{1}{T}}$ where k = reaction rate constant, A = Arrhenius constant, E_a = activation energy, R = gas constant, and T = temperature.

Reference	Ea (kJ/mol)	ln(A) (L/(mol/s))	Method
Eldridge and Piret (1950) [23]	43.2	7.5	Titration
Cleland and Wilhelm (1956) [24]	44.4	7.8	Titration
Glasser and Williams (1971) [25]	45.3	8.0	Calorimetry
Asprey et al. (1996) [26]	45.6 ± 0.6	7.7 ± 0.4	Calorimetry
Hirota et al. (2010) [27]	67.2 ± 1.4 *	15.6 ± 0.6 *	Calorimetry
Torraga et al. (2019) [28]	82.4 ± 12.4 *	22.4 ± 4.8 *	Calorimetry, NIR, and UV-vis
Garcia et al. (2021) [29]	64.4 ± 1.5	14.9 ± 0.6	Calorimetry
Average (±95% interval)	56.1 ± 36.3	12.0 ± 13.7

* Values displayed are averages ± SD of different values found by the authors using different methods.

.

The literature values are plotted linearized as ln k as a function of 1/T in Figure 4. We evaluated the reported values and found quite large deviations between experimentally determined values for the activation energy and the Arrhenius constant. We calculated the 95% confidence intervals; this results for the activation energy Ea in 56.1 ± 36.3 and for the ln(A) in 12.0 ± 13.7.

We validated the accuracy of our generalized model by plotting the measured relation between the reaction rate and temperature and comparing the reported values for activation energy Ea and ln(A) to the 95% confidence interval.

3.7. Experimental Validation of Acetic Anhydride Hydrolysis

A rotor–stator spinning disc reactor “Spinpro R10” supplied by Flowid, with a volume of 19 mL was used. This reactor has three rotors with a diameter of 70 mm in the stators with a diameter of 72 mm. The thickness of the two top rotors was 5 mm in a stator cavity of 6 mm; the third rotor had a thickness of 20.5 mm in a stator cavity of 21.5 mm because this rotor contained a magnet for propulsion. This means that the gap between the rotor and stator was 0.5 mm (G = 0.5/35 = 0.014), and the space at the rim of the rotor was 1 mm.

Acetic anhydride was added to the reactor at one of the entrances at the top with a syringe pump, and demineralized water was added to the other entrance at the top of the reactor with a small gear pump. This resulted in an initial concentration of acetic anhydride in the reaction mixture of 0.34 ± 0.03 M. The total liquid flow through the reactor was 11.7 ± 0.2 mL/min, resulting in a residence time of 100 ± 2 s inside the reactor itself. At the end of the reactor, the product exited the reactor through a metal tube with a residence time of 15 s before entering a sample tube. This sample tube contained 20 mL saturated aniline solution to stop the reaction immediately: aniline reacts very quickly with acetic anhydride, producing one molecule of acetic acid instead of two in the hydrolysis reaction, as is shown in Figure 3. A maximum of 6 mL sample per sample tube was collected to make sure that aniline was in excess. For every temperature, samples were taken in triplicate when a steady state was assumed to be reached (three residence times after the required temperature was reached). Temperatures of the reactor and the exit tube were monitored, and both were used for conversion calculations.

Experiments were conducted between 40 and 80 °C at 2000 rpm, resulting in a Reynolds and C_w number varying between respectively 3.5 × 10⁵ < Re < 7.5 × 10⁵ and 8 < C_w < 16. Above 80 °C, conversion after total residence time was so high (>90%) that a reliable determination of the reaction rate constant was not possible because minor differences in conversion have a large effect on the outcome. Below 40 °C, the measured conversion was significantly higher than the literature values, possibly because of the bad solubility of acetic anhydride at low temperatures. In that case, our assumption that the reaction is not mass-transfer limited is not valid anymore and, therefore, not suitable for comparison to literature values.

Samples were titrated with 1.0 M NaOH to determine the concentration of the product acetic acid. For every temperature, a sample without aniline was taken and titrated to calculate the initial acetic anhydride concentration. With the known initial concentration of acetic anhydride, the conversion of the hydrolysis reaction can be calculated from the final concentration of acetic acid.

4. Results

In the previous chapter, we described our aim to generalize existing models describing residence time distributions in spinning disc reactors. Since this approach inevitably will induce inaccuracies compared to the currently used practice, we developed a four-step approach (Section 3.2, Section 3.3, Section 3.4 and Section 3.5) to determine the overall error induced by this method on the accuracy of kinetic parameter determination. The results of these four steps are presented in the following paragraphs.

4.1. Modeling Reactions in PFR and 3 CSTRs in Series

Since we want to determine the effect of simplifying residence time distribution from a certain percentage of a PFR in combination with three CSTRs in series, we calculate the conversion for 100% PFR and 100% three CSTRs in series. Conversion was calculated for a wide range of reaction types (see reactions 1–6 in Table 1) using the same (randomly chosen) variable values for start concentrations, reaction constants, and residence time. The same randomly chosen values were used for calculating conversion in 100% PFR and 100% three CSTRs. For each reaction, 500 data points were calculated. The difference in calculated conversion between three CSTRs compared to the calculated conversion in PFR is plotted against the calculated conversion in a PFR in Figure 5. For reactions 1 (first-order reaction) and 2 (equimolar second-order reaction), the difference in calculated conversions found formed a continuous line; for the other reactions, the results are shown as a cloud of points. All data points and calculation methods are shown in Table S1 in the Supplementary Materials.

The calculated conversion rates for reaction 3 (general second-order reaction) lay in between the calculated conversion rates for reactions 1 and 2. This is to be expected because it is a second-order reaction that can act as a pseudo-first-order reaction when one of the concentrations of the reactants is much higher than the concentration of the other reactant. For the equilibrium reaction (no. 4), differences are equal to reaction 3 or lower because close to equilibrium, there is no difference between CSTRs and PFR. For reactions 1–4, PFR always gives the highest conversion.

For the reaction in series (no. 5), three CSTRs shows higher conversions when conversion is low. This is because the unwanted side reaction that leads to low conversion is promoted in a PFR. For the parallel reaction (no. 6), PFR has higher conversions, but in some cases, the side reaction is promoted in PFR, leading to lower conversion.

For a certain reaction, the reaction rate constants are assumed to be a function of temperature only, while residence time and initial concentrations may vary (as shown in Section 2.4 the assumption that reaction rate is a function of temperature alone requires sufficient micro-mixing; see also the Discussion Section).

Since all variables were chosen randomly not all shown data points are in line with realistic production circumstances. For commercial production only the conversion of over 80% is relevant.

The highest difference found between a PFR and three CSTRs in series is indicated with a vertical arrow in Figure 5. The reaction variables at maximum conversion were: k = 0.88 min⁻¹ mol⁻¹ m³, t = 2.21 min, C_A0 = 0.88 mol m⁻³, and C_B0 = 1.61 mol m⁻³ for a second-order reaction (reaction 3 in Table 1). With these values of the variables, a conversion of about 88% is reached in PFR, while in three CSTRs only, 79% conversion is reached. This means that the conversion in three CSTRs is 10% lower than in a PFR, and that a mistake of a maximum of 10% is made when conversion is calculated in three CSTRs when a PFR comes closer to reality, and vice versa. A generalized model using either three CSTRs or PFR for the calculation of reaction kinetics will, therefore, have a maximum inaccuracy of 10%, provided that the precondition of sufficient micro-mixing is met.

4.2. Effect of Sequence of CSTR and PFR on Conversion Rate

As described by Levenspiel (1999) [18], there will be a difference in calculated conversion when you change the sequence of a a PFR and CSTRs that are put in series. So, in principle, it does matter if you use a PFR first and then a CSTR or the other way around. Therefore, we compare conversion rates found for different sequences and percentages of CSTR and PFR (see

Table 3. Calculated conversion for the combination of three CSTRs and a PFR in different sequences. For the initial concentrations, the reaction rate constant and residence time, reaction variables were used that led to the highest difference in conversion between PFR and CSTR, as listed in the caption of Figure 5. These conversions are 0.876 for a single PFR and 0.790 for three CSTRs. Colored bars are a visualization of the reactors used in the calculations, PFR in blue and CSTR in orange. The maximum error is equal to the largest difference between two sequences divided by the average value.

	25% CSTR, 75% PFR	50% CSTR, 50% PFR	75% CSTR, 25% PFR
PFR → CSTR	0.872	0.859	0.835
CSTR→PFR	0.868	0.851	0.825
One-disc model	0.870	0.856	0.831
Three-disc model	0.869	0.854	0.829
Average ± stdev	0.870 ± 0.002	0.855 ± 0.003	0.830 ± 0.004
Maximum error	0.4%	1.0%	1.2%

).

The calculated conversions in Table 3 show that the effect of the sequence for a given percentage of the CSTRs and a PFR is small, 1.2% or less of the found average conversion. Therefore, we can conclude that in practice, the sequence of the CSTRs and a PFR in this model is not important because the differences can be neglected compared to the influence of the percentages of three CSTRs and a PFR, which was shown to be up to 10% of conversion.

4.3. Literature Data on Residence Time Distribution Versus Expected Conversion

The results on residence time distribution of De Beer et al. (2014) [4], previously shown in Figure 2, were converted to conversion by combining the percentage of a PFR with three CSTRs assuming a second-order reaction and the initial concentration, reaction rate constant and residence time as listed in the caption of Figure 5. The results are shown in Figure 6.

The conversion data in Figure 6 show the same pattern as the original data in Figure 2. This is to be expected because conversion is higher in PFR, and the data will, therefore, show a higher conversion when the percentage of the PFR is higher. Lines in the graph indicate the conversion for 100% PFR and 100% three CSTRs. Calculated conversion varies between 0.79 and 0.85, which is within 7% of the average of 0.82.

Visscher et al. (2013) [19] also modeled the residence time as a PFR and CSTRs in a spinning disc reactor, in this case, a spinning disc reactor with only one disc. Results of the calculations of conversions from the percentage of a PFR and percentage and number of CSTRs are shown in Figure 7. Figure 7 shows that also in the case of a single disc spinning disc reactor, predicted conversions are quite close to each other and are between values for a PFR and three CSTRs, although 1.9–2.4× CSTRs were used in the calculations. Values for conversion vary between 0.78 and 0.86, which is within 10% of the average value of 0.82. Only for the gap of 3 mm (G = 0.023), the conversion drops below the three CSTRs line. This shows that the boundaries of three CSTRs and a PFR also hold for a one-disc spinning disc reactor when the gap size is small.

So while the results of De Beer et al. (2014) [4] show a range of PFR volume of 0–40% combined with three CSTRs, and the results of Visscher et al. (2013) [19] show a range of 7–68% for a PFR volume combined with 1.9–2.4 CSTRs, the calculated conversions for all these data points are within 10% for the used reaction kinetics. This shows that the mixing and flow regimes only account for 10% of the variation in conversion calculated with these models if the precondition of sufficient micro-mixing is met.

4.4. Validation of the Model by Acetic Anhydride Hydrolysis

To validate the model, a well-studied reaction was chosen to see if the conversions found in a spinning disc reactor were indeed between PFR and three CSTRs. The hydrolysis of acetic anhydride has been studied extensively, and kinetic data from the literature are presented in Table 2 and Figure 4. At the temperatures we conducted this reaction, there were no mass- or heat-transfer issues reported. There is no reason to expect different values in a spinning disc reactor than in any other laboratory set up to measure kinetic properties. Figure 8 shows the conversions found in the spinning disc reactor at different temperatures after a total residence time of 115 ± 2 s and an ingoing concentration of 0.33 ± 0.01 M. Raw data of this experiment can be found in Table S2 in the Supplementary Data. The temperature on the x-axis is the temperature in the reactor itself, but for the calculations, the lower temperature in the exiting tube was considered.

For all literature data, predicted conversion in a PFR (solid line) and in three CSTRs (dotted line) were calculated using the kinetic parameters reported in the literature. This shows that for the reported literature data, the difference between predicted conversion in a PFR is similar to the prediction of the conversion in three CSTRs. The difference between the different literature data is much larger than the difference between a PFR and three CSTRs. The conversion determined in our experiment was within the range of data reported in the literature.

From our experimental data, we calculated the reaction constant (k) for every temperature, both for the assumption that the reactor can be described as a PFR and as three CSTRs in series. This results in two values for the reaction rate constant per temperature used in the experiments. These values are shown in Figure 9, together with the literature data shown before in Figure 4. Figure 9 shows that the difference in the outcome of three CSTRs or PFR described a negligible reaction at low temperatures, while the difference increases with increasing temperature. This is to be expected since the value of k calculated with three CSTRs increases stronger with temperature than calculated with a PFR.

For both calculation methods, the Arrhenius constants were calculated and are shown in

Table 4. Arrhenius parameters calculated from the reaction rate constants in Figure 9. Calculated parameter values are ±95% confidence interval.

Calculation Method	Activation Energy Ea (kJ/mol)	ln(A) (L/(mol/s))
Three CSTRs	41.3 ± 17.5	6.3 ± 6.3
PFR	35.2 ± 14.6	3.9 ± 5.3
Average	38.3 ± 10.3	5.1 ± 3.7
Average of the literature values (Table 2)	56.1 ± 36.3	12.0 ± 13.7

. The values are different but vary significantly less than the difference between the other values reported in the literature. This shows that the error induced by the simplifications in the residence time distribution model is relatively minor compared to the error in the experimental determination of the reaction rate constant as a function of temperature. When fitting the experimental data to the linearized Arrhenius equation, both ways to describe the residence time distribution (using either three CSTRs or a PFR) give a relation for k with values for the activation energy Ea and ln(A) that fall within the 95% confidence interval of the literature values.

5. Discussion

Before we continue to discuss these findings, we must stress that this model aims to identify those reactions that in the current industrial settings seem slow but that can be sped up considerably using process intensification principles. The intended outcome is to make it possible to raise some awareness of what a mini reactor could mean, not only in terms of smaller reactor volume but also in terms of decreased separation requirements and improved heat integration potential. This method is certainly not fit for reactions that are known to be fast and of which one wants to accurately determine how fast and selective they can be made, in that the other methods described in the literature study are more suitable.

Going back to our results, we found that the error induced by this approach is 10%. We also showed that the variance in reported values in the literature is much larger than this 10%, which makes it sufficiently accurate for a (pre-)feasibility study. However, both findings are only valid if the used equations apply. Therefore, we will have a more careful look at these preconditions:

• The residence time distributions described in Section 2.4 are determined on a single-phase fluid;
• The value of the reaction rate constant is a function of temperature only; thus, no limitations in mass- and or heat-transfer.

In the case of a single-phase flow, the boundary conditions are met if limitations in the heat- and mass-transfer become negligible by micro-mixing. As shown by Manzano Martinez (2017 and 2022) [10,11], micro-mixing occurs when sufficient shear is applied (ε_loc/v > 1 × 10⁸ of Re > 1 × 10⁵). The question remains, under what conditions applies for a model derived for the single-phase flow to multiphase flow? The issue here is the bubble size. If the bubble size is so small that diffusion out of the gas bubble requires significantly less time than the residence time in the dispersed phase, the gas–liquid mixture may be treated as a pseudo-single phase. The work by De Beer et al. (2016) [12] shows that with sufficient shear, all bubble diameters become smaller than the height between the rotor and the stator. The combination of gas flow and Reynolds number causing the transition from heterogeneous to homogeneous flow can be experimentally determined by measuring the pressure waves caused by gas bubbles with a diameter higher than this height.

6. Conclusions

We developed a model describing residence time distribution in spinning disc reactors that is generally applicable to all spinning disc reactors, allowing us to determine the kinetic parameters of reactions without first experimentally determining the residence time distribution of the specific spinning disc reactor that is used for the experimental determination of the kinematic parameters.

To meet this goal, we simplified the available models describing residence time distribution in such reactors. We determined what the maximum error in reaction rate would be if we assumed 100% PFR or 100% three CSTRs in series instead of a tailor-made description in the percentage of PFR and three CSTRs. We found the accuracy for the prediction of conversion to be 90% or more compared to the more accurate models. Provided that the spinning disc is operated at the optimal settings, i.e., turbulent flow in case of single-phase flow [11], or heterogenous regime in case of multiphase flow [12].

Validation with the hydrolysis of acetic anhydride shows that the 10% inaccuracy falls within the accuracy of the value of the experimental determination of the reaction rate constant.