Kinetics Analysis of Consecutive Reactions Using Nelder-Mead

21 Kinetics Analysis of Consecutive Reactions Using Nelder-Mead Simplex Optimization Gary M . Carlson and Theodore Provder 1

2

Downloaded by COLUMBIA UNIV on June 19, 2013 | http://pubs.acs.org Publication Date: June 27, 1986 | doi: 10.1021/bk-1986-0313.ch021

Glidden Coatings and Resins, S C M Corporation, Strongsville, O H 44136

The Nelder-Mead simplex minimization routine is used to determine kinetics constants for a reaction system involving a first order generation of reactive intermediate B followed by a second order consumption of B by coreaction with species C. Experimentally, the concentration of B is monitored as a function of time. Theoretical curves describing the concentration of B are generated using an Euler-Romberg solution of the differential equations. The Nelder-Mead simplex minimization routine optimizes the values of k1 and k2 to obtain the best agreement between theoretical curves and actual experimental data. The application of this technique to the isothermal deblocking and cure of blocked isocyanate containing materials monitored by Fourier transform infrared spectroscopy is demonstrated.

Crosslinking of many polymers occurs through a complex combination of consecutive and p a r a l l e l reactions. For those cases i n which the chemistry i s well understood i t i s possible to define the general reaction scheme and thus derive the appropriate d i f f e r e n t i a l equations describing the cure k i n e t i c s . A n a l y t i c a l s o l u t i o n s have been found for some of these systems of d i f f e r e n t i a l equations permitting accurate experimental determination of the i n d i v i d u a l rate constants. Often a n a l y t i c a l solutions do not e x i s t , or have not been determined. I t i s sometimes p o s s i b l e , i n these cases, to a r r i v e a t approximate expressions or to define conditions such that rate constants or functions of the rate constants can be determined experimentally. 'Current address: Research and Development, Ashland Chemical, Columbus, O H 43216 Author to whom correspondence should be addressed

2

0097-6156/86/0313-0241$06.00/0 © 1986 American Chemical Society

In Computer Applications in the Polymer Laboratory; Provder, T.; ACS Symposium Series; American Chemical Society: Washington, DC, 1986.

242

COMPUTER APPLICATIONS IN THE POLYMER LABORATORY

The use of blocked isocyanates to cure hydroxyl containing coatings i s an example of a complex system having many p r a c t i c a l applications. The chemistry of blocked isocyanates has been reviewed previously by Wicks (J.,.2). The cure reaction proceeds v i a consecutive f i r s t and second o7der reactions: H

k

1

R-NCO-B

>

k


R-NCO

+ ROH

R-NCO + BH

(1)

2 >

R-NCOOR

(2)

Deblocking l i b e r a t e s free isocyanate which reacts with hydroxyl f u n c t i o n a l i t y i n a second order reaction to form the actual c r o s s l i n k s . Anagnastou and Jaul have studied the deblocking o f MDI derivatives using d i f f e r e n t i a l scanning calorimetry ( 3 K The deblocking o f coatings systems containing the blocked isocyanates has been investigated previously by Kordomenos and coworkers (4) using thermogravimetric analysis. Previous investigations (5-8) using Fourier transform infrared analysis to determine the chemical changes which occur during l i n e a r heating supported the proposed scheme and showed that, i n the absence o f c a t a l y s t s , substantial amounts of free isocyanate were generated. Other investigators have also observed similar behavior in similar coatings systems (9, 10). The e f f e c t i v e development o f blocked isocyanate based coatings requires a complete understanding of the cure chemistry. Many materials have been i d e n t i f i e d which w i l l improve cure response. Often these are reported as deblocking catalysts even though no d i r e c t evidence exists to support t h i s claim. Most of these materials are well known catalysts o f the cure reaction between hydroxyl and isocyanate and may be improving cure response s o l e l y by catalyzing t h i s reaction. E f f e c t i v e development o f catalyst systems requires a better understanding o f the e f f e c t o f catalysts on cure. Modification of the FT-IR analysis techniques to analyze coatings under isothermal cure conditions provides the data needed to determine the rate constants for each reaction. An e f f e c t i v e method to generate the rate constants from the experimental data has been found and w i l l be described. Using the previously proposed chemical scheme (Equation 1-2) and assuming that there i s only one isocyanate species, one hydroxyl species, and that the deblocking i s i r r e v e r s i b l e ( i . e . , the blocking agent leaves the f i l m immediately) then the concentration o f each species as a function of time i s given by the following d i f f e r e n t i a l equations: -d[R-NC0(B)]/dt

=

k CR-NC0(B)]

