Open Access Article
This Open Access Article is licensed under a
Creative Commons Attribution 3.0 Unported Licence

Rapid multistep kinetic model generation from transient flow data

Christopher A. Hone a, Nicholas Holmes a, Geoffrey R. Akien ab, Richard A. Bourne *a and Frans L. Muller *a
aInstitute of Process Research and Development, School of Chemistry and School of Chemical and Process Engineering, University of Leeds, LS2 9JT, UK. E-mail: R.A.Bourne@leeds.ac.uk; F.L.Muller@leeds.ac.uk
bDepartment of Chemistry, Lancaster University, Lancaster, LA1 4YB, UK

Received 26th May 2016 , Accepted 23rd September 2016

First published on the web 3rd October 2016


Today, the generation of kinetic models is still seen as a resource intensive and specialised activity. We report an efficient method of generating reaction profiles from transient flows using a state-of-the-art continuous-flow platform. Experimental data for multistep aromatic nucleophilic substitution reactions are collected from an automated linear gradient flow ramp with online HPLC at the reactor outlet. Using this approach, we generated 16 profiles, at 3 different inlet concentrations and 4 temperatures, in less than 3 hours run time. The kinetic parameters, 4 rate constants and 4 activation energies were fitted with less than 4% uncertainty. We derived an expression for the error in the observed rate constants due to dispersion and showed that such error is 5% or lower. The large range of operational conditions prevented the need to isolate individual reaction steps. Our approach enables early identification of the sensitivity of product quality to parameter changes and early use of unit operation models to identify optimal process-equipment combinations in silico, greatly reducing scale up risks.


Over the past decade flow technologies have become established for the discovery and manufacture of active pharmaceutical ingredients.1,2 A plethora of synthetic organic methods have been developed using flow reactors.3,4 However the development of new strategies that give fundamental understanding of complex reaction systems using flow reactors has not received as much attention.5,6 In this communication, we demonstrate the rapid collection of experimental data using an automated flow reactor for the generation of kinetic models.

Reaction profiling for reaction kinetics is easily achieved under batch conditions due to the ability to collect multiple time points within a single experiment, often with only milligrams of material used.7 Flow reactors give precise control over the reaction parameters. The key benefit of flow is that the system reaches steady-state thus providing a consistent output quality. There are many examples of kinetic investigations in small-scale flow systems.8 For instance Reizman et al. used steady state flow measurements to fit four rate constants and four activation energies from a small number of flow profiles.9 The authors concluded that high confidence parameter estimates require isolation of the individual pathways by measurements starting from each intermediate.

Kinetic profiling in flow systems suffers from three significant issues: (i) significant time is required to get to steady state; (ii) flow experiments require unnecessary material usage, approx. 1.5 reactor volumes per measurement, due to the transient period prior to reaching steady-state,10 and (iii) dispersion in flow systems can influence the outlet concentration resulting in errors in the derived rate constants. In this paper we address these three issues by use of online measurements in a system with transient flow.

Several authors have generated kinetic data based on transient flows. Mozharov et al. applied a step change in flow rate to study a Knoevenagel condensation reaction using inline non-invasive Raman spectroscopy.11 The response analysis was hindered by convolution of the step change by ‘real system’ response times, and the approach is infeasible for analysis techniques with longer acquisition times, such as HPLC and GC. Moore et al. studied a Paal–Knorr pyrrole synthesis using an exponential flow ramp with online infrared spectroscopy.12 The effect of residence time and temperature were investigated, at constant inlet concentrations. The experimental data were successfully fitted to a priori rate expressions. Schaber et al. extended the approach in terms of model discrimination.13

This communication demonstrates a fully automated approach for collecting experimental data for kinetic model generation using transient flow data. The approach is illustrated using a nucleophilic aromatic substitution (SNAr) of 2,4-difluoronitrobenzene 1 with pyrrolidine 2 in ethanol to give a mixture: desired product ortho-substituted 3, para-substituted 4 and bis-adduct 5 as side products (Scheme 1).


