Accounting for Spatial Variations during Photopolymerization of 1,6-hexane-diol Diacrylate in the Presence of Oxygen

A dynamic model is proposed for photopolymerization of 1,6-hexane-diol diacrylate (HDDA) with bifunctional initiator bis-acylphosphine oxide (BAPO) in the presence of oxygen. This partial-diﬀerential-equation (PDE) model predicts time-and spatially-varying vinyl-group conversion as well as concentrations of monomer, initiator, oxygen, and seven types of radicals. Experiments to obtain diﬀusivities of oxygen, BAPO and HDDA are reported. Oxygen-related parameters are estimated using real-time Fourier-transform infrared (FTIR) conversion data. FTIR experiments were conducted using


Introduction
Free-radical photopolymerization has been increasing in popularity, especially in the production of films, coatings, and a variety of optical, dental, and microelectronic applications [1][2][3].This article focuses on modeling the photopolymerization of the divinyl monomer 1,6-hexanediol diacrylate (HDDA), which forms a cross-linked network of polymer chains.
The HDDA photopolymerization process involves a complex set of chemical reactions, especially when a bifunctional photo-initiator such as bis-acylphosphine oxide (BAPO) is used [4].Several mathematical models have been developed to provide insights into the photopolymerization kinetics of HDDA and other acrylate monomers [4][5][6][7][8][9][10][11][12][13][14][15].These models have been used to study the influences of experimental factors such as film thickness, light intensity, initiator concentration and oxygen contents [4][5][6][7][8][9][10][11][12][13][14][15].For example, Iedema et al. developed an early model to describe the photopolymerization of HDDA in an oxygen-free environment [11].This model provided predictions of vinyl-group conversion vs. time behavior for HDDA photopolymerization experiments with different light intensities, using a simplified kinetic scheme with only one type of polymeric free radical and one type of vinyl group [11].The most comprehensive oxygen-free kinetic model for HDDA photopolymerization accounts for initiation, propagation, branching, backbiting, cyclization, and termination reactions resulting in four different types of free-radical end groups on the polymer chains and two types of vinyl groups (i.e., vinyl groups on the HDDA monomer and pendant vinyl groups on the polymer) [4].Abdi et al. used conversion vs. time data and weighted least squares (WLS) parameter estimation to fit 30 kinetic parameters in this model.
The experiments used for parameter estimation were conducted using a range of different film thicknesses, light intensities, and initial BAPO concentrations [4,5].Abdi's model produces reliable predictions of vinyl-group conversions but is unable to account for oxygen presence [4].
Oxygen presents an important complication during photopolymerization processes that are in contact with air.Oxygen diffuses into the monomer/polymer mixture and reacts with initiator radicals and carbon-centered polymeric radicals to form peroxy radicals which are much less reactive, thereby inhibiting the polymerization [2,16].Several early models were developed to account for oxygen inhibition in photopolymerization processes [5,[12][13][14][17][18][19][20].One of the earliest was developed by Decker et al. to study oxygen inhibition during the photopolymerization of multi-acrylate systems [19].O'Brien et al. worked on modeling oxygen effects on photopolymerization kinetics while accounting for oxygen, initiator and monomer mass-transfer effects [17].Goodner and Bowman also accounted for simultaneous mass transfer and polymerization, using a free-volume approach to predict changes in diffusivities and rate constants as the polymerization proceeds [20].Andrzejewska et al. developed a model to study the influence of oxygen on termination kinetics, accounting for two kinds of radicals with different propensities to diffuse and react [18].All these models provided new insights into the influence of oxygen inhibition on photopolymerization processes; however, they relied on highly simplified reaction schemes [12,[17][18][19][20]  does not account explicitly for spatial gradients in oxygen and other species in the film as the photopolymerization progresses [5].
The influence of oxygen-induced spatial concentration gradients during photopolymerization of HDDA and other acrylate monomers has been explored by several research groups [12,17,[20][21][22][23][24][25][26][27][28].In thicker films (~100 ), Goodner and Bowman [20] predicted that a thin layer at the top surface has lower monomer conversion due to ongoing oxygen diffusion into the film as polymerization proceeds.However, with thinner films (~10 ), the oxygen diffuses through the entire depth of the film, impeding conversion to various extents at different depths [12].As such, oxygen inhibition may cause polymer films to have varying mechanical properties over their thicknesses [25,26].The model of Iedema et al. provides reliable fits for conversion vs. time data in thin films with low oxygen levels, but overpredicts conversion in an oxygen-rich environment.
Difficulties in predicting experimental data over a wider range of conditions may be due to the simplified reaction mechanism and inaccurate parameter values used in Iedema's model [12].The aim of the current modeling research is to develop an improved model using the comprehensive reaction scheme of Vo et al. with explicit modeling of diffusion within the HDDA polymeric film [5].A further objective is to provide new diffusivity data for oxygen, BAPO and HDDA so that improved estimates can be obtained for model parameters related to mass-transfer and oxygen inhibition.The resulting improved HDDA photopolymerization model will be helpful for selecting appropriate recipes and operating conditions to achieve high-quality rapid printing processes [5].
The remainder of the article is organized as follows.First, new experimental results and statistical analysis are discussed.Second, an updated model is proposed which relies on Vo's reaction mechanism, but explicitly accounts for oxygen, initiator, and monomer concentration gradients in the film.This updated model, which relies on an assumption of constant film thickness, uses partial differential equations (PDEs) to describe concentrations of reacting species and functional groups, which vary with time and position in the film.Six oxygen-related and diffusion-related parameters are then estimated based on conversion versus time data from thirty-two experimental runs.The predictive ability of the model is assessed using additional data not utilized for parameter estimation.Finally, recommendations are made regarding a future model that will account for changes in density and film thickness during the photopolymerization process.

Experimental Data and Statistical Analysis 2.1 Measurement of the diffusion coefficient of oxygen in HDDA
Pure HDDA films with a thickness of 1.5  were used to determine the oxygen solubility and oxygen diffusion coefficient in HDDA monomer.The film was first exposed to a nitrogen gas flow for approximately 90 minutes so that all atmospheric oxygen was purged from the system.After 90 minutes, air was fed to the chamber where the HDDA film is located.Oxygen concentration at the bottom of the film was measured using an optical oxygen sensor in a custom-built enclosure, which is shown schematically in Fig. S1, in the Supplementary Information.
A plot of the oxygen concentration vs time is shown in Fig. 1.The data shown in Figure 1 are responses from two replicate runs.Whiskers on the error bars correspond to one standard deviation, which was computed from a pooled variance estimate.The curve shown in Figure 1 is a plot of: which is the analytical solution of the oxygen material balance PDE, were obtained from the data in Fig. 1 using least squares regression.

Measurement of self diffusion coefficient for HDDA and diffusion coeffiecient of BAPO in HDDA
Diffusion coefficients of HDDA and BAPO were measured using diffusion-ordered NMR spectroscopy (DOSY).DOSY NMR works by analyzing the results of a pulsed field gradient spin echo experiment to obtain the diffusion coefficients of individual signals in a spectrum [30].In the experimental setup, a solution containing 1 % BAPO was placed in a 400 MHz NMR spectrometer (magnetic field B0) and exposed to a linear magnetic field gradient ( → ) which spatially encodes the HDDA and BAPO molecules.After a diffusion delay the molecules are similarly decoded.The diffusion coefficient can then be calculated from the signal attenuation of each molecule.An encoding/decoding step with duration  = 3000  and a diffusion step with duration ∆ = 0.2  were used to obtain NMR peaks, which resulted in calculated diffusion coefficients   = 1.59 × 10 −10  2  and   = 8.78 × 10 −11  2  .

Measurement of overall vinyl-group conversion in HDDA polymer films
Overall vinyl-group conversion vs. time data were collected using Fourier-transform infrared reflection (FTIR) spectroscopy.The setup of the experiments used to obtain these data consists of a pure HDDA film with added BAPO initiator and Tego Rad 2250 surfactant (2 %).These films were exposed to the air so that the initial oxygen concentration was uniform throughout.The The additional 27 runs used in the current study contain new information.

Reaction mechanism
The current modeling study relies on the reaction mechanism of Abdi et al. and the extension proposed by Vo et al. to account for oxygen inhibition [4,5].A detailed list of reactions is provided in Table S1 in the Supplementary Information.The first step in the mechanism is light-induced decomposition of one of the carbon-phosphorus bonds of the BAPO initiator (see Fig. 2).
Decomposition of BAPO results in two types of radicals, a carbon-centered  • radical and a phosphorous-centered  ̃• radical.The phosphorous-centered radical has a remaining carbonphosphorous bond that can be decomposed in a later reaction to produce   • (see Fig. 2b).Initiator radicals can react with vinyl groups on the HDDA monomer generating  =• .This macroradical end group  =• has a free radical and an unreacted vinyl group (see Fig. 3b).We assume that the life of  =• is too short for the vinyl groups in  =• to be consumed by branching before the corresponding free radical reacts (assumption 4 and assumption 7 in Table S1 in the Supplementary Information).
The end group  =• can propagate with HDDA monomer (Fig. 3a).When this reaction happens, the vinyl group from the  =• whose radical is consumed becomes a pendant vinyl group   .End groups of type  =• can also undergo cyclization reactions to produce the cyclic radical  • (Fig. 3b).Free radical end groups such as  =• and  • can also participate in branching reactions by propagating with pendant vinyl groups   (Fig. 3c) to produce a branch-point radical  • .
Backbiting reactions (Fig. 3d) produce a tertiary radical  • .In addition to these reactions, all radical species can react with oxygen to produce peroxidic radical end groups  • (Fig. 3e).
Peroxidic radicals undergo very slow propagation, branching and back-biting reactions.Finally, Peroxidic radicals can terminate with other radicals in the reaction mixture.Consumption of  =• ,  • ,  • and  • radicals via termination with  • radicals lead to an important reduction in the rate of consumption of vinyl groups.

Model equations
The proposed model accounts for spatial variation in concentrations of chemical species and end groups, leading to 13 partial differential equations (PDEs) as shown in Table 3   )] [5] Table 3: PDEs based on material balances for chemical species and end-groups, including steady state hypothesis where it applies.

