Pseudo-Homopolymerization Approach to Predict ... - ACS Publications

or the discrete Galerkin approximation of the MWD by Chebyshev ...... function of the degree of polymerization (DPn), by applying the Chebyshev inequa...
1 downloads 0 Views 2MB Size
Subscriber access provided by UNIV OF DURHAM

Kinetics, Catalysis, and Reaction Engineering

Pseudo-Homopolymerization Approach to Predict the Molecular Weight Distribution in the Copolymerization via Activator Regenerated by Electron Transfer Atom Transfer Radical Polymerization Ivan Zapata-González, Enrique Saldivar-Guerra, and Jesus Ruiz-Villegas Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b02123 • Publication Date (Web): 07 Aug 2018 Downloaded from http://pubs.acs.org on August 7, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Pseudo-Homopolymerization Approach to Predict the

Molecular

Weight

Distribution

in the

Copolymerization via Activator Regenerated by Electron

Transfer

Atom

Transfer

Radical

Polymerization AUTHOR NAMES Iván Zapata-González,§ *, Enrique Saldívar-Guerra †*, Jesús Ruiz-Villegas ‡ AUTHOR ADDRESS § Catedrático CONACYT, Tecnológico Nacional de México-I.T. de Tijuana, Centro de Graduados e Investigación en Química, Tijuana, B.C., México. † Centro de Investigación en Química Aplicada, Saltillo, Coah. México. ‡ Centro de Graduados e Investigación en Química, Instituto Tecnológico de Tijuana, Tijuana, B.C., México KEYWORDS Kinetic model; Control radical copolymerization; molecular weight distribution; QuasiSteady-State Approximation; ATRP modeling

1 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 52

ABSTRACT An efficient kinetic model, based on the Pseudo-Homopolymerization (PHP) and the RSQSSA methodology, is presented to compute the molecular weight distribution (MWD) in the Activator Regenerated by Electron Transfer via Atom Transfer Radical Copolymerization. The copolymerization of butyl methacrylate and butyl acrylate is taken as basis of the study, the validation with the experimental data results in good agreement and a deviation between the prediction of the Mayo-Lewis equation and the model estimation is presented. The MWD is dependent on the fed monomer composition. A semibatch process with an extra addition of reduced agent resulted in an increased polymerization rate, but this affects the process control and a broader MWD is obtained. The computational time expended to compute the MWD is in the order of seconds or minutes, which can be reduced by using a new expression to calculate the bound of the MWD.

1. INTRODUCTION One of the most successful techniques to control the monomer addition (namely Reversible Deactivation Radical Polymerization, RDRP) in the polymer chains and to design complex molecular architectures is the Atom Transfer Radical Polymerization (ATRP),1 in which a narrow and unimodal Molecular Weight Distribution (MWD) is a primary target as a result of a well-controlled process, avoiding the dead chains obtained by termination and chain transfer events. It is well-known that the application properties of the polymer are directly correlated with the shape of MWD,2,3 and its prediction is a challenge that has been 2 ACS Paragon Plus Environment

Page 3 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

addressed in ATRP by several polymer engineers in order to produce numerical tools for analysis, control, optimization and scaling-up of the process, generating macromolecules with target properties. While polymer materials produced by homopolymerization are limited in their applications, Atom Transfer Radical copolymerization (ATRcoP) expands the range of the physical,4 chemical5 and biological properties6 of the products that can be produced, so the polymer practitioners are able to synthesize designed molecular architectures with specific copolymer composition and microstructure. However, the prediction of the detailed polymer structure, specifically in the ATRcoP, becomes more complex as the number of components is increased in the process, since the kinetic expressions to describe the copolymerization rapidly increase in complexity. A traditional treatment

for

the

modeling

of

free-radical

copolymerization

is

the

Pseudo-

Homopolymerization (PHP) method,7 by which the description of a copolymerization system of any number of monomers is treated as a homopolymerization system using average rate constants. Even though the PHP method has become a standard tool in Polymer Reaction Engineering, its application in RDRP (or controlled radical polymerization), in order to obtain averages of the MWD, has been only recently used essentially by one group,8,9,10,11,12 and only in conjunction with the method of moments (please see more on this issue in the sections where the PHP model is discussed); to the extent of our knowledge it has not been used to obtain the full MWD for RDRP copolymerization. In this work we develop a kinetic model to predict the MWD in a copolymerization via ARGET ATRP using a hybrid methodology that combines the PHP technique and the Reduced Stiffness by Quasi-Steady State Approximation (RSQSSA) method, which involves a change from an ODE system to a differential-algebraic system and its solution, using direct integration with a sequential strategy to solve the algebraic 3 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 52

states,13,14,15,16. It is important to mention that the method presented here for the calculation of the full MWD provides an easy-to-understand and easy-to-apply, intuitive methodology compared to other techniques which are available for the solution of the full MWD but are more mathematically involved, such as the probability generating functions with numerical inversion,17,18 or the discrete Galerkin approximation of the MWD by Chebyshev polynomials of a discrete variable,19 provided by the popular commercial package Predici, which requires a license fee for its use.

1.1 Activator Regenerated by Electron Transfer ATRP Before explaining the methodology which is the main subject of this work, the system that will be used for illustration is introduced and discussed in detail. The example system selected is the ARGET ATRcoP (Activator Regenerated by Electron Transfer via Atom Transfer Radical Copolymerization) of butyl methacrylate (BMA) and butyl acrylate (BA), An extensive review covering the kinetic basis of the ATRP process has been recently reported by Krys and Matyjaszewski.20 A brief description of the ATRP mechanism is shown in Figure 1, an activation/ deactivation cycle of the growing chains (Pr) controls the growth of the backbone chain, this cycle is generated involving a metal transition complex (Cu is generally used) which changes the oxidation state to lead the formation of activator (XCuI/L) and deactivator (CuIIX2/L) complexes. During this step, a halide atom is exchanged and it caps the propagating chains under more stable species (XPr), named dormant chains.

4 ACS Paragon Plus Environment

Page 5 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 1. Kinetic mechanism of normal ATRP. The parameters ka, kd, kp and kt denote the activation, deactivation, propagation and termination rate coefficients. However, a great disadvantage of the normal ATRP is the high concentration of Cu required to maintain a good development during the polymerization, since the equilibrium of the reversible reaction shown in Figure 2 is affected by bimolecular termination and the persistent radical effect (PRE).21 Therefore, the equilibrium is displaced to the formation of dormant species giving rise to a strong decrement in the polymerization rate. Techniques such as ARGET ATRP22 and Initiators for Continuous Activator Regeneration (ICAR ATRP)23 only need a low amount of Cu species (parts per million) to reach good control in the process, resulting in low dispersity values (Đ), narrow MWD, and high end-group functionality (EGF). In ARGET ATRP a reducing agent (ReAr, where r is the oxidation number) is supplied in the polymerization in order to reduce the oxidation state of the deactivator complex and to regenerate XCuI/L, resulting in a mitigation of the PRE.

5 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 52

Figure 2. Kinetic mechanism of ARGET ATRP. Butyl Acrylate (BA) and Butyl Methacrylate (BMA) polymerizations and (BA-co-BMA) copolymerizations via ARGET ATRP have been studied by Payne et al.24,25,26 using copper(II) bromide (Cu2Br), ethyl 2-bromoisobutyrate (EBIB), tin(II) 2-ethylhexanoate Sn(EH)2 and tris[(2-pyridyl)methyl]amine (TPMA), as deactivator complex, ATRP initiator, reducing agent, and ligand respectively. A low polymerization rate produced by the slow deactivation of EBIB was seen in the early stages of the reaction in the homopolymerization of BMA at 70°C,

24,25

affecting considerably the initiator efficiency

(fint=Mn,theo / Mn,exp where Mn,theo and Mn,exp denote the theoretical and the experimental number-average molecular weight, respectively). The MWD exhibited a low molecular weight tail produced by the slow initiation of the alkyl halide initiator. The BA homopolymerization at 90 °C resulted in a value of fint almost constant throughout the reaction time and a low molecular weight tail was also seen when a low Cu amount was added in the recipe.24 Copolymerization of BMA and BA at 70 °C was also carried out, a faster activation of EBIB in BA homopolymerizations, with respect to BMA, only was seen when feed compositions above 70 wt % of BA were used. However, the addition of a low 6 ACS Paragon Plus Environment

Page 7 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

amount of BA in the system gives rise to a remarkable decreasing in the Mn with higher dispersity values.26 Recently, Ribelli et al.27 reported the contribution of the bimolecular reaction termination and the catalytic radical termination (CRT) in the polymerization of methyl acrylate (MA) via

ATRP,

using

Cu2Br,

(dimethylamino)ethyl]amine

(pMA-Br)

with

(Me6TREN)

several and

ligands:

TPMA,

