Photopolymerized additive manufacturing materials: Modeling of the printing process, mechanical behavior, and sensitivity analysis

The physical – chemical processes involved in light-induced polymerization (photopolymerization) are widely exploited in additive manufacturing (AM) technologies such as Stereolithography and Digital Light Processing. The influence of the AM process parameters on the physical properties of manufactured components has been often investigated through empirical methods based on the trial and error approach, that is, by collecting and interpreting a large amount of experimental data. However, when desired physical properties are required, accurate modeling of the liquid – solid conversion is necessary. In this work, in order to determine the properties of the resulting material according to the adopted process setup, we present a multi-physics approach to model the physical – chemical transformation taking place in photopolymerization. The role played on the final mechanical properties by the laser light intensity and by its moving speed is considered. Further, the influence of the uncertainty of the process parameters is investigated through a sensitivity analysis. The proposed approach is suitable for investigating the reliability of additively manufactured components as well as for their design according to an optimum printing strategy. From the perspective of making innovative functional materials, the proposed multi-physics model allows tuning the printing process in order to get the desired distribution of mechanical properties.

subsequently developed, such as the Digital Light Processing (DLP), the Projection Micro-Stereolithography, and the Microstereolithograpy. 1,2 All these technologies enable producing components to be used in a wide range of applications, such as in medicine, bio-engineering, soft-robotics, automotive, aerospace, and nanotechnology.
As a matter of fact, the photopolymerization process allows getting polymeric components characterized by intricate geometries even with micro-size dimensions, a result almost impossible to be achieved with standard technologies. Moreover, differently from standard technologies based on a subtractive manufacturing, it enables minimizing the waste material. 3 An open issue still existing for all the AM technologies is represented by the possibility to tune the mechanical features of the printed material not only on the basis of the raw material used (consisting of a liquid monomer in the case of photopolymerization) but also on the specific way the solidification process, driven by the AM process parameters, is conducted. 4 In fact, the same component can show very different properties according to the printing setup adopted to fabricate it. 1,2 The existence of these "degrees of freedom" in producing a photopolymerized material, introduces new complexities that need to be fully controlled in order to fulfill the requirements related to the application in turn.
Traditionally, the influence of the printing setup on the mechanical behavior of AM materials has been investigated by means of empirical approaches, mostly based on the collection of a large amount of experimental data, [5][6][7] aimed at correlating the final mechanical properties to the AM parameters used. Despite these approaches can help understanding the relation between AM process parameters and achieved mechanical features, their use is limited to the specified material studied and the obtained relations are often limited to specific conditions, thus do not have the required generality.
In order to face this problem, understanding and modeling the chemical-physical phenomena in AM processes is crucial for the use of the AM production in advanced applications requiring a precise control of the final materials' properties. This approach, besides being useful in predicting the mechanical behavior of AM components, more intriguingly allows designing the printing process strategy according to the desired final properties and their distribution within the component.
The paper is organized as follows: in Section 2 we present the multiphysics model for the photopolymerization process, illustrating in Section 2.1 the kinetics model for the chemical species and mechanical properties evolution, in Section 2.2 the model of the light (triggering the polymerization) diffusion in a continuum and in Section 2.3 the micromechanical model used to describe the mechanical behavior of photopolymerized materials. Then, in Section 3-based on the FEM modeling of the involved processes-we show and discuss the results related to the sensitivity analysis conducted, useful to understand the relationship between the photopolymerization features and the mechanics. Finally, Section 4 closes the paper with some concluding remarks and future perspectives.

| MULTIPHYSICS MODEL OF PHOTOPOLYMERIZATION
Photopolymerization is a chemical-physical process requiring a multidisciplinary approach. The process takes place according to the following phenomena: (1) the light propagation process in the medium being solidified, (2) the kinetics of the chemical species, and (3) the evolution of the mechanical properties. The light intensity distribution within the domain (induced by the moving laser source, typically used in Stereolithography and in other AM technologies based on the same principle) represents the stimulus responsible for the polymerization process to occur, consisting in initiation, propagation and termination of polymer chains, see Figure 1.