𝜕
which becomes, after applying SSH: which becomes, after applying SSH:

Parameter estimation
The proposed model requires 39 parameters to simulate the photopolymerization process.
Fortunately, 30 of these parameters were previously estimated by Abdi et al. using their oxygenfree model (see Table 4) [4].Three of the 9 new parameters are diffusivities ( 2− ,  − and  − ) which are measured directly in the current study (see Table 5).
The remaining six parameters which require estimation are free-volume parameters and oxygenrelated kinetic parameters whose initial guesses are provided in Table 6.Initial guesses for and are values obtained from Vo's parameter estimation study [5].They are re-estimated in the current study to account for the improved model structure and information in the new data.Upper and lower parameter bounds in Table 6 are specified to ensure that all estimated parameters are physically realistic.
The least-squares objective function used for parameter estimation is: where    , is the i th measured value of overall vinyl-group conversion,  , is the corresponding model prediction and    = 1,655 is the number of measured values used for parameter estimation (obtained from 32 experimental runs).Prior to conducting the parameter estimation, a formal estimability analysis was performed to confirm that all six parameters could be estimated uniquely using the available data [34,35].
The resulting parameter estimates are shown in the fourth column in Table 6.Notice that the estimate of the initiator free-volume parameter   is zero.The optimizer moved this parameter to its lower bound because accounting for the influence of conversion on initiator diffusion did not improve the model predictions.This result is not surprising because the initiator is consumed rapidly, so there is not enough time for considerable initiator diffusion to occur.As a result, the diffusivity of the initiator has little influence on the model predictions.Also notice that the estimated oxygen free-volume parameter   2 is lower than the estimated HDDA free-volume parameter   .This result makes physical sense because oxygen molecules are smaller than HDDA and are less impeded by reductions in free volume as polymerization proceeds.