tris[2-

tris(3,5-dimethyl-4-methoxy-2-

pyridylmethyl) (TPMA*3). By using size exclusion chromatography (SEC), the MWD of the polymers obtained from the activation of the pMA-Br macroinitiator in more conventional radical conditions exhibited a bimodality, concluding that the higher molecular weight population corresponds to dead chains predominantly obtained by termination by combination for acrylate systems. However, a unimodal MWD was generated when controlled initial conditions were induced, predominating the CRT process for the three systems, which was confirmed by computational simulations. 1.2 Computer Simulations in ATRcoP Several simulation works in homopolymerization via ATRP have been reported in the literature in order to study complex kinetic aspects,25,26 multifunctional ATRP initiators,28,29 high pressure conditions,30 semicontinuous reactors,31,32 continuous reactors,33 reagent feed policy, diffusion controlled termination,34,35 prediction of the MWD,29,36 optimization conditions, etc. Only a few works have reported the study of the kinetic behavior in transition metal complex mediated copolymerization. Muller et al.37 simulated with the Predici ® software package the chain extension of pBA-Br using MA or MA/Styrene (10 mol%) via ICAR

7 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 52

ATRP, in which the monomer conversion was notably increased with the addition of a low amount of St, due to the fact that equilibrium ratios were modified. Payne et al.26 used a kinetic Monte Carlo (kMC) method in the study of BA and BMA in ARGET ATRcoP; the modeling results helped to understand the influence of the BA addition on the kinetics of the copolymerization, especially on the fint. The copolymer composition did not suffer a deviation from the Mayo-Lewis equation in the initial period, as seen in other controlled radical copolymerizations. An extra addition of reducing agent, at approximately 20% of conversion, increased the polymerization rate and the fint, although at the expense of lower EGF and Mn. kMC was also used by Payne et al.38 to analyze the influence of the penultimate monomer unit (PMU) on the copolymerization of BA-MMA via ICAR ATRP with feeding policy. As expected, the terminal model showed differences with the solution by the PMU model. Wang et al.39 used the method of moments to analyze the molecular properties (Mn, Mw and Đ) of methyl methacrylate and 2-(trimethylsilyl) ethyl methacrylate (MMA-co-HEMATSM) via normal ATRcoP generated in three reactors (batch, CSTR and FPR). While in the batch reactor polymers with high Mn and low dispersities can be produced, in semicontinuous reactors the operation strongly controls the copolymer composition. When the process is carried out in a continuous reactor, the fundamental advantages of this reactor configuration with respect to the discontinuous reactors are present, such as a higher product quality (homogeneity of the materials) and production efficiency. However, as the final properties depend on the residence time, the continuous operation could lead to disadvantages such as higher dispersities at long residence time or low conversion at short residence times. 8 ACS Paragon Plus Environment

Page 9 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Steenberge et al.40 developed a model based in kMC to study the microstructural characteristics, such as the linear gradient quality and control over the chain length and EGF, in ATRcoP of the systems BA/MMA, BA/Styrene, and MMA/Styrene. Apparent rate coefficients were used in all the systems to include the diffusion-controlled termination rate, using the well-known composite model of RAFT chain length dependent termination (RAFT-CLD-T). Steenberge et al. concluded that the copolymerization of BA/MMA could exhibit a good linear gradient and quality and low Đ if the reaction is carried out in a batch reactor at high monomer conversion, and an appropriate catalyst system is used. 1.3 Pseudo-Homopolymerization (PHP) Approach The prediction of the MWD generated by a copolymerization process in a rigorous development implies the distinction of three chemical characteristics: radical type i (where i = 1 or 2), n number of monomeric units type 1 added in the chain and m number of monomeric units type 2 added in the chain. The overall representation for a radical polymer chain is Pi ,n,m (where i=1 or 2; n= 1, 2, 3,… ; m= 1,2,3,…). Due to the very large dimension of the mathematical system, the integration of the previous population balance equations (PBE) could be restricted to low monomer conversion, or otherwise result in very long computational times, depending on the kinetics, CPU characteristics of the computer used, etc. Another kinetic modeling approach was reported many years ago by several groups almost simultaneously in conventional free radical copolymerizations; this method implies the treatment of the PBE for a copolymer system as that of a homopolymer41 and has received different names: pseudokinetic rate constants method,7 apparent rate constants,42 and pseudo-homopolymer approach,43 which can be used in copolymerizations described by first-order Markov chain statistics (terminal model), as well as higher order statistics 9 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 52

(such as the penultimate model).7 Some assumptions are mandatory in order to derive the PHP methodology. First, the long chain hypothesis (LCH), which proposes that the main monomer consumption is due to the propagation step, and therefore the contribution given by initiation and termination (chain transfer) steps are negligible. In addition, the QuasiSteady-State Approximation (QSSA) must hold true, since the reaction rate of the propagating chain radicals is approximately zero (generation rate ≅ consumption rate) in comparison with the overall rate of slow dynamic species. In summary, the PHP approach in conventional free radical copolymerizations consists in the evaluation of the fraction of chain radical φi, and the remaining monomer molar fraction fi or polymeric species, according to the following definitions: ∞

φi =

∑  P  i n

2

n =1 ∞

∑∑  P

n

j =1 n =1

2

with

∑φ

i

(1) j



=1

i =1

fi =

[M i ] 2

∑  M j =1

j

(2) 

where the  Pni  denotes the concentration of a growing chain radical with ultimate monomer unit i and total chain length n, and [ M i ] is the concentration of remaining monomer type i. Notice that, in the case of RDRP, it is also necessary to deal with the calculation of the probabilities of the dormant species ending with a specific type of monomer unit as the last 10 ACS Paragon Plus Environment

Page 11 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

potentially active unit since, upon re-activation, the radical type will be determined by this unit. This issue is addressed below in the PHP implementation section. Afterwards, the apparent rate constants (using the definition given in Eq. 3-5 for initiation, propagation and termination, respectively) are computed from the chemical rate constant (ki, kpij and ktij being the initiation, propagation and termination rate coefficients, respectively) and the previous fraction values. 2

ki = ∑ f i k I , i

(3)

i =1

2

2

k p = ∑∑ φi f j k pij

(4)

i =1 j =1

2

2

kt = ∑∑ φiφ j ktij

(5)

i =1 j =1

2

2

ktc = ∑∑ ktcijφiφ j

(5a)

i =1 j =1

2

2

ktd = ∑∑ ktd ijφiφ j

(5b)

i =1 j =1

The next step is to write a PBE, considering a total chain length n (without the distinction of the composition of the polymer chain); for instance, Equation (6) describes the PBE for growing chain radicals in a traditional free radical copolymerization, with initiation propagation and termination steps.

11 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

d  Pni  dt

2

2

2

j =1

j =1

j =1

Page 12 of 52



2

= ∑ kI ,i [ R *][ M i ] + ∑ k pji  Pnj−1  [ M i ] −  Pni  ∑ k pij  M j  −  Pni  ∑∑ ktij  Pmj  m =1 j =1

(6) i = 1,2; j = 1,2; n= 2,3,…

Additionally, the concentration of the growing chain radicals can be described by  Pni  = [ Pn ]φi

(7)

2

where

[ Pn ] = ∑  Pni  i =1

The right-hand side of Equation (7) is substituted in Equation (6) and after the summation over i = 1 and 2, the replacement of [Mi] by fi [M], (where [M] = [M1] +[M2]) and the substitution of the apparent rate constants, the following expression is obtained: ∞ d [ Pn ] = k I [ M ]  R*  + k p [ Pn −1 ][ M ] − k p [ Pn ][ M ] − kt [ Pn ] ∑ [ Pm ] dt m =1

(8)

n=2,3,…

Note that the previous expression is analogous to that equation obtained by the derivation of the mass balance in a homopolymerization for a growing polymer chain. Also, a large reduction of the dimension of the ordinary differential equation system is possible with the application of the PHP method. Recently, the PHP method has been used by Vivaldo et al. in the modeling of copolymerization

with

crosslinking

agents

via

RDRP,

in

Nitroxide-Mediated

Polymerization,12 and Reversible Addition Fragmentation chain Transfer (RAFT),8,10,11 and

12 ACS Paragon Plus Environment

Page 13 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

ATRP,9 but only limited to the method of moments for computation of the averages of the MWD. Unfortunately, not enough details about the use of the PHP method were exposed in the reports. 2. MODEL DEVELOPMENT 2.1 Kinetic mechanism The set of reactions and the kinetic rate constant considered for the copolymerization of BA and BMA are presented in Table 1 at 70 °C, using EBIB, CuBr2 / TPMA and Sn (EH)2 as ATRP initiator, deactivator complex and reducing agent respectively, and is based in that proposed by Payne et al.26 The mechanism includes the conventional steps of a free radical copolymerization, such as initiation (reactions 7 and 8) of polymer chains, homo and crosspropagation (reactions 9-12) and termination (reactions 15-20). Additionally, the molecular control is obtained by the three reversible reactions: i) Alkyl halide initiator (R-X) equilibrium (reaction 1 and 2), ii) BA ATRP equilibrium (reaction 3 and 4), and iii) BMA ATRP equilibrium (reaction 5 and 6). The process to regenerate the activator complex through the reducing agent in this system (reactions 13 and 14) is in agreement with the two reduction steps suggested by Jakubowski et al.44 In this work, the ATRP copolymerization model is carried out under the assumptions of terminal model kinetics, isothermal process, perfect mixing in a batch reactor, and intrinsic, non-chain-length dependent, values of the kinetic rate constants. Intramolecular and intermolecular chain transfer to polymer mechanisms are important phenomena in the acrylate polymerizations,45,46 at least at temperatures at or above 80 °C, and mid-chain radicals formed by these mechanisms could present lower polymerization rate; for simplicity, for a first application of the methodology,

13 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 52

a rough representation of this mechanism is taken here into account via the decrement of the propagation rate coefficient, as proposed by Fu et al.31 (see note a in Table 1). 2.2 Population Balance Equations The mass balance equations are derived from the kinetic mechanism, Table 1. The mathematical system is a set of stiff ODEs with initial conditions. PBEs corresponding with small species and polymer species are presented in the following Equations (9-24), where in the case of polymeric chain species the subscript denotes the chain length and the superscript the ultimate monomer unit. Primary radicals: d  R* 

dt

= − k p1  R*  [ M 1 ] − k p 2  R*  [ M 2 ] + ka 0 [ R − X ] Cu I X / L  − kd 0  R*  Cu II X 2 / L  (9)

Monomers (M1 = BMA and M2 = BA) ∞ ∞ d [ M1 ] = −k p1  R*  [ M 1 ] − k p11 [ M 1 ] ∑  Pn1  − k p 21 [ M 1 ] ∑  Pn2  dt n =1 n =1