| Kinetics of the chemical species modeling and mechanical properties evolution
The light-induced polymerization of a liquid monomer resin (photoplymerization) can be chemically described by a set of chemical reactions (see Figure 1) describing the evolution of the chemical species involved, namely, (1) the monomer molecules (M), (2) the photo-initiators molecules (Ph I ), and (3) the free radicals (R), whose concentrations in the initial reference configuration at a generic time are C M ðX, t), C I ðX, tÞ and C R ðX, tÞ, respectively. The evolution of these species is responsible for the initiation, propagation and termination of polymer chains in the continuum and consequently for the solidification process (also known as curing): the light irradiation induces the conversion (and therefore the consumption) of Ph I into R; the last ones, reacting with M, induce the activation of functional groups P that grow within the material by reacting with other monomers M and forming the polymer chains. The chains grow within the medium until another chain or a free radical is encountered (end of propagation, P end ). The above-illustrated phenomena are described by a set of partial differential equations 8-10 : together with the related initial conditions: C I ðX, t ¼ 0Þ ¼ C I0 , C R ðX, t ¼ 0Þ ¼ C R0 and C M ðX,t ¼ 0Þ ¼ C M0 : In Equations 1-3, the superposed dot indicates the time derivative of the variable in turn while k pr , k p , and k t are the photodecomposition, the propagation and the termination rate constants (hereafter assumed do be independent by the curing evolution), respectively; finally, m is the number of radicals generated by the decomposition of each photo-initiator. The solution of this system of differential equations requires the knowledge of the light intensity I at each point of the domain (see Section 2.2). Finally, the so-called degree of cure (ϱ), a fundamental quantity quantifying the degree of polymer chains development in the material, can be obtained. This parameter affects the polymer chains concentration (crosslinked chains) c a , thus representing an internal variable related to the stiffness of the material; it can be simply evaluated as ϱðX, tÞ ¼ 1 À C M ðX, tÞ=C M0 ðXÞ. 8,11,12 The polymer chain concentration can be finally evaluated as follows: where H is the Heaviside function defined as H ¼ 0 for ϱðX, tÞ < ϱ gel or H ¼ 1 for ϱðX, tÞ ≥ ϱ gel . ϱ gel represents the degree of conversion of the material at the gelation point, w is a model parameter regulating the rate of the material stiffness increase after the gelation point, while μ is the shear modulus of the fully-cured polymer (i.e., the maximum shear modulus achievable with the specific resin being printed).