image file: c6re00109b-s1.tif
Scheme 1 SNAr reaction of 2,4-difluoronitrobenzene 1 with pyrrolidine 2.

In this study an automated mesoscale flow reactor system (5 mL reactor volume, 0.79 mm internal diameter, 1.58 mm outer diameter) with online HPLC, see Fig. 1, was used. A linear gradient flow ramp, see Fig. 2, which allows for investigating a complete reaction profile from a single transient experiment was developed. A transient profile measurement starts by setting the maximum flow for all pumps. The flow rate ratio of P3 (pyrrolidine in EtOH) and P2 (EtOH) is varied at the beginning of each ramp to obtain different molar equivalents (1.5, 4 and 7) of pyrrolidine 2. As a sudden increase in flow rate produces unstable flow, pump flow rates were kept constant long enough to establish steady-state. Subsequently the flow rate, Q, is slowed at a constant rate, α (0.836 mL s−2), so as to increase the residence time (τres) over time t, whilst maintaining constant inlet concentrations. Reactor effluent is injected into the HPLC at 2 min time intervals. The pump flow ramps, reactor temperature and sample loop injection were automatically controlled by a MATLAB based computer program. Rather than isolate the individual pathways we extended the range of reaction conditions so as to over- and under-react significantly. Application of the full operational range allowed by the equipment e.g. from the mildest (e.g. dilute, low temperature) to the harshest (e.g. concentrated, high reagent to substrate ratios, high temperature) results in higher concentrations of intermediates and by-products. This represents a richer dataset that increases the confidence in parameters of the fitted kinetic motifs.14 Hessel et al. introduced the concept of novel process windows (high temperature, high pressure and high concentration) for accessing conditions in flow which are not typically accessed in conventional practice for organic synthesis in batch.15 Thus linear gradient flow ramps (Fig. 2) were developed to explore from the mildest (2 = 1.5 mol eq., = 0.5 min and T = 30 °C) to the harshest conditions (2 = 7 mol eq., = 2 min and T = 120 °C). Residence time points less than 0.5 min were not collected as transient effects dominate.16 The approach focuses on using reaction conditions that maximise the confidence in the kinetic parameter estimates.


image file: c6re00109b-f1.tif
Fig. 1 Automated continuous-flow reactor system equipped with software for controlled flow ramps. The reactor comprises of three Jasco PU-980 dual piston HPLC pumps which feed into Swagelok tee-pieces for mixing. A polar bear flow synthesiser (Cambridge Reactor Design, UK) was used for heating and cooling PTFE reactor tubing (5 mL internal volume, 0.79 mm internal diameter). A VICI Valco internal sample injector (SL) extracted aliquots of neat reaction for HPLC analysis. The composition of the reactor outlet was determined by online HPLC. Pressure control was achieved with a back-pressure regulator at the outlet of the system. The internal reactor temperature was measured by a thermocouple inserted into the centre of the reactor.

image file: c6re00109b-f2.tif
Fig. 2 Changes in volumetric flow rate over time, where QP1, QP2, QP3 and Qtotal were for image file: c6re00109b-u1.tif pump 1 (Ar in EtOH), image file: c6re00109b-u2.tif pump 2 (EtOH), image file: c6re00109b-u3.tif pump 3 (pyrrolidine in EtOH), image file: c6re00109b-u4.tif total volumetric flow rate, and ◆ HPLC injection respectively. The linear flow ramps correspond to pyrrolidine 2 to 2,4-difluoronitrobenzene 1 molar ratios: (i) 1.5[thin space (1/6-em)]:[thin space (1/6-em)]1, (ii) 4[thin space (1/6-em)]:[thin space (1/6-em)]1 and (iii) 7[thin space (1/6-em)]:[thin space (1/6-em)]1 using Qtotal from 10 to 1.5 mL min−1.

The transformation of the data is shown in Fig. 3. For a linear flow ramp the velocity in the coil is described by:

 
u(t) = β(T)(uoαt)(1)
where uo is the initial velocity, α the rate of change and β(T) = 1 + αv(T1T0) is the correction for the thermal expansion of EtOH due to the temperature rise from ambient to the reactor temperature. A fluid element leaving the coil at time t will have travelled the length of the tube (L) and resided a total time τres so:
 