(10)

∞ ∞ d [M 2 ] = −k p 2  R*  [ M 2 ] − k p12 [ M 2 ] ∑  Pn1  − k p 22 [ M 2 ] ∑  Pn2  dt n =1 n =1

(11)

Activator and deactivator complexes d Cu II X 2 / L  dt





= ka 0 Cu I X / L  [ R − X ] + k a1 Cu I X / L  ∑  Pn1 − X  + k a 2 Cu I X / L  ∑  Pn2 − X  n =1

n =1





n =1

n =1

−kd 0 Cu II X 2 / L   R*  − kd 1 Cu II X 2 / L  ∑  Pn1  − kd 2 Cu II X 2 / L  ∑  Pn2  − kr1 Cu II X 2 / L   Sn II L2  −kr 2 Cu X 2 / L   Sn L2 X  II

III

(12)

14 ACS Paragon Plus Environment

Page 15 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Table 1. Set of reactions and kinetic parameters for the copolymerization of BMA=1, BA = 2 in an ARGET ATRcoP of BA and BMA at 70°C. Reaction number

Reaction

Activation / Deactivation 1 R - X + Cu I X / L → R * + Cu II X 2 / L 2 R * + Cu II X 2 / L → R - X + Cu I X / L 3 Pn1 - X + Cu I X / L → Pn1 + Cu II X 2 / L 4 Pn1 + Cu II X 2 / L → Pn1 - X + Cu I X / L 5 Pn2 - X + Cu I X / L → Pn2 + Cu II X 2 / L 6 Pn2 + Cu II X 2 / L → Pn2 - X + Cu I X / L Initiation * 7 R + M 1 → P11 8 R * + M 2 → P12 Propagation 1 9 Pn + M 1 → Pn1+1 10 Pn1 + M 2 → Pn2+1 11 Pn2 + M 1 → Pn1+1 12 Pn2 + M 2 → Pn2+1 Reduction II II 13 Cu X 2 / L + Sn L2 → Cu I X / L + Sn III L2 X 14 Cu II X 2 / L + Sn III L2 X → Cu I X / L + Sn IV L2 X 2

15 16 17 18 19 20 a. b.

Kinetic Value Coef. (L mol-1 s-1) ka0 kd0 ka1 kd1 ka2 kd2

3.4

25,26

2.4x107

26, 47

2.4x102

26

1.2x107

25, 26

8

b

2x107

b

kp1 kp2

1.23x103

26, 48

4.32x104

49

kp11 kp12 kp21 kp22a

1.23x103

48, 26

6.54x102

26

2.06x104

b

6.71x103

b

3x10-1

25, 26

1

25, 26

kr1 kr2

Termination ktc11 P + Pm1 → Dn + m 2 ktc12 P + Pm → Dn + m 2 ktc22 P + Pm → D n + m ktd11 Pn1 + Pm1 → D n + D m 1 2 ktd12 Pn + Pm → D n + D m ktd22 Pn2 + Pm2 → D n + D m 4 -1 -1 26 Original values used by Payne et al. kp22 = 4.32x10 L mol s , using a factor of 1/7 times, proposed by Fu et al.31 This work 1 n 1 n 2 n

Ref.

9.0x106

36

3.23x108

c

1.16x1010

50

8.1x107

36

3.23x108

c

1.29x109

50

which is modified

15 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 52

ktc12 = ktc11ktc 22 and ktd 12 = ktd 11ktd 22

c.

d CuI X / L dt





= −ka0 CuI X / L [ R − X ] − ka1 CuI X / L ∑ Pn1 − X  −ka2 CuI X / L ∑ Pn2 − X  n=1

n=1





+ kd 0 CuII X 2 / L R*  + kd1 Cu II X 2 / L ∑ Pn1  + kd 2 CuII X 2 / L ∑  Pn2  n=1

n=1

+kr1 Cu X 2 / L Sn L2  + kr 2 Cu X 2 / L Sn L2 X  II

II

II

III

(13) Reducing agent in three oxidation states d  Sn II L2 

= −kr1  Sn II L2  Cu II X 2 / L 

dt

d  Sn III L2 X 

dt d  Sn IV L2 X 

dt

(14)

= kr1  Sn II L2  Cu II X 2 / L  − kr 2  Sn III L2 X  Cu II X 2 / L 

(15)

= kr 2  Sn III L2 X  Cu II X 2 / L 

(16)

Alkyl halide species d [R − X ]

dt

= − ka 0 [ R − X ] Cu I X / L  + kd 0  R*  Cu II X 2 / L 

(17)

Dormant polymer species d  Pn1 − X 

dt d  Pn2 − X 

dt

= − ka1  Pn1 − X  Cu I X / L  + kd 1  Pn1  Cu II X 2 / L 

(18)

= − ka 2  Pn2 − X  Cu I X / L  + kd 2  Pn2  Cu II X 2 / L 

(19)

Growing polymer chains with n= 1

16 ACS Paragon Plus Environment

Page 17 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

d  P11  dt

= k p1  R*  [ M 1 ] − k p11 [ M 1 ]  P11  − k p12 [ M 2 ]  P11  +ka1 Cu I X / L   P11 − X  ∞





n =1

n =1

n =1

− kd 1 Cu II X 2 / L   P11  − ktc11  P11  ∑  Pn1  − ktc12  P11  ∑  Pn2  − ktd 11  P11  ∑  Pn1  ∞

− ktd 12  P11  ∑  Pn2  n =1

(20) d  P12  dt

= k p 2  R*  [ M 2 ] − k p 21 [ M 1 ]  P12  − k p 22 [ M 2 ]  P12  +ka 2 Cu I X / L   P12 − X  ∞





n =1

n =1

n =1

− kd 2 Cu II X 2 / L   P12  − ktc 21  P12  ∑  Pn1  − ktc 22  P12  ∑  Pn2  − ktd 12  P12  ∑  Pn1  ∞

− ktd 22  P12  ∑  Pn2  n =1

(21) Growing chain species with chain length n˃1 d  Pn1  = k p 11  Pn1−1  [ M 1 ] + k p 21  Pn2−1  [ M 1 ] − k p 11  Pn1  [ M 1 ] − k p 12  Pn1  [ M 2 ] dt ∞

+ k a 1  Cu I X / L   Pn1 − X  − k d 1  Cu II X 2 / L   Pn1  − k tc11  Pn1  ∑  Pm1  m =1







m =1

n =1

n =1

− k tc12  Pn1  ∑  Pm2  − k td 11  Pn1  ∑  Pm1  − k td 12  Pn1  ∑  Pm2  (22)

d  Pn2  = k p 22  Pn2−1  [ M 2 ] + k p 12  Pn1−1  [ M 2 ] − k p 21  Pn2  [ M 1 ] − k p 22  Pn2  [ M 2 ] dt ∞

+ k a 2  C u I X / L   Pn2 − X  − k d 2  C u II X 2 / L   Pn2  − k tc 12  Pn2  ∑  Pm1  m =1







m =1

n =1

n =1

− k tc 22  Pn2  ∑  Pm2  − k td 21  Pn2  ∑  Pm1  − k td 22  Pn2  ∑  Pm2  (23)

17 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 52

Dead polymer chains n˃1 d [D n ] dt

= +

n −1 n −1 1 1 k tc 1 1 ∑  Ps1   Pn1− s  + k tc 1 2 ∑  Ps1   Pn2− s  + 2 2 s =1 s =1 n −1 ∞ ∞ 1 k tc 2 2 ∑  Ps 2   Pn2− s  + k td 1 1  Pn1  ∑  Ps1  + k td 1 2  Pn1  ∑  Ps2  2 s =1 s =1 s =1



+ k td 2 2  Pn2  ∑  Ps 2  s =1

(24)

2.3 Method of Moments The definitions of the moments k (where k = 0, 1 and 2) with radical type i (where i = 1 or 2), for the growing chain species (α), dormant species (γ) and dead chain species (β), are shown in the Equation (25-27). ∞

α ki = ∑ n k  Pni 

(25)

s =1



γ ki = ∑ n k  Pni − X 

(26)

s =1



β k = ∑ n k [ Dn ]

(27)

n =1

As in previous studies, zeroth, first and second moments are obtained for each polymeric species, such as growing, dormant or dead polymer species, and also for each ultimate

18 ACS Paragon Plus Environment

Page 19 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

monomeric unit. The resulting expressions are shown in the Appendix A (Supporting Information), Equation (A1-A15). Number- and weight-average molecular weight, Mn and Mw, and their ratio in the dispersity index (Đ), are estimated by the well-known expression, given in the Equations (28-30). These values are directly compared with the experimental molar mass obtained by SEC measurements of absolute or relative nature with the correction of the MW using the Mark– Houwink–Sakurada constants. End-group functionality is calculated by Equation (31).

Mn =

2

2

2

k =1 2

k =1 2

k =1 2

∑ α1k + ∑ β1k + ∑ γ 1k ∑α + ∑ β + ∑ γ k 0

k 0

k =1

k =1

2

2

∑α + ∑ β + ∑ γ

k 2

∑α + ∑ β + ∑ γ

k 1

k 2

k =1 2

k =1 2

k 1

Đ=

k =1 2

k 1

k =1

k =1

( MW1Facc1 + MW2 Facc 2 )

(29)

k =1

Mw Mn

(30) 2

EGF =

(28)

k =1

2

k 2

Mw =

( MW1 Facc1 + MW2 Facc 2 )

k 0

2

∑ α 0k + ∑ β0k k =1

2

k =1

2

k =1

k 0

k =1

(31)

2

∑α + ∑ β + ∑ γ k 0

k 0

k =1

Where MW1 and MW2 are the molecular weight of BA and BMA, respectively; and Facc1 and Facc2 represent the cumulative copolymer composition of BA and BMA, respectively, which are estimated with Equation (32a and 32b):