(3)

-d[R-NC0]/dt

=

K [R-NC0][0H] - kjCR-NCOCB)]

(4)

-d[0H]/dt

=

k CR-NC0][0H]

(5)

=

k [R-NC0][0H]

(6)

d[R-NC00R]/dt

1

2

2

2



21. CARLSON AND PROVDER

Kinetics Analysis of Consecutive Reactions

243

In the absence of simplifying assumptions a non-closed form a n a l y t i c a l solution of these equations giving the concentration o f each species as a function of time has been accomplished through the use of Bessel functions (11). Simplifying assumptions, such as a large hydroxyl excess, are acceptable t h e o r e t i c a l l y but are not h e l p f u l i n determining the rate constants of coatings materials with standard compositions. As seen i n Equation 4, the concentration of free isocyanate as a function of time i s dependent on both rate constants. Measuring the concentration of the isocyanate as a function of time during isothermal cure provides a convenient method to determine both rate constants. Many p r a c t i c a l applications of cure characterization involve samples f o r which the data required to convert isocyanate absorbance to concentration i s unavailable. The emphasis i s often placed on rapid analysis of many samples rather than an exhaustive characterization of a single sample. I t i s p a r t i c u l a r l y desirable to develop a procedure which can determine the rate constants describing the cure reaction without converting the infrared absorbance curve to concentration. This has been accomplished by normalizing the data i n such a way that the rate constants are determined from the shape of the cure curve. The experimental technique used to obtain the cure curves and a detailed examination of the procedure used to generate the rate constants, an i t e r a t i v e procedure to determine the best values of k^ and kp based on the set of d i f f e r e n t i a l equations describing the cure process, for these systems w i l l be discussed. Experimental Spectra were obtained using a Digilab FTS-15E Fourier Transform Spectrophotometer. A NaCl c r y s t a l mounted i n a heated c e l l (Model #018-5322 Foxboro/Analabs, N. Haven, Ct.) was placed i n the infrared beam and the chamber allowed to purge f o r several minutes while the c e l l was brought to the desired temperature. The temperature of the c e l l was controlled using a DuPont 900 D i f f e r e n t i a l Thermal Analyzer interfaced to the spectrometer c e l l . A chlorobenzene solution (ca. 10% by wt.) of the sample was then applied to the c r y s t a l using cotton tipped wood s p l i n t . Immediately after applying the sample, spectral data acquisition was started. The small size of the sample, i n r e l a t i o n to the heated c e l l assembly, assures nearly instantaneous e q u i l i b r a t i o n at the desired temperature. Twenty scans were co-added to produce one spectrum at a resolution of 4 cm" every 30 seconds. The reaction was followed by monitoring the absorbance of isocyanate as a function of time. Film thickness changes were compensated f o r by normalizing the isocyanate absorbance to the 1446 cm" band which remains constant i n absorbance during the reaction. The infrared absorbance of the free isocyanate i s converted to concentration by comparison to the absorbance observed i n a similar sample which has no functional groups with which the isocyanate can react.


244



Results and Discussion The concentration o f each chemical species, as a function o f time, during cure can be calculated numerically from Equations 3-6 using the Euler-Romberg integration method i f the i n i t i a l concentrations of blocked isocyanate and hydroxyl f u n c t i o n a l i t y are known. I t i s a s e l f - s t a r t i n g technique and i s generally well behaved under a wide v a r i e t y of conditions. D e t a i l s of t h i s numerical procedure are given by McCalla 0 2 ) . In the case under study, the rate constants are unknown but a t r i a l and error approach can be used t o determine the rate constants which give the best agreement t o the experimental data. The NelderMead simplex algorithm i s used to choose the t r i a l and error values (13)* This algorithm has been shown t o be very u s e f u l i n other r e l a t e d experiments ( J 4 ) . The procedure s t a r t s using an i n i t i a l guess for the i n d i v i d u a l rate constants k« and kp. The Euler-Romberg integration method i s then used t o generate concentration curves based on the t r i a l rate constants. An objective function describing the difference between the experimental and t r i a l curves i s calculated and returned t o the Nelder-Mead simplex routine. This objective function i s then used, i n combination with the other previous best guesses t o a r r i v e at a new p a i r of rate constants. The process i s repeated u n t i l the values of k- and kp which best describe the experimental data have been i d e n t i f i e d . The complete process i s i l l u s t r a t e d i n Figure 1. The simplest case for determining rate constants occurs when the s t a r t i n g concentrations o f blocked isocyanate and hydroxyl f u n c t i o n a l i t y are accurately known. In t h i s case the objective function F, i s defined as