| Model of the light diffusion in a continuum
When the boundary of the domain is irradiated by a laser source, the field of the light intensity inside the material can be determined; starting from the lighted surface, its intensity decreases going down in the domain because of the chemical properties of the material (material influence) and because of the distance from the laser source (geometrical F I G U R E 1 Schematic view of photopolymerization: liquid resin doped with photo-initiator molecules at the initial stage before light irradiation (A) and polymer chains evolution due to light spread within the domain (B); the chemical reactions of the photopolymerization process are also reported influence). The reduction of the light intensity within a semi-transparent continuum, whose domain is here denoted with Ω, with respect to the light intensity irradiated on the boundary ∂Ω, depends on the so-called material absorbance A(X,t), being X Ω the position vector of a point within the domain and t the time.
The light irradiating the boundary of a semi-transparent continuum induces an instantaneous light intensity field within the continuum; this field can be evaluated by solving the well-known Beer-Lambert law, 13 whose expression in a 3D domain is: lðX, tÞ Á r X IðX, tÞ ¼ ÀAðX, tÞIðX, tÞ for X Ω , where 5 is the gradient operator while l is the versor defining the orientation of the incoming laser beam; the Dirichlet boundary condition, namely I X,t ð Þ¼ I X,t ð Þ¼ for X Ω 0 (where I m is the peak intensity of the laser light, v is the laser speed and c is a parameter regulating the width of the light irradiated zone) are required to solve the problem, see Figure 1B.
It is worth mentioning that the Beer-Lambert law in 1D problems with a constant material's absorbance can be solved analytically; in a general 2D or 3D setting, even if the absorbance is assumed to be constant, the closed-form solution of the Beer-Lambert equation becomes almost impossible because of the generic shape of the domain and/or of the boundary conditions. However, the use of standard computational techniques, such as the finite element (FE) one, allows efficiently solving the multi-physics problem even for a generic problem's domain with the related boundary conditions.
In our approach we have reformulated the light diffusion problem by using a weak form, suitable to be readily implemented in a FE framework. After discretizing the domain with FE, the above equations are expressed by the following system of algebraic equations: where Q is the vector of the FE nodal values of the incoming light intensity (assigned on the top boundary of the domain), I is the vector of the unknown light intensities to be evaluated (i.e., the vector containing the light intensities at each node of the discretized domain), while the matrices E , A and D represent the light gradient matrix, the depletion matrix and the stabilization matrix of the overall discretized problem, respectively. The above matrices can be evaluated by assembling the corresponding matrices related to each FE, that is, and [B] the standard shape functions and compatibility matrix, respectively, while Q > 1 and h e represent a coefficient and a characteristic size of the FE, respectively. It is worth mentioning that the material absorbance used to define the depletion matrix is related to the physical-chemical evolution of the material being solidified. It is usually expressed as AðX, tÞ ¼ θC I ðX, tÞ þ A pol ϱðX,tÞ þ A mon ð1Þ À ϱðX, tÞÞ, where θ is the photo-initiator molar absorbance while A pol and A mon are the absorbance of the polymer and of the monomer molecules, respectively. 8 In our model, we consider that the photopolymerization process does not induce any temperature increase of the material, since all the energy provided by the light source is assumed to be used for inducing the polymerization reactions and does not heat the material. Moreover, the conversion of photoinitiators into free radicals and the reaction of the latter with monomer units activating the functional groups responsible for the chain to grow, is assumed to be an isothermal process. This hypothesis allows neglecting the thermal problem in the photopolymerization phenomenon; it is worth mentioning that the temperature increase during the photopolymerization could affect the kinetics rate constants (see Equations 1-3) and in turn the rate of monomer conversion. It has been shown 14 that if the thermal diffusion takes place faster than the rate of polymerization, the system evolves as if it were isothermal; on the other hand, when thermal diffusion is slower, an accumulation of heat behind the polymerization front takes place, resulting in an increase of the degree of conversion corresponding to an acceleration of the growth of the solid.

| Mechanics of polymer's network
The molecular structure of a polymer suggests describing its mechanical behavior on the basis of the entropic energy, usually adopted to characterize disordered systems. Within this framework, a well-established approach is commonly used: the state of a polymer chain (whose conformation is often described by the so-called freely jointed chain model, FJC) is assumed to depend only on its end-to-end vector r identified by the extremities of the chain itself. In the reference configuration Ω 0 and in the corresponding deformed one Ω such a vector is identified by r 0 and r ¼ Fr 0 , respectively, where F represents the deformation gradient tensor. 15 The above expression for r is based on the so-called affine deformation hypothesis, usually adopted in the mechanics of polymers, that allows us to identify the microscopic deformation with the macroscopic one. 16 The mechanical energy stored in a single chain ψ is usually assumed to depend on the end-to-end vector; the force existing in the chain is thus given by the following derivative 16 : We assume the time t of the mechanical problem to be measured from the end of the polymerization process, that is, t ≥ t c . By considering a single point of a polymeric material, once the distribution of the end-to-end vectors is known at that point its state is completely defined; to this end, let us introduce the distribution function of the end-to-end vectors ρ(r) providing the number of chains, whose end-to-end vector falls within the sphere centered at r and with radius dr j j, per unit volume of material. 17 From a mechanical perspective, such a function accounts only for the mechanically active chains, that is, those connected at both ends to the network; the amount of active chains per unit volume, or chain concentration, is hereafter indicated with c a , so the following equality holds: where the integral is evaluated over the chain configuration space. 18 The energy per unit volume is given by The derivative of the energy function Equation 8 with respect to the deformation gradient allows determining the nominal stress tensor: P , p being the first Piola stress tensor and the hydrostatic pressure enforcing the isochoric deformation, respectively, while J = detF is the volumetric deformation. 19 Correspondingly, the Cauchy stress is given by where it is J = 1 because of the assumed isochoric deformation. From the above discussion, it appears that the stress state of the network is known once the current distribution function ρ(r,t) and the corresponding one in the stress-free reference state ρ 0 are known. In a deformed material, the current ρ(r,t) can be obtained starting from ρðr, t ¼ 0Þ ¼ ρ 0 ðrÞ as ρ r, t ð Þ¼ρ 0 r ð Þþ Ð t 0 _ ρ r, τ ð Þdτ; the time rate _ ρ r, t ð Þ required in the previous integral is given by 18 : where L is the velocity gradient tensor, L ¼ _ FF À1 . The mechanical behavior of a polymer network is thus conveniently describable by the above illustrated approach; in a photopolymerized polymer the chain density c a is evaluated by solving the kinetic Equations 1-4 and by adopting the standard Gaussian distribution for ρ 0 (r). 18 The present multi-physics approach has been implemented in ad-hoc in house FE code where the governing equations (and the related boundary conditions) of the light diffusion, chemical kinetics and mechanics have been coded in nonlinear solution algorithms.