19 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

M 10 − M 1 ( M 10 + M 20 ) − ( M 1 + M 2 )

Facc1 =

Page 20 of 52

(32a)

Facc 2 = 1 − Facc1

(32b)

M10 and M20 are the initial concentrations of BA and BMA, respectively.

2.4 Pseudo-Homopolymerization Model As previously described, the first step of the Pseudo-Homopolymer approach is to establish the fraction of radicals type i, φi, and the fraction of remaining monomer type i, fi, as shown in the Equation (1 and 2). Additionally, in ARGET ATRP an extra definition of the fraction for dormant species with terminal unit i is necessary, as:

vi =

 Pi − X    2

∑ j =1

P j − X   

2

with

∑v

i

 Pi − X   =  P − X 

(33)

= 1 , [P-X] being the total concentration of dormant species.

i =1

Notice that in previous implementations of the PHP in other RDPR chemistries, done by Vivaldo et al for the calculation of moments of the MWD (not the full MWD)8,9,10,11,12 the authors unfortunately do not describe how they calculate these important probabilities, so it is not possible to know what underlying assumptions they considered. This is important and not a trivial issue, since these calculations can be done in different ways and, depending on the underlying assumptions, it is anticipated that the results can be different. Given the importance of this issue, it will be pursued in a future publication of our group but, for the time being, the most rigorous way of calculating these probabilities has been adopted in this

20 ACS Paragon Plus Environment

Page 21 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

work. This emphasizes one of the main contributions of this work: the detailed description for the application of the PHP in RDRP copolymerizations for the calculation of the full MWD. Following the methodology described and reported previously, and after some simple algebraic manipulation, the resulting equation system to solve the MWD for the ATRP copolymerization is reduced to a PHP model, Equation (34-43). First radicals

d  R* 

dt

= -k I  R*  [ M ] + ka0 [ R - X ] Cu I X / L  - kd0  R*  Cu II X 2 / L  (34)

Monomer ∞ d [M ] = −k I  R*  [ M ] − k p [ M ] ∑ [ Pn ] dt n =1

(35) Activator and Deactivator complexes ∞ 2 d Cu I X / L  = − k a 0 Cu I X / L  [ R − X ] − ka Cu II X / L  ∑ [ Pn − X ]+ kd Cu II X 2 / L  ∑ [ Pn ] dt n =1 n =1

+ k d 0 Cu II X 2 / L   R *  + k r1 Cu II X 2 / L   Sn II L2  + kr 2 Cu II X 2 / L   Sn III L2 X 

(36) d Cu II X 2 / L  dt



2

= k a 0 Cu I X / L  [ R − X ] + ka Cu II X / L  ∑ [ Pn − X ] − k d Cu II X 2 / L  ∑ [ Pn ] n =1

n =1

− k d 0 Cu X 2 / L   R  − k r1 Cu X 2 / L   Sn L2  − k r 2 Cu X 2 / L   Sn L2 X  II

*

II

II

II

III

21 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 52

(37) Growing polymer species with chain length n=1

d [ P1 ]

= kl  R*  [ M ] − k p [ M ][ P1 ] +ka Cu I X / L  [ P1 − X ] − kd Cu II X 2 / L  [ P1 ]

dt





m =1

m =1

− [ P1 ] ∑ ktc [ Pm ] − [ P1 ] ∑ ktd [ Pm ]

(38) Growing polymer species with chain length n > 1

d [ Pn ]

= k p [ M ][ Pn−1 ] − k p [ M ][ Pn ] + ka Cu I X / L  [ Pn − X ] − kd Cu II X 2 / L  [ Pn ]

dt





m =1

m =1

− [ Pn ] ∑ ktc [ Pm ] − [ Pn ] ∑ ktd [ Pm ]

(39) Dormant polymer species

d [ Pn − X ]

dt

= − ka Cu I X / L  [ Pn − X ] + kd Cu II X 2 / L  [ Pn ]

(40)

Dead polymer species

d [ Dn ] dt

=

∞ 1 n k P P P ktd [ Ps ] + [ ][ ] [ ] ∑ tc s n−s n ∑ 2 s =1 s =1

(41)

Note that the reducing agent with two oxidation states and the alkyl halide do not suffer any changes and their mass balances remain identical to Equation (14-16) and (17), respectively. The definition of the initiation, propagation, termination and ATRP

22 ACS Paragon Plus Environment

Page 23 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

deactivation/ activation apparent rate constants are described in the Equation (3-5) and (42 and 43), respectively. The final ODE system includes the Equation (14-17, 34-41). 2

ka = ∑ kai vi

(42)

i =1

2

kd = ∑ kd iφi

(43)

i =1

2.5 RSQSSA Method The reduction of the PBE produces a system analog to the ARGET ATRP homopolymerization model reported by our group.36 The total concentration of the three polymer populations is needed during the integration with a sequential calculation for the live polymer population and a stable explicit numeric method. Throughout the text the superscript Tot means the addition of both polymer chain types (Tot =1+2) giving rise to total concentration of the polymer species. For instance, the total concentration for growing species is derived from the PHP model, using the following Equation (44): ∞

α kTot = ∑ s k [ Ps ] = α k1 + α k2

(44)

s =1

Zeroth moments for dormant species type 1 ( γ 01 ) and 2 ( γ 12 ) are preserved in the final model, with the purpose of the evaluation of

vi , thus Equation (A7 and A10) are contained

within the model, and the zero moments for the type i growing radicals are computed with the expression α 0i = φiα 0Tot . As previously seen,36 the total concentration of growing radicals 23 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 52

could be described under QSSA, since the processes leading to its formation are essentially the same as those consuming it. Setting to zero the derivative of the total concentration of

 dα kTot  = 0  and after algebraic procedures, a quadratic equation is growing radicals   dt  obtained, Equation (45).

α

Tot 0

=

−b ± ( b 2 − 4ac )

(45)

2a

Where:

a = −(ktc + ktd ) , b = − k d Cu II X 2 / L 

c = ki [ M ]  R*  + ka Cu I X / L   PTot − X 

In order to take advantage of the dynamic differences exhibited by the live (fast dynamics) polymer species and the rest of the species (slow dynamics), the RSQSSA method is applied, considering as unstable species (fast) the first radicals and the growing polymer species. Therefore, following the RSQSSA methodology, the overall rates for R* and Pr are negligible with respect to those exhibited by slow species, such as M1, Pn-X or Dn, which

 d [ R *] d [ Pn ]  = = 0  . The resulting Equation (46) allows to set to zero their derivatives  dt  dt  (48) are derived by setting to zero the derivatives of the differential equations (34), (38) and (39), yielding respectively: 24 ACS Paragon Plus Environment

Page 25 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

k a 0 [ R − X ] Cu I X / L   R  = ki [ M ] + k d 0 Cu II X 2 / L  *

(46)

ki [ M ][ R *] + k a Cu I X / L  [ Pr − X ] [ P1 ] = k p [ M ] + k d Cu II X 2 / L  + ktd + ktc α 0Tot

(47)

k p [ M ][ Pn −1 ] + ka Cu I X / L  [ Pr − X ] k p [ M ] + kd Cu II X 2 / L  + ktd + ktc α 0Tot

(48)

(

[ Pn ] =

(

)

)

As shown in Table 1, the chain length dependence (CLD) of the termination rate coefficient has been neglected in this model, using only intrinsic rate coefficients. However, an extension of the RSQSSA model to compute the apparent nature of the CLD termination rate coefficient has been reported by De Rybel et al.51 That procedure avoids the computation of the moments by the direct solution of the ODE system, since the moments method is only valid for constant values of the rate coefficients to obtain a “closed” system. Additionally, number average chain length-dependent termination rate coefficients, such as the expression proposed by Griffiths et al.52 for BMA and MMA, can be easily incorporated in the model without any strong impact in computation times.

2.6 Numerical Aspects To assess the validity of the QSSA using only moment equations, the ODE system of the exact moments model, Equation (9-17 and A1-A15) was programmed in Fortran language, and the solution was integrated using DDASSL, which is considered as a robust algorithm available to solve differential-algebraic equation systems or stiff differential equation systems. It employs an implicit technique with adaptive step size and a collection of integrators of the same type as those employed in Gear's method.

25 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 52

On the other hand, a second order explicit technique, specifically Heun, was used in the solution of the algebraic-differential system to predict the MWD in the PHP model. The algebraic part (growing polymer population) was solved at each time step using sequential explicit expressions (Equation 45-48). The model comprises the Equation (14-17, 35-37, 40, 41, and 45-48). A previously reported algorithm was adapted to include the equations of the PHP model.36 Their moments (k) were obtained by direct summation using the moment definitions shown in Equation (25-27). The flow diagram of the algorithm used to describe the ARGET ATRcoP kinetics, programmed also in Fortran code, is shown in Figure 3, where N must be the number of iterations to reach a final condition (tfin or conversion) with the specified time step. A maximum chain length (CLmax) is an input to the program as initial guess; however, the CLmax could be either kept constant or it can be modified as a function of the degree of polymerization (DPn), by applying the Chebyshev inequality.53

CLmax ≤ DPn + K × ( DPmed − DPn 2 )

1/ 2

(49)

Where DPmed = µ2 / µ1, µ2 and µ1 is the second and first total moment, respectively; K is a safety factor (times the standard deviation) in order to capture most of the information. Also, in this work the following empirical approximation is used. CLmax = DPn + 30 ( DPw − DPn )

12

(50)

As initial condition, CLmax = 100 in order to ascertain that all chain lengths are taken into account in the computation of the MWD. For the simulations we used a conventional computer with 8 GB RAM memory and Intel Core i7 7th Gen processor running at 2.7 GHz.