Simulation results and model assessment
Fig. 4 compares the model predictions with the experimental results for the two runs that had the best fit to the data (in green) and the worst fit to the data (in red), respectively, according to their contributions to their least-squares objective function.The best-fit simulation results in green correspond to a relatively thick film (~17 ), a medium light intensity (~1180 / 2 ) and a high BAPO level (4 %).The worst-fit simulation results in red are obtained for a run with similar film thickness (~20 ) and light intensity (~1200 / 2 ), but with a low BAPO level (1 %).Plots showing the fit to the data for several additional runs used for parameter estimation are provided in the Supplementary Information.In general, plots for runs with low BAPO (e.g.see Fig. S2) reveal that the model tends to underpredict the overall conversion when initial BAPO concentration is low.This result is not surprising because most of the kinetic parameters used in the simulations were estimated by Abdi et al., using a data set where all runs were conducted using a high BAPO concentration of 4 %.Abdi noted that there is considerable uncertainty associated with many of his model parameters and recommended future parameter estimation studies using experiments conducted over a wider operating range.Fortunately, the runs shown in Fig. 4 confirm that the model gives reliable qualitative predictions (i.e., both the model predictions and the data are in agreement that lower BAPO levels lead to lower vinyl-group conversion).Simulation and experimental results in Fig. S3, which compare the results for two runs with different film thicknesses, confirm that thicker films tend to result in higher conversions because they experience less oxygen penetration during the runs.
Fig. 5 shows the model's predictive ability using two of the three runs saved for validation.As expected, the model predicts correctly that as light intensity increases the overall vinyl-group conversion also increases.Unfortunately, the model overpredicts the conversions in the two runs shown in Fig. 5. Fig. S4 shows that the model slightly overpredicts the conversion for the third validation run, which has lower film thickness and lower light intensity than the runs in Fig. 5.
In general, the proposed model performs better at fitting the objective function than predictions from Vo's ODE oxygen-inhibition model (i.e., the value of the objective function in Equation ( 2) is 49% lower).This is not surprising, since this model uses more data to get more accurate parameter estimates and uses PDEs to better account for oxygen and monomer diffusion.when all the reactions cease, the oxygen concentration becomes uniform within the film.Fig. 7 shows the total vinyl-group concentration changing with time at different depths within the film.
It confirms that at deeper levels within the film, the total vinyl-group concentration is lower than the concentration near the surface, because less oxygen inhibition occurs at the bottom of the film.In future, it will be important to re-estimate some of the kinetic parameters obtained by Abdi et al which were held constant during the current modelling study.We are hopeful that improved parameter estimates will lead to better predictions for experiments with low BAPO levels and light