| SENSITIVITY ANALYSIS
A sensitivity analysis is aimed at assessing the variation of the response of an engineering system because of small variations of the input quantities, and is typically employed to quantify the effect of the uncertainty of the input data. In this work, the sensitivity analysis is performed to quantify how the uncertainties of the speed and light peak intensity of the laser beam influence the mechanical properties of the final printed component.
The engineering system here considered is represented by a rectangular domain of liquid monomer resin undergoing photopolymerization; the process leads to a simple structural element, namely a single-printed layer beam with rectangular cross-section.
We simulate the curing process of a liquid monomer resin contained in a 2D vat having depth s = 10 mm and width w = 30 mm (single-layer). The initial properties of the system used in Equations 1-3 are assumed to be C M0 = 3000 mol/ m 3 , C I0 = 20 mol/m 3 , C R0 = 0 mol/m 3 . The chemical rate constants are assumed as follows: k pr = 8 Á 10 À4 s 2 /kg, k p = k t = 0.21 m 3 Á s Á mol À1 . The shear modulus of the fully-cured polymer is assumed to be μ ¼ 267 MPa while we adopt w ¼ 3. The sensitivity analysis has been performed by simulating several curing processes by assuming both a constant material absorbance AðX,tÞ ¼ A L ¼ 600 m À1 (representative of the liquid resin), and a variable absorbance evolving during the curing process as shown in Section 2.2, by assuming θ ¼ 30m 2 =mol, A pol ¼ 2400 m À1 and A mon ¼ 0 m À1 . In the case of a variable absorbance evolving during the curing process, when the material attains the solid state all the available photo-intiators are consumed, i.e. C I ðX, tÞ ! 0, and we get AðX, tÞ ¼ A S ¼ 2400 m À1 . After simulating the printing process, the mechanical response of a printed cantilever beam under a transversal load is considered.
The sensitivity analysis outlines the role played by the input uncertainties on the mechanical properties of the printed material (see Section 2.1) and on its mechanical response evaluated by using the micromechanical model of Section 2.3.
The dimensionless shear modulus plotted along a vertical section placed in the middle of the printed beam is displayed in Figure 2: in order to quantify the effect of uncertainties of the laser intensity, the reference laser velocity has been assumed equal to v = 1 mm/s, while the reference maximum laser intensity of I m = 100 W/m 2 has been varied by F I G U R E 2 Dimensionless shear modulus along a vertical line placed in the middle of beam element, whose photopolymerization process has been simulated by (i) varying the light intensity while keeping a constant laser speed of v = 1 mm/s (A, B), and (ii) varying the laser speed for a fixed laser light intensity of I m = 100 W/m 2 (C, D), by assuming a constant (A, C) and a variable (B, D) material absorbance, being A S = 2400 mm À1 . Dimensionless depth η = d/s (±5 % , ± 10 % , ± 50%); a constant and a variable material absorbance has been adopted, Figure 2A,B. Similarly, in Figure 2C,D the effect of the uncertainties of the laser speed (variations of ±5 % , ± 10 % , ± 50% with respect to the reference value v = 1 mm/s have been adopted), for a fixed maximum laser light intensity I m = 100 W/m 2 is displayed. In all the investigated cases, small variations of the input quantities (±5% and ±10%) do not provide noticeable effects on the mechanical quality of the material, quantified through its dimensionless shear modulus. On the other hand, greater variations of the input quantities (±50%) provide a more evident changes of the mechanical property of the material. Moreover, as can be noticed from Figure 2, the final material seems to be more sensitive to input variations when a constant material absorbance is assumed, see Figures 2A,B and 2C,D. Figure 3 displays the load-displacement curves of the beam element considering the same input variations described above; a slight variation of the photopolymerization parameters (±5% and ±10%) is responsible for a limited change of the load-displacement response, while greater variations (±50%) induces a more severe variation of its response. Interestingly, the variation of the maximum light intensity provides load-displacement curves close to the original one, either assuming a constant or a variable material absorbance, Figure 3A,B. On the other hand, the variation of the laser speed strongly affects the load-displacement curves, see Figure 3C,D. According to the above-discussed results, the investigated photopolymerized structural element results to be more sensitive to the variation of the laser speed rather than to that of the maximum light intensity.
Moreover, it can be noticed that the stiffness of the structural element (as well as its deformed shape at a certain load level, see Figure 3C, points P 1 and P 2 and the corresponding deformed shapes plotted in the same figure) can be tuned by properly setting the printing parameters. From a functional material perspective, starting from a given raw material, this offers the possibility to tailor the mechanical behavior of the printed element by acting on the values of the printing parameter.
In Figure 4 the dimensionless displacement attained for a load level of P = 20 N (normalized with respect to the beam made of the fully cured material, that is, the one cured for an infinite time) for a fixed laser speed and by varying the maximum light intensity ( Figure 4A) or by fixing the maximum light intensity and varying the laser speed ( Figure 4B), is shown for both a constant and a variable material absorbance. As already mentioned, the response of F I G U R E 3 Load-displacement curves of the beam element obtained according to the corresponding photopolymerization parameters reported in Figure 2 the system is more sensible to the laser speed variations. On the other hand, irrespectively of the input parameters, by considering the printing process to be characterized by the lowest value of the absorbance (that of the liquid monomer, black curves in Figure 4A,B) the model overestimates the stiffness of the element compared with the case of a variable absorbance (red curves in Figure 4A,B). The simulation of the photopolymerization process by adopting a variable absorbance results to be more realistic, since the shielding effect produced by the solidified part of the resin in limiting the penetration of light in the material, so hindering its curing, is accounted for.