26 ACS Paragon Plus Environment

Page 27 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 3. Flowchart of the algorithm used in the numerical solution of the PHP model for ARGET ATRcoP.

3. RESULTS AND DISCUSSION 3.1 Validation of the QSSA and the Equilibrium The first three moments of the PHP/QSSA model obtained by the integration of the differential part and involving the QSSA hypothesis, are compared with those generated by the direct integration of the moments model without the QSSA assumption, Figure 4. The direct integration of the moments model without the QSSA (exact model) requires the use of an integration algorithm for stiff systems (DDASSL in our case) as a result of the presence of the time derivatives for the moments of the live polymer, as explained before. 27 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 52

The validity of the QSSA assumption on R* and moments of the Pn is confirmed, since for the three moments of each polymer species both solutions are overlapped. Hence, the fast dynamic behavior exhibited by short-lived species (R*, Pn1 and Pn2) can be well-described by algebraic expressions and the corresponding time derivatives vanish from the ODE system. Practically, the QSSA hypothesis is valid from the beginning of the process and throughout the full course of the polymerization, the transient period to reach the QSSA condition being negligible (in the order of 10-2 s). The concentrations of the small species (such as R*, R-X, CuIIX2/L, CuIX/L) are also exactly computed by the PHP/QSSA model in comparison with the exact moment solution. Another important aspect is the establishment of the chemical equilibria in the activation/deactivation processes, Reaction (1-6) in Table 1: i) equilibrium between small species involving R-X (Reaction 1 and 2), ii) equilibrium formed between the growing radicals with terminal monomer unit of BMA (α01) and the deactivator complex (Reaction 3 and 4), and iii) equilibrium formed between the growing radicals with terminal monomer unit of BA (α02) and the deactivator complex. Chemical equilibrium is defined as the

ultimate point in where the rate in both directions are identical for a reversible process, so that the system gives the appearance of having a static composition at which the Gibbs energy is minimum, and the equilibrium constant, K, is given by the mass-law effect.54 So, on a mathematical basis, the derivatives in time of the involved components in an equilibrium process will be exactly zero (RA/D,i = 0, i= 0, 1, 2).

28 ACS Paragon Plus Environment

Page 29 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 4. Comparison of the first three moments obtained by the model of moments (lines) and the PHP/QSSA model (marks) in an ARGET ATRcoP. a) Evolution of the zeroth moment, b) evolution of the first and c) second moment. The moments are denoted as αkTot, βkTot and γkTot, Initial molar ratio [M]:[R-X]:[CuIIX2/L]:[SnIIL2]=39:1:0.00125:0.0125 at 70 ªC, with 70 wt % BA in the feed monomer composition. Kinetic parameters are shown in Table 1. The rates RA/D,0, RA/D,BMA and RA/D,BA for both models (exact and PHP/QSSA) corresponding with Equation (17, A7 and A10), respectively, are shown in Figure 5a. A strong deviation from the equilibrium conditions is seen in the first stages of the copolymerization, approximately until 1000 min. During this pre-equilibrium period, the deactivation rate is much higher than the activation of R-X, since the production of activation complex from the two reduction steps is very low. However, higher initial concentration levels of the SnIIL2 29 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 52

lead to a negative effect in the conservation of the EGF, affecting the copolymer composition in the chains. Both rates of polymer species RA/D,BMA and RA/D, BA are positive due to the high concentration of deactivator complex and the very high deactivation rate constant. BMA reaches its equilibrium condition much before BA, since the former exhibits a monomer feed composition of only fBMA = 0.3. As previously reported,55,56 the knowledge of the characteristics of the pre-equilibrium stage in the process is deemed as an important information since it is a condition to assess the validity of the Mayo-Lewis (ML) equation in a specific system. Additionally, as predicted by Payne et al., the instantaneous copolymer composition computed by their model does not present strong deviations with respect to the ML equation in the BA/BMA copolymerization at high conversion, due to the fact that activation/deactivation rate coefficients are not very different (Table 1). Note that in that work, the instantaneous composition is only obtained for a specific conversion value, but the first stage of the copolymerization analysis was not analyzed. However, an imbalance in the rates captured by the PHP/QSSA model, can be observed below 60% of conversion, in both instantaneous and cumulative copolymer composition, see Figure 5b and 5c respectively. The discrepancies of the PHP/QSSA model with respect to the ML equation prediction can be attributable to the differences in the deactivation/activation rates and the reactivity ratios (Ri

= kii/ kij), as explained before,55,

56

also known as competitive equilibria and processes.57

Additionally, Fins,1 and Fcum,1 obtained by the exact moments model overlap with the PHP/QSSA model simulations, confirming these results.

30 ACS Paragon Plus Environment

Page 31 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 5. a) Comparison of the individual activation/ deactivation rates on time for each dormant type, model of moments (lines) and PHP model (marks). Labels: RA/D, 0: terms of LHS Equation 17; RA/D, BA: terms of LHS Equation A7; RA/D, BMA: terms of LHS Equation A10. b) Evolution of Fins1 on overall conversion. c) Evolution of Fcum1 on overall conversion. Same recipe and conditions as those used in Figure 4; the kinetic parameters are shown in Table 1.

3.2 Validation of the PHP Model and MWD Prediction The solution of the PHP model (lines) is compared with experimental data reported by Payne et al.26 (marks), showing a very good fitting in the evolution of conversion, Mn, dispersity and initiator efficiency, Figure 6a, b, c and d respectively. Most of the kinetic parameter values were taken from literature sources (see Table 1); the only parameters 31 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 52

fitted were those of activation and deactivation of polymeric units ended with monomer 2 (BA), as well as the homopropagation kinetic constant of BA, which turned out to be an important parameter to fit the experimental conversion. As previously seen by Schroeder et al.31 in the polymerization of BA via normal ATRP at 90°C using N,N,N’,N’’,N’’pentamethyldiethylenetriamine (PMDETA) as ligand, the polymerization rate can be perturbed by side reactions not taken into account in the kinetic mechanism shown in Table 1. According to that report, the slow polymerization rate exhibited could be attributed to either an intramolecular chain transfer producing a midchain radical with lower propagation rates, or a chain transfer between the ligand and the acrylate radical, as recently reported.27 Figure 6a also shows the evolution of each monomer in time, M1 (BMA) being consumed faster than M2 (BA) due to the asymmetry of initial concentrations and the kinetic parameters. On the other hand, the strong deviation between Mn,exp (also Mn obtained by PHP model) and the Mn,theo is due to the poor efficiency of the ATRP initiator producing, from the beginning of the process, a low amount of polymer chains with longer chain length than desired, as previously reported.25,36,58 Although this process exhibits characteristics of a well-controlled copolymerization, the dispersity index is close to 1.6 throughout the reaction.

32 ACS Paragon Plus Environment

Page 33 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 6. Comparison of experimental (marks) and PHP/QSSA solution data (lines) of the ARGET ATRP copolymerization of BMA and BA with 70 wt % BA in the feed monomer composition. a) Evolution of conversion in time, b) Mn, c) Dispersity and d) Initiator efficiency versus time. Same recipe and conditions as those used in Figure 4; the kinetic parameters are shown in Table 1. In the previous section the results shown of the PHP/QSSA model are related to average quantities, but its main feature is that it provides the full MWD. The molecular weight distribution obtained at 10, 30, 50 and 70% of overall conversion is shown in Figure 7a, a large tail of low molecular weight is produced by the slow initiation process of the alkyl halide. The generation of dead polymer chains is avoided, its concentration being around two orders of magnitude lower than that of the dormant species and therefore the contribution of the dead chains in the full-MWD being negligible. The PHP/QSSA model is

33 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 52

able to capture the evolution of the partial MWD corresponding to each species population (growing, dormant and dead chains) in the conversion, as shown in Figure 7b. Dormant chains are the preponderant species, their concentration being around two orders of magnitude higher than for the dead polymer, which is increased throughout the copolymerization. Maximum concentrations for both dormant and dead chains are obtained in the final conversion (at 70%) of the simulation, but for the growing chains a maximum is found at 27% of conversion, afterward, the concentration falls attributable to more termination events. The evolution of the MWD with conversion is illustrated in Figure 7c, in which the effect of the slow initiation period is clearly seen in the first stage of the process. Additionally, propagating radicals are kept at very low levels of concentration, as seen in Figure 7b; however, it is important to maintain a suitable concentration of them in order to keep a good polymerization rate. In that way, it has been seen that increasing the content of BA in the feed monomer composition, the copolymerization rate increases due to the higher propagation rate coefficient and the faster activation of the ATRP initiator. As a result of the increment of BA content, the ATRP initiator is activated faster in early stages of the copolymerization and a high number of chains is produced, which are perfectly in dormant states, leading to a lower Mn than that obtained by only BMA, Figure 7d. In fact, the results are qualitatively in accordance with experimental data obtained by Payne et al.26 and theoretical studies about the initiation process exhibited in ARGET ATRP.36

34 ACS Paragon Plus Environment

Page 35 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 7. Molecular weight distribution predicted for the ARGET ATRP copolymerization of BMA and BA. a) Evolution of the distribution at 10, 30, 50 and 70 % (labels) of overall conversion. b) Partial MWD for each type of polymer species (Pn-X = dormant, Dn = dead, and Pn = propagating radical). c) MWD for all polymer species with 70 wt % BA in the feed monomer composition d) Effect of variation of the initial monomer composition BA: BMA on the MWD at 70% of conversion.

35 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 52