image file: c6re00109b-t1.tif(2)


image file: c6re00109b-f3.tif
Fig. 3 i) Total volumetric flow rate (Qtotal) change throughout the duration of a ramp, ♦ HPLC injection; (ii) residence time, τ, as a function of operation time, t; (iii) resulting concentration-residence time profile generated from a typical multi-step reaction scenario, ■ 2,4-DF 1, ● ortho-3, ● para-4, ▼ bis-5.

Integration leads to a second order equation from which the residence time τres of a sample taken at time t can be worked out as:

 
image file: c6re00109b-t2.tif(3)

The method was applied to generate a series of concentration-time profiles at four temperature levels (30, 60, 90 and 120 °C) giving a total of 12 profiles (Fig. 5). The 72 data points shown in Fig. 5 were collected in less than 3 hours collection time. The reaction conditions were selected to give a wide range of conversions values for the different reaction components.

The data were then fitted to the kinetic scheme shown in Scheme 1, with different orders assessed with respect to the aromatic components and pyrrolidine (see ESI). The kinetic motif in which all steps are second order gave the best fit. The rate constants were initially fitted at isothermal conditions (90 °C) using the Levenburg–Marquardt algorithm, a non-linear least squares algorithm, in DynoChem software (Scale-up Systems). Subsequently all the experimental data (72 data points) were simultaneously fitted to give all the kinetic parameters (Table 1). The kinetic model very closely corresponded to the experimental data, with R2 = 0.9995 (Fig. 4). The rate constant k1 for ortho-3 is 20 times larger than k2 of the para-4 product formation and their activation energies are 33.3 ± 0.3 kJ mol−1 and 35.3 ± 0.5 kJ mol−1 respectively; thus temperature influences the rate, but not selectivity. Parameter uncertainties were all less than 4%. Even the rate parameters for the overreaction pathways identified with minimal uncertainty. Unsteady-state results were compared to results at steady-state, the two methods gave statistically similar results to data collected at steady-state conditions (ESI). More aggressive conditions e.g. high molar equivalents of pyrrolidine and high temperature result in elevated formation of the bis adduct 5 product to which k3 and k4 were fitted with high confidence, thus preventing the need to synthesise and isolate each reaction component. Our approach offers significant time savings and minimises material consumption compared to a steady-state approach, addressing shortcomings (i) time required to reached steady-state for each measurement and (ii) material wastage reaching steady-state. Further material savings could be made using a flow system with a small internal volume, such as on a microliter scale as reported by McMullen et al.17

Table 1 Kinetic parameter estimates and standard errors (SE) from the fitting based on 95% confidence level. Rate constants, k, are given at Tref = 90 °C
  k ± SE (10−2 M−1 s−1) E a ± SE (kJ mol−1)
Step 1 57.9 ± 0.7 33.3 ± 0.3
Step 2 2.70 ± 0.06 35.3 ± 0.5
Step 3 0.865 ± 0.004 38.9 ± 1.5
Step 4 1.63 ± 0.11 44.8 ± 1.8



image file: c6re00109b-f4.tif
Fig. 4 The dispersion ratio κ vs. DeSc0.5 for various literature data. Where possible data have or correlations have been taken from the papers. If the data could not be retrieved and no correlation was given the data has been fitted with a line that visually represents the data well.

Limitation (iii) concerned the influence of dispersion on the measured rate constants. If the coil reactor is described by a plug flow model the conversion for a first order reaction is:

 
image file: c6re00109b-t3.tif(4)

It is well established that dispersion can have a significant effect on the progress of a reaction.18 In the 1950s Taylor described dispersion in of a solute a straight pipe with laminar parabolic flow and this was modified by Aris,19,20 to give:

 
image file: c6re00109b-t4.tif(5)