Conclusion
intensities.The current model assumes that the film has constant density during polymerization, which is not the case.As photopolymerization proceeds, film shrinkage on the order of 10% can occur because HDDA polymer is denser than the monomer.An updated model that accounts for shrinkage is currently under development and will be used in future parameter estimation studies.
The measured diffusivities and parameter estimate from the current study will be used as initial guesses for this ongoing research, which will produce a more accurate model that can be used reliably over a wider range of conditions of industrial interest.

Fig. 1 .
Fig. 1.Measured oxygen concentration (▲) and model predictions (-) using Equation(1) experiments are divided into two sets.The first data set arises from 12 experimental runs that had a film thickness of ~12 , an initial BAPO concentration of 4 %, and light intensities ranging from 200 / 2 to 6000 / 2 .The second data set arose from an attempted factorial design involving three factors with two or three replicates for every combination, leading to 23 experimental runs.The factors film thickness (8 or 16 ), initial BAPO concentration (1 or 4 %), and light intensity (1000 or 5000 / 2 ).All experiments were conducted in the presence of atmospheric oxygen during polymerization.Thirty-two out of the 35 experimental runs are used for estimating the model parameters in this study, while the remaining 3 are saved for model validation.Details about the photopolymerization method and data collection are provided elsewhere [5, 12, 31].Note that, Vo et al. used 8 out of the 35 runs for parameter estimation in their previous modelling study (i.e., runs involving thin films with relatively low light intensities).

Fig. 4 .
Fig. 4. Model predictions (curved lines) and the corresponding FTIR measured values for the experiments of those runs.

Fig. 6
Fig. 6 shows predicted changes of oxygen concentrations with time at different depths in the film for a run conducted with a thick film (~18 ), high light intensity (5916 / 2 ) and a high

Fig. 6 .Fig. 7 .
Fig. 6.Oxygen concentration levels with respect to time at different film depths

Table 1 -
Initial and boundary conditions for solving Fick's second law.

Table 2 )
. Equation 3.1 is a material balance on the unreacted initiator where [] is the initiator concentration in   ⁄ ,  is the time in , and  is the distance from the top surface of the film in .The first term on the right-.Initial and boundary conditions for Equation 3.1 and the other PDEs in Table3are provided in TableS3in the Supplementary Information.Equations 3.2 to 3.10 are dynamic material balances on unreacted initiator fragments  ̃, monomer , pendant vinyl groups   , dissolved oxygen, and free radicals of types  =• ,  • ,  • ,  • and  • , respectively.Notice that the balances on monomer, initiator and oxygen contain diffusion terms, but balances on the polymeric species (i.e.,  ̃,   ,  =• ,  • ,  • ,  • and  • ) do not (see assumption 3.1 and 3.8 are provided in the Supplementary Information.The other PDEs in Table 3 are derived in a similar fashion.
[12] side accounts for initiator diffusion within the film, and the second term accounts for lightinduced initiator consumption.The diffusivity   depends on the free volume   within the film, which depends on the monomer conversion  [12].Algebraic equations for   and other variables are provided in Table2.Notice that the second term in Equation3.1 contains a factor of two since BAPO is a bifunctional initiator.The initiator decomposition rate constant   depends on light intensity   which decreases gradually with the film depth  (see Equations 2.1 and 2.2 in

Table 4 :
[4]t of parameters and their estimates obtained from Abdi's model[4]

Table 5 :
List of experimentally obtained parameters.

Table 6 :
List of estimated parameters