As it is known, the experimental measurement of the SEC MWD for copolymers is a complex issue, compared with the homopolymerization case. However, good estimations can be obtained by using the concept of the universal calibration curve and adequate MarkHouwink coefficients. Hence, the validation of the kinetic mechanism via the direct comparison between the theoretical MWD and the SEC measurement is in principle possible. Furthermore, semi-batch reactor strategies have been studied in which an additional solution of reducing agent is fed in the system; if a low concentration solution is added in the reactor, the effect over the molecular structure properties (Mn and Đ) is almost negligible. In order to delve deeper into this interesting behavior, we are expanding the investigation following what was initially explored by Payne et al.26 In the following simulations, a high concentration of reducing agent (0.1 M) is added after 10, 20, 30, 40, 50 and 60 % of overall conversion, assuming that the addition does not change the density of the reaction mixture. As expected, a strong increment in the polymerization rate is obtained after the supplementary charge of reducing agent occurs, Figure 8a; since a higher concentration of CuIX/L complex is regenerated and therefore the equilibrium ratio is modified to produce more growing chains. Note that the slope of the monomer consumption in a semi-batch process is the same, regardless of the injection time for the extra charge. Moreover, the activation/ deactivation cycles become faster, which leads to an increase in the number of propagation events, resulting in polymers with longer chains and a broader MWD at 70% of conversion, Figure 8d, and higher Mn, Figure 8b. A minimum in the dispersity is predicted after the addition of reducing agent with respect to the batch system, Figure 8c, which is due to a sudden increase in the degree of polymerization, but

36 ACS Paragon Plus Environment

Page 37 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

the Mw values gradually turn higher, even surpassing the dispersity value obtained by a batch system at 70% of conversion. The addition of extra reducing agent has a significant effect on the EGF, Figure 8e; as a consequence of the higher concentration of growing chains, there are probably more bimolecular termination events. When a supplement of reducing agent is added in the first stages of the reaction, the equilibrium ratio changes to build up dead polymer chains, producing a severe decrement in the EGF. Otherwise, if the injection of a reducing agent is done during the last stages, it also results in a high quantity of dead polymer, but the reaction is stopped before a strong reduction of the EGF occurs. Generally, the initiator efficiency is only slightly decreased with the incorporation of extra reducing agent, Figure 8f; for instance, the value of fin in a batch system is 0.83, decreasing at 0.73 when an extra addition of reducing agent is injected at 10% of conversion. A solution of lower concentration results in an insignificant effect on the values of fin, as seen in the experimental data reported previously.26

37 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 52

Figure 8. Effect of the variation on the addition of reducing agent at different conversions (10, 20, 30, 40, 50 and 60%) in a semi-batch system in ARGET ATRcoP, with 70 wt % BA in the feed monomer composition. Evolution of a) conversion on time, b) Mn and c) dispersity on conversion, d) chain length distribution, e) EGF and f) Initiator efficiency on conversion. Same recipe as in Figure 4 and the kinetic parameters are shown in Table 1.

38 ACS Paragon Plus Environment

Page 39 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

3.3 Computation Times and Numerical Accuracy A very short computation time was required in order to describe the MWD in the ARGET ATRP process by the PHP model, only 14.22 s took to compute a simulation until 70% of overall conversion with the base case, using a time step (∆t) of 0.1 s and maximum chain length of 300, see Figure 9a. As ∆t is increased to 0.4 s the computation time is reduced to a short time of 3.38 s. The integration efficiency is improved due to the diminution of the dimension in the ODE system, from 1800 equations of the original system to 600 equations of the PHP model (considering a maximum chain length of 300), but more important is the fact of the stiffness reduction. Additionally, considering CLmax as the function of average properties of the total polymer product, such as DPn and DPw, Equation (50), the running time is decreased around a half, even the time expended is lower than that obtained using the Chebyshev inequality with K=10 proposed by Wulkow.53 The safety factor could be reduced to K = 6, lower values do not capture the full MWD information, as explained below. Additionally, Figure 9b shows a contour plot of the MWD obtained by the PHP model and the bounds calculated by the four expressions of CLmax. When a constant value of CLmax = 300 is used, half of the calculations are significant for the description of the MWD and the rest of the information is negligible. The bound obtained by the Chebyshev inequality, using K = 10, reaches too far from the more dense region of the population, but decreasing the safety factor to K= 6, the calculated bound is closer. However, caution is recommended to use these values since with an initial CLmax lower than 50, it produces an incomplete distribution at low conversion (˂ 20%), obtaining incorrect predictions at high conversions. The empirical expression, Equation (50), leads to a much more efficient adaptive bound, which follows the shape of the level contours and runs parallel to them in 39 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 40 of 52

the MWD, only those chain lengths with significant values of polymer concentrations are computed. More studies about this empirical expression are being done. Figure 9c shows the relative error (ε) of the total concentration of the growing, dormant and dead chains, α0Tot, β0Tot and γ0Tot respectively, versus ∆t used during the simulation. The relative error is calculated by the following equation:

ε=

g0 − g g

(51)

Where g denotes any total moment obtained by exact moment equations (α0Tot, β0Tot, and

γ0Tot) and g0 is the corresponding total moment obtained by summation, Equation (25-27). The errors computed for all species remain at very low values, the error for dead chains being one order of magnitude higher than those of growing and dormant species. Even more, the relative errors of the dead chains present a minimum value around 0.2 s, but a good linear fitting in the Log-Log plot is presented for both growing and dormant chains. Figure 9d illustrates the MWD curves obtained by the variation of ∆t; values in the range 0.01 ≤ ∆t ≤ 0.4 s produce practically the same MWD shape. However, when ∆t ˃ 0.4 s a slight deviation of the curves is generated. Real-time optimization and control applications demand kinetic algorithms with reliability, high efficiency and short computation times. The model developed here complies satisfactorily with those characteristics, using special strategies (PHP and RSQSSA methodology) in order to keep the computational effort low.

40 ACS Paragon Plus Environment

Page 41 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

41 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 42 of 52

Figure 9. Plots of computation time and accuracy of the PHP model. a) Relation between computation time and the time step ∆t in the PHP model; b) Contour plot of the MWD using CLmax=300 and the bound obtained by: CLmax = Eq. (49) denoting the Chebyshev inequality, with K= 6 and 10; and CLmax denoting the empirical expression Eq. (50). c) Error between the method of moments and the integration of the PBE, for each polymeric species (growing, dormant and dead chains, sum of type 1 and 2) using a Heun method for the time integration. d) MWD produced by the variation of time step in the computational time, using constant chain length (CLmax) obtained by the PHP model, at 0.01 0.1, 0.4. a 0.5 s. e) Relation between simulation time and the maximum chain length (CLmax) obtained by the PHP model, using a ∆t = 0.4 s. Same recipe as in Figure 4 and the kinetic parameters are shown in Table 1.

4. CONCLUSIONS A new methodology, using both the pseudo-homopolymerization and RSQSSA techniques, has been developed to compute the molecular weight distribution in the copolymerization via RDRP. It has been demonstrated the viability of the new methodology to reliably compute, from a complex kinetic mechanism, the average polymerization characteristics such as partial and overall monomer conversion, Mn, Mw, Đ, EGF and fin. The methodology could also be applicable to other copolymerization systems with minimum or no changes, and also adapted to other controlled techniques. More importantly, the complete MWD and the contribution of each chain species were well-described as a result of the efficient algorithm. A rigorous validation of the QSSA has demonstrated that the fast dynamics exhibited by the growing chains, can be described by algebraic equations, reducing the 42 ACS Paragon Plus Environment

Page 43 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

mathematical stiffness of the ODE, the dimension of the original system and the computation time. Simulations in the copolymerization of BMA y BA via ARGET ATRcoP show a long preequilibrium period and large differences in the reactivity ratios, resulting in a prediction of copolymer composition with differences with respect to the ML equation, these discrepancies being higher in the first stage of the process (until 60% of conversion). Also, the solution of the PHP model exhibited a good fitting with the experimental data previously presented with minimal parameter fitting. The most critical parameters to fit the available experimental data were found to be kp22, kact2 and kdeact2, since the other values have been taken from the literature. The model was able to describe the MWD of the process, which was wide at low conversion due to the slow activation of the ATRP initiator. The evolution of MWD in time, corresponding to each polymeric species, offers valuable information about the global behavior of the system; for instance, a maximum concentration of the growing radicals was found at 25 % of conversion, after that the growing radicals gradually fall. In order to increase the polymerization rate, the effect of a previous semi-batch strategy in the average polymerization characteristics and the MWD was analyzed. In this analysis, it was found that the best time to charge the extra addition of the reduced agent to obtain a balance between fast polymerization rate and a preservation of the EGF, is approximately at half of the final conversion. The code developed in this work can be shared with interested research groups to analyze other ARGET copolymerization systems on the basis of a collaboration. The numerical accuracy and the computation times generated by the PHP model turned out to be good, with times in the scales of seconds and minutes. Different computation bounds 43 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 44 of 52

for the maximum chain length of the MWD were studied in the algorithm in order to decrease the total number of equations to solve. The Chebyshev inequality is highly dependent of the safety factor K, producing unnecessary computation of very long chains without significant impact in the MWD when 10 times the standard deviation is used, K= 6 being the best value to minimize the computational time. An empirical expression reported here guarantees the computation of those chain length values with significant polymer concentrations.

AUTHOR INFORMATION Corresponding Author:

Iván Zapata-González Address: Centro de Graduados e Investigación en Química, Tecnológico Nacional de México/ Instituto Tecnológico de Tijuana, Blvd. Industrial and Ave. ITT de Tijuana s/n, Otay, CP 22430 Tijuana, B.C., México. email: [email protected], [email protected] Phone: +52 (664) 6233772, ORCID: 0000-0002-0912-3627

Enrique Saldívar-Guerra Address: Blvd. Enrique Reyna Hermosillo No.140 C.P. 25294 Saltillo, Coahuila México. email: [email protected], Phone: +52 (844) 438-9830, ORCID: 0000-00019623-575X

Author Contributions

44 ACS Paragon Plus Environment