With DS the dispersion coefficient, Dm the diffusion coefficient, dt the tube diameter and u the mean velocity in the tube. In liquid systems the term Dm is negligible, and κ represents the ratio of dispersion in a conduit to dispersion in a straight tube with similar diameter; for a straight cylindrical tube under laminar flow κ = 1. The dispersion ratio κ is dependent on the channel geometry and the flow regime (e.g. laminar or turbulent). Laminar flow in coiled tubes deviates from the parabolic velocity profile as a result from centrifugal forces. So called Dean vortices form which introduce radial flow that reduces dispersion. The intensity of the radial flow is characterised by the Dean number:

 
image file: c6re00109b-t5.tif(6)

Re is the Reynolds number (ρudt/μ) and dc the coil diameter (0.79 mm in our case). In the 1970s dispersion in coils under laminar conditions was studied experimentally in the group of Vasudeva in wide bore tube (4.4–20 mm),21–23 and by Van den Berg and Kockmann in 0.5 and 1 mm ID tubes.24,25 They found that the dispersion ratio κ reduces significantly below 1 (Fig. 4). Theoretical work by Janssen and later by Johnson showed dispersion ratio in coils may be correlated by DeSc0.5, where the Schmidt number (Sc) is defined as μ/ρDm.26,27 In the 1960s the effect of dispersion was coupled to a reaction system by Wehner28 to give a general equation for the effect of dispersion on the conversion of a solute due to a first order reaction with rate constant k:

 
image file: c6re00109b-t6.tif(7)

In which dispersion is characterised with the dimensionless Péclet number Pe = D/uL. Eqn (7) can be simplified for values of Pe < 0.05 by applying a second order Taylor series expansion for a:

 
a ≈ 1 + 2resPe − 2(resPe)2(8)
where resPe = kD/u2 is assumed to be small, therefore the dispersion is small and peaks remain symmetrical. Combining eqn (5) and (7), and the observation that (i) in liquid systems the term Dm becomes negligible compared to Ds and (ii) τres = L/u, it follows that:
 
image file: c6re00109b-t7.tif(9)

This shows that in a coil the observed rate constant kobs may be given as:

 
image file: c6re00109b-t8.tif(10)
where the Damköhler number for radial diffusion may be defined as Dar = kdt2/Dm. The deviation εDC in the observed rate constant kobs obtained from the profiles measured in continuous flow in a coil with respect to the true rate constant k will be:
 
image file: c6re00109b-t9.tif(11)

Eqn (11) shows the kinetic constants will always be underestimated if dispersion is significant. The above result shows that for a tube with 1 mm internal diameter and a typical value for Dm = 0.8 × 10−9 m2 s−1 in EtOH (ref. 29) the reduction in rate constant is εDC ≈ −500% κk. For instance, Durant measured the rate constant for the thermolysis of 1,3-dioxin-4-ones to be <0.015 s−1 in a 1 mm ID coil at DeSc0.5 = 13.8, 138 and 1380.30 This corresponds to a deviation εDC ≤ −5.6%, −1.5%, 0.38% respectively (κ from Fig. 4), which corresponds well to their observation of a negligible effect of dispersion. In a 1 mm coil system with κ < 0.1 rates as fast as 0.1 s−1 can be measured accurately. Smaller tube diameters will reduce the error, as κDar reduces.

To assess the impact of dispersion on the second order rate constants for conversion of 2,4-difluoronitrobenzene 1, at a particular excess of pyrrolidine the first order rate constant k may be approximated by:

 
k ≈ (k1 + k2)C2,0(12)

In our system we measured an F curve at 6 min residence time giving Dar = 788 and DeSc0.5 = 65 (ESI). This resulted in a dispersion ratio of κ = 0.31 which corresponds well with the data from Van den Berg.24 Using the correlation from Van den Berg for κ in coiled tubes we estimate that for the experimental conditions used to generate the data on which Table 1 is based κ, is in the range of 0.05 to 0.15 (Fig. 4). Using eqn (11) the extent to which k1 + k2 is underestimated can be evaluated for profiles where the starting material concentration deviates significantly from zero e.g. plots i, ii and v for 1 in Fig. 5 (for the profiles where the reactant has fully converted or no stating material is available at all residence times no sensible estimate dispersion effect on k can be made). To be conservative we assumed a value of κ = 0.15 for all errors in Table 2. The unshaded cells represent profiles with significant concentrations of reacting components. The fitted motif is thus fit for purpose with respect to evaluate parameter sensitivity, alternate unit operations and scale-up.