F

=Z

f

( [NC03

-

e

[NC03 )

(7)

c

where t i s the i n i t i a l time, t i s the f i n a l time, [NCO] i s the experimental isocyanate concentration at time t and [NCO] i s the calculated concentration at time t . The effectiveness of the above described computational procedure was tested by generating an a n a l y t i c a l ("ideal") data curve by c a l c u l a t i n g the Isocyanate concentration as a function of time assuming rate constants of k = 1.0/min and kg = 1.0 i/mol/min and i n i t i a l concentrations of blocked isocyanate and hydroxyl o f 1.0 M. The objective function, for various values of k- and kp, for t h i s " i d e a l " data was calculated and a contour plot f o r constant values of F was generated and i s shown i n Figure 2. This plot i s a map of the surface over which the search f o r the best values occurs. The solution i s located i n the objective function 'well' at k. = 1.0/min and kp = 1.0 £/mol/min. The path taken by the program when the i n i t i a l guesses of k = 1.5/min and = 1.5 fc/mol/min also i s shown i n Figure 2. This s t a r t i n g point i s located on a s l i g h t slope and begins i t s descent towards the minimum value with only one major deviation, occcurring a f t e r the eighth i t e r a t i o n . The contour p l o t also demonstrates that no l o c a l minima occur over the region i n which reasonable answers occur. f

1

1



CARLSON AND PROVDER

INPUT DATA 1

>

SEL ECT INITIAL k k 1t

s

NELDER MEAO Tl^LE)TMTNli^Z^fibNl


CALCULATE TRIAL CURE CURVE

FIGURE 1

KINETICS ANALYSIS FLOW CHART

1.0

FIGURE 2 CONTOUR PLOT AND SEARCH TRAJECTORY FOR TEST DATA (INITIAL ESTIMATE: k = 1.5 /min, k = 1.5 A/mol/min) 1

2



246


The path taken when the search i s started at k. = 0.5/min and kps0.5 A/mol/min i s shown i n Figure 3. This search begins on a steeper slope and the function, i n t h i s case, s p i r a l s i n on the minimum value. Again the correct values of k- = 1.0/min and k^ = 1*0 £/mol/min are found. Total computer time ror t h i s analysis averaged 3-5 minutes with a D i g i t a l Equipment Corporation PDP-11/44 minicomputer. An experimentally determined curve describing the concentration of free isocyanate during isothermal curing of a blocked isocyanate containing coating i s shown i n Figure 4. From the minimization procedure the rate constants were determined to be k^ = 1.86 /min and kp s 0.35 Jt/mol/min. The predicted concentration i s shown by the s o l i d l i n e labeled by squares i n Figure 4. I t i s also possible knowing the rate constants to calculate the concentration of blocked isocyanate and the concentration of crosslinks as a function of time. Interestingly, i t i s seen i n Figure 4 that the maximum concentrat i o n of free isocyanate observed during the deblocking corresponds to 69.2 percent of the t o t a l isocyanate. The blocked isocyanate has decreased to a negligible l e v e l after only 3 minutes but at the time of complete deblocking only about one-half of the t o t a l potential crosslinks have formed. In a coating such as t h i s a better catalyst for deblocking w i l l have l i t t l e e f f e c t on the o v e r a l l cure rate of the coating. A more e f f e c t i v e choice for catalyst would be a catalyst for the cure reaction. The contour plot for t h i s experiment i s shown i n Figure 5. The shape of the surface i s somewhat d i f f e r e n t than for the 'ideal* data. This i s partly a result of the experimental error which i s present i n the data and raises the value of the objective function for any t r i a l set of rate constants. There also i s a difference i n the shape of the surface. The contours for the ' i d e a l ' data, shown i n Figures 2 and 3, described a deep c i r c u l a r well with nearly concentric c i r c l e s surrounding i t . The surface of the data i s more elongated i n the k. axis. Thus a wide range of k. values w i l l permit the answer to f a l l within a given objective function value. The rate constant for cure has a smaller range of acceptable answers. This gives an estimate of the margin of error for each determined rate constant. In t h i s data, kp i s known, to a better degree than k-. The objective function, F, as a function of i t e r a t i o n , during the search i s shown i n Figure 6. The objective function decreases rapidly and has neared the minimum value after approximately 20 i t e r a t i o n s . The t r i a l values for the rate constants are shown i n Figure 7. The values for kp l e v e l o f f after about 30 i t e r a t i o n s , while k«j, which has less e f f e c t on the value of the objective function, does not l e v e l o f f u n t i l about 60 i t e r a t i o n s . The convergence c r i t e r i a i s not met u n t i l almost 120 i t e r a t i o n s . Since the rate constants are e s s e n t i a l l y determined after 60 i t e r a t i o n s , the computer time could be decreased without a f f e c t i n g the accuracy of the f i n a l answer by using a less stringent c r i t e r i a . The residual plot of the difference between the best computed curve and the experimental data i s shown i n Figure 8. The largest variance i s observed i n the v i c i n i t y of the maximum isocyanate absorbance and probably arises because of the large changes i n concentration which arise during the early stages of cure. Using the previously described procedure, determination of the rate constants can be made when an accurate curve describing the




CARLSON AND PROVDER

FIGURE 3 CONTOUR PLOT AND SEARCH TRAJECTORY FOR TEST DATA (INITIAL ESTIMATE: k = 0.5 /min, k = 0.5 A/mol/min) 1

2

10H

z g

< 0C

§0.5-1

o o < -I

o

2 0.0H 10

TIME (min)

—i— 15

20

FIGURE 4 EXPERIMENTAL DATA AND CALCULATED ISOTHERMAL CURE CURVES (A PREDICTEDD [BLOCKED ISOCYANATE], 0 EXPERIMENTAL [NCO], •PREDICTED [NCO], + PREDICTED [CROSSLINKS])

American Chemical Society Library 1155 16th St., N.W. Washington, D.& 20036




1.4

1.8

FIGURE 5

2.2

2.6

CONTOUR PLOTS FOR EXPERIMENTAL DATA

0.3n

40 FIGURE 6

80 ITERATIONS

120

OBJECTIVE FUNCTION TRAJECTORY DURING MINIMIZATION



CARLSON AND PROVDER


40 80 ITERATIONS FIGURE 7

120

RATE CONSTANT TRAJECTORY DURING MINIMIZATION

TIME (min) FIGURE 8

RESIDUAL FOR FITTED CURVE




250

isocyanate concentration as a function of time can be obtained experimentally. This requires either detailed knowledge of the resin system, including the molar absorbtivity of the isocyanate group, required to convert the absorbance readings to actual concentration* It i s also necessary to know the i n i t i a l hydroxyl and blocked isocyanate concentrations. In many cases, t h i s data i s unavailable. It would be p a r t i c u l a r l y useful to determine the rate constants for each step using the isocyanate absorbance versus time rather than concentration as a function of time. Indeed, since the absorbance i s e s s e n t i a l l y an a r b i t r a r y measure, the curve can be measured i n any convenient units such as FT-IR chart height. E s s e n t i a l l y , the procedure bases i t s f i t on the shape of the curve rather than on the actual magnitude of the absorbance. The key to using absorbance rather than concentration i s the normalization of the t r i a l curve and the experimental curve. This i s done by d i v i d i n g each point of the t r i a l curve by an a r b i t r a r y value. This a r b i t r a r y value i s the measured height on the t r i a l curve of the isocyanate absorbance corresponding to the time when the experimental data reaches i t s maximum. In t h i s way a curve i s obtained which has a maximum value of 1 . 0 i n dimensionless u n i t s . Since the units of the measurement cancel, the measurement can be done using any convenient u n i t s . Each t r i a l curve generated by the Euler-Romberg integration method i s s i m i l a r l y normalized by dividing each of i t s points by a value corresponding to the time at which the experimental curve reached i t s maximum. f(t)

new

=

f(t> old

/f(tmax)

(8)

where f ( t ) i s the normalized value of the calculated curve and f(tmax) i s Ine height of the t r i a l curve at the time when the experimental data had i t s maximum value. This process forces each curve to have an a r b i t r a r y concentration value of 1 . 0 at the time i n which the maximum experimental absorbance was measured. Figure 9 shows a t h e o r e t i c a l 'data' curve for the case where k. s 1 . 0 /min and L s 1 . 0 #mol/min. A t r i a l curve also i s depicted where k^ = 1 . 5 /min and kp = 1 . 5 Mnol/min. After normalization the two curves are transformed to those shown i n Figure 10. The ' t r i a l * curve now has an a r b i t r a r y maximum greater than 1 . 0 but passes through the point at which the 'data' curve was maximum. The e f f e c t of t h i s normalization procedure can be seen i n the contour plot of Figure 1 1 . The minimum, rather than being a 'well' as i n the procedure based on concentration now i s more of a 'valley' in which a wide range of values of k and kp w i l l provide reasonable solutions to the equation. Values for k or from . 8 to 1 . 3 /min and for k of from . 5 to 1 . 5 Vmol/min can result i n answers with F = 0.005. The trajectory of the minimization procedure i s shown i n Figure 1 1 . The function rapidly finds the valley f l o o r and then travels through the valley u n t i l i t reaches the minimum. A similar trajectory i s shown i n Figure 1 2 i n which the search i s started from a d i f f e r e n t point. In the case of " i d e a l " data the procedure w i l l s t i l l find the minimum along the valley f l o o r . A better understanding of the nature of t h i s 'valley' i s obtained by looking at other points which l i e within the 'valley'. The rate constants k- s 1 . 2 5 /min and k- s 0 . 7 Mnol/min also l i e within the 1

1

2


CARLSON AND PROVDER



•48r

8 16 TIME (MIN.)

FIGURE 9 TYPICAL TRIAL CURVE BEFORE NORMALIZATION ( D k = 1.5/min, k = 1.5 £/mol/min TRIAL CURVE; 0 k = 1.0/min,V 1.0*/mol/mXn DATA CURVE) 1

p

1

1.2 O

TIME (MIN.)

FIGURE 10 NORMALIZED TRIAL CURVE ( O k - = 1.5/min, k = 1.5*/mol/min TRIAL CURVE; 0 L = 1.0/min, k = 1.0£/m5l/min DATA CURVE) p

9

d

d




FIGURE 11 CONTOUR PLOT AND SEARCH TRAJECTORY USING NORMALIZATION (INITIAL ESTIMATE: k = 0.5 /min. kg = 1.5 j/mol/min) 1

FIGURE 12 CONTOUR PLOT AND SEARCH TRAJECTORY USING NORMALIZATION (INITIAL ESTIMATE: k = 1.5 /min, kg = 0.5 $/mol/min) 1




21. CARLSON AND PROVDER

253

v a l l e y where the objective function F

Kinetics Analysis of Consecutive Reactions Using Nelder-Mead

Recommend Documents