Page 45 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

The manuscript was written through contributions of all authors. All authors have given approval to the final version of the manuscript. §†‡These authors contributed equally.

Notes The authors declare no competing financial interest

ACKNOWLEDGMENTS I. Zapata-González and E. Saldívar-Guerra thank Prof. Robin Hutchinson for his valuable comments and discussions on this work. Also, Ivan Zapata would like to acknowledge generous financial support from CONACYT (Mexico National Council of Science and Technology) through: Cátedras CONACYT project N. 707- 2015, and INFR- 2690072016. E. Saldívar-Guerra also thanks CONACYT for generous grant 256358 of the Basic Science Fund.

NOMENCLATURE Abbreviations ARGET, Activator Regenerated by Electron Transfer; ATRcoP, Atom Transfer Radical copolymerization; ATRP, Atom Transfer Radical Polymerization; BA, butyl acrylate; BMA, butyl methacrylate; CCT, catalytic radical termination; EGF, end-group functionality; ICAR, Initiators for Continuous Activator Regeneration; kMC, kinetic Monte Carlo; LCH, Long Chain Hypothesis; MA, methyl acrylate; MMA, methyl methacrylate; MWD, molecular weight distribution; PHP, Pseudo-Homopolymerization; QSSA, QuasiSteady State Approximation; RDRP, Reversible Deactivation Radical Polymerization; 45 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 46 of 52

RSQSSA, Reduced Stiffness by Quasi-Steady State Approximation; SEC, size exclusion chromatography.

Symbology

φi , fraction of chain radical type i; fi, remaining monomer molar fraction;

vi fraction for dormant species with terminal unit i; Pni , growing chain radicals with terminal unit i and total chain length n; Pni − X , dormant species with terminal unit i and total chain length n;

Dn , dead polymer with total chain length n; M i , monomer i; R − X , Alkyl halide; R * , primary radicals

XCuI/L, activator complex; CuIIX2/L, deactivator complex; Sn r L2 , reducing agent, where r is the oxidation number;

α ki , Moment k of growing chain radical with terminal unit i;

46 ACS Paragon Plus Environment

Page 47 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

γ ki , Moment k of dormant chain with terminal unit i;

βk , Moment k of dead polymer; fint, initiator efficiency; Mn, number-average molecular weight Mn,theo, theoretical number-average molecular weight Mn,exp, experimental number-average molecular weight Mw, weight-average molecular weight Facci, cumulative copolymer composition of monomer i;

Subscript i, monomer unit, 1= BMA and 2=BA

Supporting Information. Model obtained by Method of Moments for ARGET ATRcoP system, derivation of Equation (45).

REFERENCES (1)

Matyjaszewski, K. Atom Transfer Radical Polymerization (ATRP): Current Status and Future Perspectives. Macromolecules 2012, 45 (10), 4015–4039.

(2)

Rudd, J. F. The Effect of Molecular Weight Distribution on the Rheological Properties of Polystyrene. J. Polym. Sci. 1960, 44 (144), 459–474.

(3)

Read, D. J.; Auhl, D.; Das, C.; den Doelder, J.; Kapnistos, M.; Vittorias, I.; McLeish, T. C. B. Linking Models of Polymerization and Dynamics to Predict Branched Polymer Structure and Flow. Science (80-. ). 2011, 333 (6051), 1871–1874.

(4)

Chang, C. J.; Chang, S. J.; Shih, K. C.; Pan, F. U. L. Improving Mechanical Properties and Chemical Resistance of Ink-Jet Printer Color Filter by Using Diblock Polymeric Dispersants. J. Polym. Sci. Part B Polym. Phys. 2005, 43 (22), 3337– 47 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 48 of 52

3353. (5)

Fuchise, K.; Sone, M.; Miura, Y.; Sakai, R.; Narumi, A.; Sato, S. I.; Satoh, T.; Kakuchi, T. Precise Synthesis of Poly(1-Adamantyl Methacrylate) by Atom Transfer Radical Polymerization. Polym. J. 2010, 42 (8), 626–631.

(6)

Boyer, C.; Corrigan, N. A.; Jung, K.; Nguyen, D.; Nguyen, T. K.; Adnan, N. N. M.; Oliver, S.; Shanmugam, S.; Yeow, J. Copper-Mediated Living Radical Polymerization (Atom Transfer Radical Polymerization and Copper(0) Mediated Polymerization): From Fundamentals to Bioapplications. Chem. Rev. 2016, 116 (4), 1803–1949.

(7)

Tobita, H.; Hamielec, A. E. Kinetics of Free-Radical Copolymerization: The PseudoKinetic Rate Constant Method. Polymer (Guildf). 1991, 32 (14), 2641–2647.

(8)

Hernández-Ortiz, J. C.; Vivaldo-Lima, E.; Penlidis, A. Modeling of Network Formation in Nitroxide-Mediated Radical Copolymerization of Vinyl/Divinyl Monomers Using a Multifunctional Polymer Molecule Approach. Macromol. Theory Simulations 2012, 21 (5), 302–321.

(9)

Hernández-Ortiz, J. C.; Vivaldo-Lima, E.; Dubé, M. A.; Penlidis, A. Modeling of Network Formation in Reversible Addition-Fragmentation Transfer (RAFT) Copolymerization of Vinyl/Divinyl Monomers Using a Multifunctional Polymer Molecule Approach. Macromol. Theory Simulations 2014, 23 (3), 147–169.

(10)

López-Domínguez, P.; Hernández-Ortiz, J. C.; Barlow, K. J.; Vivaldo-Lima, E.; Moad, G. Modeling the Kinetics of Monolith Formation by RAFT Copolymerization of Styrene and Divinylbenzene. Macromol. React. Eng. 2014, 8 (10), 706–722.

(11)

Espinosa-Pérez, L.; Hernández-Ortiz, J. C.; López-Domínguez, P.; Jaramillo-Soto, G.; Vivaldo-Lima, E.; Pérez-Salinas, P.; Rosas-Aburto, A.; Licea-Claverie, A.; Vázquez-Torres, H.; Bernad-Bernad, M. J. Modeling of the Production of Hydrogels from Hydroxyethyl Methacrylate and (Di) Ethylene Glycol Dimethacrylate in the Presence of RAFT Agents. Macromol. React. Eng. 2014, 8 (8), 564–579.

(12)

Scott, A. J.; Nabifar, A.; Hernández-Ortiz, J. C.; McManus, N. T.; Vivaldo-Lima, E.; Penlidis, A. Crosslinking Nitroxide-Mediated Radical Copolymerization of Styrene with Divinylbenzene. Eur. Polym. J. 2014, 51 (1), 87–111.

(13)

Saldívar-Guerra, E.; Infante-Martínez, R.; Vivaldo-Lima, E.; Flores-Tlacuahuac, A. Returning to Basics : Direct Integration of the Full Molecular-Weight Distribution Equations in Addition Polymerization. Macromol. Theory Simul. 2010, 19, 151–157.

(14)

Zapata-González, I.; Saldívar-Guerra, E.; Ortiz-Cisneros, J. Full Molecular Weight Distribution in RAFT Polymerization. New Mechanistic Insight by Direct Integration of the Equations. Macromol. Theory Simulations 2011, 20 (6), 370–388.

(15)

Zapata-González, I.; Saldívar-Guerra, E.; Flores-Tlacuahuac, A.; Vivaldo-Lima, E.; Ortiz-Cisneros, J. Efficient Numerical Integration of Stiff Differential Equations in Polymerisation Reaction Engineering: Computational Aspects and Applications. Can. J. Chem. Eng. 2012, 90 (4), 804–823. 48 ACS Paragon Plus Environment

Page 49 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(16)

Zapata-Gonzalez, I.; Saldivar-Guerra, E.; Licea-Claverie, A. Kinetic Modeling of RAFT Polymerization via Dithiobenzoate Agents Considering the Missing Step Theory. Chem. Eng. J. 2017, 326, 1242–1254.

(17)

Gonçalves, M. A. D.; Trigo, I. M. R.; Dias, R. C.; Costa, M. R. P. F. N. Kinetic Modeling of the Molecular Architecture of Cross-Linked Copolymers Synthesized by Controlled Radical Polymerization Techniques. Macromol. Symp. 2010, 291–292 (1), 239–250.

(18)

Sarmoria, C.; Asteasuain, M.; Brandolin, A. Prediction of Molecular Weight Distributions Functions. Can. J. Chem. Eng. 2012, 90, 263–273.

(19)

Wulkow, M. Computer Aided Modeling of Polymer Reaction Engineering- The Status of Predici, 1.- Simulation. Macromol. React. Eng. 2008, 2, 461–494.

(20)

Krys, P.; Matyjaszewski, K. Kinetics of Atom Transfer Radical Polymerization. Eur. Polym. J. 2017, 89 (January), 482–523.

(21)

Fischer, H. The Persistent Radical Effect In “Living” Radical Polymerization. Macromolecules 1997, 30 (19), 5666–5672.

(22)

Jakubowski, W.; Min, K.; Matyjaszewski, K. Activators Regenerated by Electron Transfer for Atom Transfer Radical Polymerization of Styrene. Macromolecules 2006, 39 (1), 39–45.

(23)

Matyjaszewski, K.; Jakubowski, W.; Min, K.; Tang, W.; Huang, J.; Braunecker, W. a; Tsarevsky, N. V. Diminishing Catalyst Concentration in Atom Transfer Radical Polymerization with Reducing Agents. P. Natl. Acad. Sci. USA 2006, 103 (42), 15309–15314.

(24)

Payne, K. A.; Cunningham, M. F.; Hutchinson, R. A. ARGET ATRP of BMA and BA: Exploring Limitations at Low Copper Levels. ACS Symp. Ser. 2012, 1100 (Mmd), 183–202.