image file: c6re00109b-f5.tif
Fig. 5 Concentration-time profiles from simultaneous parameter fitting, points = experiments ■ 2,4-DF 1, ● ortho-3, ● para-4, ▼ bis-5, lines = model using Table 1 kinetic parameter estimates (i–xii).
Table 2 Estimated predicted % effect in rate constants k1 + k2 and k3 + k4 from dispersion. The cells highlighted in grey show where the deviation in concentration is too low to significantly affect the fitting of the rate constants

image file: c6re00109b-u5.tif

Conclusions

An automated continuous-flow platform with quantifiable online analysis has been developed as an enabling tool for the rapid and economic collection of kinetic profiles. In less than 3 hours, 12 reaction profiles were collected; sufficient to fit a kinetic motif consisting of 4 reactions, 8 fitting parameters, with less than 4% uncertainty. Furthermore, it has been demonstrated that the effect of dispersion in these system results in an underestimation of the rate constants by 5% or less. The efficiency of the linear flow ramp is reduced if the product analysis takes a long time, however; with recent advances in analytical techniques, e.g. UPLC and FlowIR many reaction systems can be analysed successfully within seconds to minute timescales. Data collected from transient flow profiles gave statistically similar results to data collected at steady-state further validating the flow ramp approach. The combination of a linear gradient flow ramp and extreme conditions far away from the preferred operating point, easily accessed in continuous-flow reactors, enables rapid data generation with a quality suitable for fitting parameters of multistep kinetic motifs. Our approach allows kinetic models to be generated much earlier in process development, allowing early estimation of the sensitivity of product quality to input parameter changes. The model can be applied in silico to simulate alternative optimisation scenarios, equipment configurations, and to achieve significant reductions in scale up risks and costs.

Acknowledgements

The authors gratefully acknowledge the generous support by AstraZeneca, EPSRC DTG funding and the University of Leeds. RAB thanks the Royal Academy of Engineering for an Industrial Secondment at AstraZeneca.