| CONCLUSIONS
A multi-physics model of the photopolymerization used in AM has been illustrated. The governing equations of the polymerization process induced by light and those related to the mechanics of the polymer network have been illustrated. The developed approach, implemented in a FE framework, has been used to study the sensitivity of the mechanical properties and of the corresponding mechanical response of the photopolymerized material to the uncertainty of the printing parameters. The numerical analyses revealed that the properties of the final material are slightly affected by the variation of the laser speed and the maximum laser light intensity when their variations are limited (±5% and ±10%). On the other hand, variations of the order of about ≥50% are responsible for a significant variation of the material's properties and of its mechanical response. Further, the effect of the material absorbance has been considered: by assuming the photopolymerization taking place at a constant value of the absorbance (that of the initial liquid monomer) overestimates the stiffness of the material, while if the evolution of the absorbance during the curing process is accounted for, more realistic results are obtained. The proposed multi-physics approach is readily extendable to the simulation of the photopolymerization taking place in multiple superposed layers as occurs in real AM printing; the forthcoming development of the proposed approach will consider such an aspect by accounting for the effective aborbance of the already solidified layers when solving the light diffusion problem.
The developed approach allows investigating and optimizing the printing strategy of additively manufactured components. Moreover, the possibility to precisely design the mechanical properties distribution, opens new paths in the optimal material design aimed at tuning the structural components' behavior, enabling the production of elements with tailored and programmable responses.

DATA AVAILABILITY STATEMENT
The data that support the findings of this study are available from the corresponding author upon reasonable request.