(25)

Payne, K. A.; D’hooge, D. R.; Van Steenberge, P. H. M.; Reyniers, M.-F.; Cunningham, M. F.; Hutchinson, R. A.; Marin, G. B. ARGET ATRP of Butyl Methacrylate: Utilizing Kinetic Modeling To Understand Experimental Trends. Macromolecules 2013, 46 (10), 3828–3840.

(26)

Payne, K. A.; Van Steenberge, P. H.; D’hooge, D. R.; Reyniers, M. F.; Marin, G. B.; Hutchinson, R. A.; Cunningham, M. F. Controlled Synthesis of Poly[(Butyl Methacrylate)-Co-(Butyl Acrylate)] via Activator Regenerated by Electron Transfer Atom Transfer Radical Polymerization: Insights and Improvement. Polym. Int. 2014, 63 (5), 848–857.

(27)

Ribelli, T. G.; Augustine, K. F.; Fantin, M.; Krys, P.; Poli, R.; Matyjaszewski, K. Disproportionation or Combination? The Termination of Acrylate Radicals in ATRP. Macromolecules 2017, 50 (20), 7920–7929.

(28)

Al-Harthi, M.; Soares, J. B. P.; Simon, L. C. Dynamic Monte Carlo Simulation of ATRP with Bifunctional Initiators. Macromol. React. Eng. 2007, 1 (1), 95–105. 49 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 50 of 52

(29)

Al-Harthi, M.; Soares, J. B. P.; Simon, L. C. Mathematical Modeling of AtomTransfer Radical Polymerization Using Bifunctional Initiators. Macromol. Theory Simulations 2006, 15 (3), 198–214.

(30)

Schroeder, H.; Buback, J.; Schrooten, J.; Buback, M.; Matyjaszewski, K. Modeling Atom-Transfer Radical Polymerization of Butyl Acrylate. Macromol. Theory Simulations 2014, 23 (4), 279–287.

(31)

Fu, Y.; Cunningham, M. F.; Hutchinson, R. a. Semibatch Atom Transfer Radical Copolymerization of Styrene and Butyl Acrylate. Macromol. Symp. 2007, 259 (1), 151–163.

(32)

Fu, Y.; Mirzaei, A.; Cunningham, M. F.; Hutchinson, R. a. Atom-Transfer Radical Batch and Semibatch Polymerization of Styrene. Macromol. React. Eng. 2007, 1, 425–439.

(33)

Zhang, M.; Ray, W. H. Modeling of “Living” Free-Radical Polymerization Processes. I. Batch, Semibatch, and Continuous Tank Reactors. J. Appl. Polym. Sci. 2002, 86 (7), 1630–1662.

(34)

Toloza Porras, C.; D’hooge, D. R.; Reyniers, M.-F.; Marin, G. B. Computer-Aided Optimization of Conditions for Fast and Controlled ICAR ATRP of n -Butyl Acrylate. Macromol. Theory Simulations 2013, 22 (2), 136–149.

(35)

D’hooge, D. R.; Reyniers, M.-F.; Marin, G. B. Methodology for Kinetic Modeling of Atom Transfer Radical Polymerization. Macromol. React. Eng. 2009, 3 (4), 185– 209.

(36)

Zapata-Gonzalez, I.; Hutchinson, R. A.; Payne, K. A.; Saldívar-Guerra, E. Mathematical Modeling of the Full Molecular Weight Distribution in ATRP Techniques. AIChE J. 2016, 62 (8), 2762–2777.

(37)

Mueller, L.; Jakubowski, W.; Tang, W.; Matyjaszewski, K. Successful Chain Extension of Polyacrylate and Polystyrene Macroinitiators with Methacrylates in an ARGET and ICAR ATRP. Macromolecules 2007, 40 (18), 6464–6472.

(38)

Fierens, S. K.; Van Steenberge, P. H. M.; Reyniers, M. F.; Marin, G. B.; D’hooge, D. R. How Penultimate Monomer Unit Effects and Initiator Influence ICAR ATRP of N-Butyl Acrylate and Methyl Methacrylate. AIChE J. 2017, 63 (11), 4971–4986.

(39)

Wang, W.; Zhou, Y.; Shi, L.; Luo, Z. Modeling of the Atom Transfer Radical Copolymerization Processes of Methyl Methacrylate and 2‑(Trimethylsilyl) Ethyl Methacrylate under Batch, Semibatch, and Continuous Feeding: A Chemical Reactor Engineering Viewpoint. Ind. Eng. Chem. Res. 2014, 53 (30), 11873–11883.

(40)

Steenberge, P. H. M. Van; Dagmar, R. D.; Wang, Y.; Zhong, M.; Konkolewicz, D.; Matyjaszewski, K.; Marin, G. B. Linear Gradient Quality of ATRP Copolymers. Macromolecules 2012, 45, 8519–8531.

(41)

Ballard, M. J.; Napper, D. H.; Gilbert, R. G. Theory of Emulsion Copolymerization Kinetics. J. Polym. Sci. Polym. Chem. Ed. 1981, 19 (4), 939–954. 50 ACS Paragon Plus Environment

Page 51 of 52 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(42)

Saldívar, E.; Dafniotis, P.; Ray, W. H. Mathematical Modeling of Emulsion Copolymerization Reactors . I . Model Formulation and Application to Reactors Operating with Micellar Nucleation. J Macromol Sci-Rev Macromol Chem Phys 1998, 38 (2), 207–325.

(43)

Storti, G.; Carrà, S.; Morbidelli, M.; Vita, G. Kinetics of Multimonomer Emulsion Polymerization. The Pseudo‑homopolymerization Approach. J. Appl. Polym. Sci. 1989, 37 (9), 2443–2467.

(44)

Jakubowski, W.; Matyjaszewski, K. Activator Generated by Electron Transfer for Atom Transfer Radical Polymerization. Macromolecules 2005, 38 (10), 4139–4146.

(45)

Guo, J.-K.; Zhou, Y.-N.; Luo, Z.-H. Kinetic Insight into Electrochemically Mediated ATRP Gained Through Modeling. AIChE J. 2015, 61 (12), 4347–4357.

(46)

Konkolewicz, D.; Sosnowski, S.; D’hooge, D. R.; Szymanski, R.; Reyniers, M.-F.; Marin, G. B.; Matyjaszewski, K. Origin of the Difference between Branching in Acrylates Polymerization under Controlled and Free Radical Conditions: A Computational Study of Competitive Processes. Macromolecules 2011, 44 (21), 8361–8373.

(47)

Tang, W.; Matyjaszewski, K. Effects of Initiator Structure on Activation Rate Constants in ATRP. Macromolecules 2007, 40 (6), 1858–1863.

(48)

Beuermann, S.; Buback, M.; Davis, T. P.; Gilbert, R. G.; Hutchinson, R. A.; Kajiwara, A.; Klumperman, B.; Russell, G. T. Critically Evaluated Rate Coefficients for Free-Radical Polymerization, 3. Propagation Rate Coefficients for Alkyl Methacrylates. Macromol. Chem. Phys. 2000, 201 (12), 1355–1364.

(49)

Wang, W.; Hutchinson, R. A. A Comprehensive Kinetic Model for HighTemperature Free Radical Production of Styrene/Methacrylate/Acrylate Resins. AIChE J. 2011, 57 (1), 227–238.

(50)

Hamzehlou, S.; Reyes, Y.; Hutchinson, R.; Leiza, J. R. Copolymerization of N-Butyl Acrylate and Styrene: Terminal Vs Penultimate Model. Macromol. Chem. Phys. 2014, 215 (17), 1668–1678.

(51)

De Rybel, N.; Van Steenberge, P. H. M.; Reyniers, M. F.; D’hooge, D. R.; Marin, G. B. How Chain Length Dependencies Interfere with the Bulk RAFT Polymerization Rate and Microstructural Control. Chem. Eng. Sci. 2018, 177, 163–179.

(52)

Griffiths, M. C.; Strauch, J.; Monteiro, M. J.; Gilbert, R. G. Measurement of Diffusion Coefficients of Oligomeric Penetrants in Rubbery Polymer Matrixes. Macromolecules 1998, 31 (22), 7835–7844.

(53)

Wulkow, M. The Simulation of Molecular Weight Distributions in Polyreaction Kinetics by Discrete Galerkin Methods. Macromol. Theory Simulations 1996, 5 (3), 393–416.

(54)

McNaught, A. D.; Wilkinson, A. IUPAC. Compendium of Chemical Terminology, 2nd ed.; Blackwell Scientific Publications: Oxford, 1997. 51 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 52 of 52

(55)

Zapata-González, I.; Hutchinson, R. A.; Matyjaszewski, K.; Saldívar-Guerra, E.; Ortiz-Cisneros, J. Copolymer Composition Deviations from Mayo-Lewis Conventional Free Radical Behavior in Nitroxide Mediated Copolymerization. Macromol. Theory Simulations 2014, 23 (4), 245–265.

(56)

Matyjaszewski, K. Factors Affecting Rates of Comonomer Consumption in Copolymerization Processes with Intermittent Activation. Macromolecules 2002, 35 (18), 6773–6781.

(57)

Konkolewicz, D.; Krys, P.; Matyjaszewski, K. Explaining Unexpected Data via Competitive Equilibria and Processes in Radical Reactions with Reversible Deactivation. Acc. Chem. Res. 2014, 47 (10), 3028–3036.

(58)

Guo, J.; Zhou, Y.; Luo, Z. Kinetic Insights into the Iron-Based Electrochemically Mediated Atom Transfer Radical Polymerization of Methyl Methacrylate. Macromolecules 2016, 49 (11), 4038–4046.

TABLE OF CONTENTS AND ABSTRACT GRAPHICS

52 ACS Paragon Plus Environment