References

  1. A. Adamo, R. L. Beingessner, M. Behnam, J. Chen, T. F. Jamison, K. F. Jensen, J.-C. M. Monbaiu, A. S. Myerson, E. Revalor, D. R. Snead, T. Stelzer, N. Weeranoppanant, S. Y. Wong and P. Zhang, Science, 2016, 352, 61–67 Search PubMed.
  2. B. Gutmann, D. Cantillo and C. O. Kappe, Angew. Chem., Int. Ed., 2015, 54, 6688–6728 Search PubMed.
  3. L. Malet-Sanz and F. Susanne, J. Med. Chem., 2012, 55, 4062–4098 Search PubMed.
  4. N. G. Anderson, Org. Process Res. Dev., 2012, 16, 852–869 Search PubMed.
  5. R. L. Hartman, J. P. McMullen and K. F. Jensen, Angew. Chem., Int. Ed., 2011, 50, 7502–7519 Search PubMed.
  6. P. Poechlauer, J. Colberg, E. Fisher, M. Jansen, M. D. Johnson, S. G. Koenig, M. Lawler, T. Laporte, J. Manley, B. Martin and A. O'Kearney-McMullan, Org. Process Res. Dev., 2013, 17, 1472–1478 CrossRef CAS.
  7. D. G. Blackmond, Angew. Chem., Int. Ed., 2005, 44, 4302–4320 Search PubMed.
  8. For selected kinetic examples see; (a) J. Weon, Z. Horváth, A. G. O'Brien, P. H. Seeberger and A. Seidel-Morgenstern, Chem. Eng. J., 2014, 251, 355–370 Search PubMed; (b) D. M. Roberge, C. Noti, E. Erle, M. Eyholzer, B. Rittiner, G. Penn, G. Sedelmeier and B. Schenkel, J. Flow Chem., 2014, 4, 26–34 Search PubMed; (c) H. Song and R. F. Ismagilov, J. Am. Chem. Soc., 2003, 125, 14613–14619 Search PubMed; (d) A. Gholamipour-Shirazi and C. Rolando, Org. Process Res. Dev., 2012, 16, 811–818 Search PubMed.
  9. B. J. Reizman and K. F. Jensen, Org. Process Res. Dev., 2012, 16, 1770–1782 Search PubMed.
  10. F. E. Valera, M. Quaranta, A. Moran, J. Blacker, A. Armstrong, J. T. Cabral and D. G. Blackmond, Angew. Chem., Int. Ed., 2010, 49, 2478–2485 Search PubMed.
  11. S. Mozharov, A. Nordon, D. Littlejohn, C. Wiles, P. Watts, P. Dallin and J. M. Girkin, J. Am. Chem. Soc., 2011, 133, 3601–3608 Search PubMed.
  12. J. S. Moore and K. F. Jensen, Angew. Chem., Int. Ed., 2014, 53, 470–473 Search PubMed.
  13. S. D. Schaber, S. C. Born, K. F. Jensen and P. I. Barton, Org. Process Res. Dev., 2014, 18, 1461–1467 Search PubMed.
  14. G. E. P. Box and W. G. Hunter, Technometrics, 1965, 7, 23–42 CrossRef.
  15. V. Hessel, D. Kralisch, N. Kockmann, T. Noël and Q. Wang, ChemSusChem, 2013, 6, 746–789 Search PubMed.
  16. K. B. Nagy, B. Shen, T. F. Jamison and K. F. Jensen, Org. Process Res. Dev., 2012, 16, 976–981 Search PubMed.
  17. J. P. McMullen and K. F. Jensen, Annu. Rev. Anal. Chem., 2010, 3, 19–42 Search PubMed.
  18. C. B. Rosas, Ind. Eng. Chem. Fundam., 1969, 8, 361–364 Search PubMed.
  19. G. Taylor, Proc. R. Soc. London, Ser. A, 1953, 219, 186–203 Search PubMed.
  20. R. Aris, Proc. R. Soc. London, Ser. A, 1956, 235, 67–77 Search PubMed.
  21. R. N. Trivedi and K. Vasudeva, Chem. Eng. Sci., 1975, 30, 317–325 Search PubMed.
  22. K. D. P. Nigam and K. Vasudeva, Chem. Eng. Sci., 1976, 31, 835–837 Search PubMed.
  23. V. D. Shetty and K. Vasudeva, Chem. Eng. Sci., 1977, 32, 782–783 Search PubMed.
  24. J. H. Van den Berg and R. S. Deelder, Chem. Eng. Sci., 1979, 34, 1345–1347 Search PubMed.
  25. M. G. Gelhausen, S. K. Kurt and N. Kockmann, Mixing and Heat Transfer in Helical Capillary Flow Reactors With Alternating Bends, 12th International Conference on Nanochannels, Microchannels and Minichannels, Chicago, Illinois, USA, 2014 Search PubMed.
  26. L. A. M. Janssen, Chem. Eng. Sci., 1976, 31, 215–218 Search PubMed.
  27. M. Johnson and R. D. Kamm, J. Fluid Mech., 1986, 172, 329–345 CrossRef CAS.
  28. J. F. Wehner and R. H. Wilhelm, Chem. Eng. Sci., 1956, 6, 89–93 Search PubMed.
  29. E. E. Hills, M. H. Abraham, A. Hersey and C. D. Bevan, Fluid Phase Equilib., 2011, 303, 45–55 CrossRef CAS.
  30. T. Durand, C. Henry, D. Bolien, D. C. Harrowven, S. Bloodworth, X. Franck and R. J. Whitby, React. Chem. Eng., 2016, 1, 82–89 Search PubMed.

Footnote

Electronic supplementary information (ESI) available. See DOI: 10.1039/c6re00109b

This journal is © The Royal Society of Chemistry 2017