Keywords: surrogate modeling, injection molding, crystallization, simulation, thermoplastic, model parameter, metamodel. v modélisation de substitution, moulage par injection, cristallisation, simulation, thermoplastique, paramètre de modèle, métamodèle

The computational costs of injection molding simulations have been increasing in the past years due to the higher complexity of the embedded models. This is especially problematic in case of using such simulation models for optimization routines or sensitivity analyses. One way to overcome this challenge is by having a surrogate model, also known as a metamodel, of these high-fidelity simulations, which provides a cheaper way to perform these types of analyses. These surrogate models can play an important role in the case of the injection molding of semi-crystalline polymers to model the flow-induced crystallization process. To date, most commercial software do not explicitly take polymer crystallization into account leading to various errors in the fill predictions as well as the calculation of warpage and shrinkage. This is mainly due to the common complexity of the models used to describe crystallization and the challenging respective model parameter identification process under injection molding conditions. To close this gap, in this thesis, the feasibility of using surrogate modeling to identify modeling parameters is first studied. This is then followed by the implementation of a thermomechanical crystallization model in order to describe the flow-induced and quiescent crystallization of an unreinforced semi-crystalline thermoplastic material during injection molding. The crystallization model is defined alongside crystallization-dependent viscosity, PVT and solidification models in the commercial software Autodesk R ⃝ Moldflow R ⃝ Insight 2021 using the Solver API feature. The model parameters are identified using a calibration scheme that employs three surrogate models representing the simulated pressure results to perform a multi-objective optimization. The fill predictions as well as the calculated pressure fields are presented using the calibrated model parameters in comparison to those measured during the actual injection molding of a polyoxymethylene part with different process conditions. The results show major improvements in the predictions of the pressure signals as well as the filling status of the produced parts and the estimated skin layer thicknesses formed under high-shear conditions. Additionally, the calibrated models are tested using various mold geometries to assess the calibrated models' performance.

Introduction

List of Tables

Introduction

Injection molding is one of the main polymer processing techniques with which more than onethird of thermoplastic materials are produced [START_REF] Osswald | Material science of polymers for engineers[END_REF]. The injection molding process is cyclic and ideal for the mass-production of complex geometries with tight tolerances [2]. In the case of part defects during production, a change in process conditions may not be sufficient and mold modifications must be applied to solve the problem [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. However, this is expensive and time-consuming. Thus, it is important to be able to foresee such issues in the design phase of the part using simulation techniques.

Therefore, to minimize such costs while guaranteeing a superior product quality, injection molding simulation is widely used to efficiently design molds and provide a tool to determine optimal process settings that overcome commonly observed defects such as warpage, shrinkage and short shots.

However, over the last years, the growing interest in highly accurate simulations resulted in the increase of the complexity of the embedded models and in turn to higher computational costs. This problem is specifically underlined when using such simulation models for optimization, sensitivity analyses or uncertainty quantification. Surrogate modeling offers a way to overcome the limitation of these high-fidelity simulations and to perform such analyses in a more cost-and time-efficient way. These models are usually generated using a relatively low number of simulations obtained by varying a specific number of input parameters and are then able to approximate the output of interest in a matter of seconds instead of minutes or hours using the high-fidelity simulation. An example application is presented in Figure 1 where a surrogate model is used to calibrate simulation model parameters by performing an efficient optimization routine. In the presented case, the surrogate model represents one output result learned from running multiple high-fidelity simulation runs and is then used independently to obtain the optimal input parameters that best fit the experimental results.

Currently, in the field of injection molding simulation, the use of surrogate modeling in the literature is mainly limited to the optimization of process parameters such as the mold and melt temperatures INTRODUCTION Figure 1 -Model parameter optimization scheme using a surrogate model [4]. [5,[START_REF] Chen | Multi-objective optimal approach for injection molding based on surrogate model and particle swarm optimization algorithm[END_REF], packing and cooling times [START_REF] Li | A modified global optimization method based on surrogate model and its application in packing profile optimization of injection molding process[END_REF][START_REF] Cheng | Multiobjective optimization of injection molding parameters based on soft computing and variable complexity method[END_REF], packing pressure [START_REF] Gao | Intelligent methods for the process parameter determination of plastic injection molding[END_REF][START_REF] Shi | A warpage optimization method for injection molding using artificial neural network with parametric sampling evaluation strategy[END_REF]... One recent publication by Ivan et al. [START_REF] Ivan | Improving Numerical Modeling Accuracy for Fiber Orientation and Mechanical Properties of Injection Molded Glass Fiber Reinforced Thermoplastics[END_REF] stands out in which the authors use a surrogate model to identify two model parameters that describe the fiber orientation. A similar approach could be utilized to calibrate other types of complex material models such as a crystallization model. Due to the complexity of experimentally identifying crystallization parameters under injection molding conditions, many commercial software do not take it into account even when dealing with the injection of semi-crystalline materials. This leads to inaccuracies in the fill predictions as well as in the estimation of warpage and shrinkage in the simulation [START_REF] Kennedy | Flow analysis of injection molds[END_REF].

In the last decades, great efforts have been made in order to model polymer crystallization while taking into account the flow history of the melt. A common approach utilized to describe flow-induced crystallization, observed during injection molding, is based on the Avrami-Kolmogorov nucleation and growth theory [START_REF] Avrami | Kinetics of Phase Change. I General Theory[END_REF][START_REF] Kolmogotov | A Statistical Theory for the Recrystallisation of Metals[END_REF] combined with Schneider's rate equations [START_REF] Schneider | Non -Isothermal Crystallization Crystallization of Polymers: System of Rate Equations[END_REF]. Such method necessitates the experimental determination of a high number of model parameters and solving of an intricate and interdependent system of differential equations [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF]. Whereas, Poitou et al. [START_REF] Poitou | Cristallisation induite par écoulement ou par déformation d'un polymère -une approche thermodynamique[END_REF][START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF] proposed a thermo-mechanical approach describing both the flow-induced and quiescent crystallization with a much smaller number of model parameters making it more attractive to use in simulation software.

The main aim of this thesis is the use of surrogate modeling in the model parameter calibration The models are implemented in 3D and calibrated using three surrogate models representing the simulated pressure at different processing conditions. The ability of the calibrated models to predict accurately the pressure, fill status and skin layer thickness is tested at different conditions using a base mold geometry and three variations of it.

INTRODUCTION

Chapter 1

Theory and State of the Art

INJECTION MOLDING PROCESS

Injection Molding Process

The injection molding process is one of the most widely used plastics manufacturing techniques in the industry such as more than one-third of all thermoplastic materials are injection molded [START_REF] Osswald | Material science of polymers for engineers[END_REF][START_REF] Zheng | Injection Molding[END_REF].

It enables the high production volumes of complex plastic parts while having short cycle times [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. A cycle constitutes mainly of six phases [START_REF] Osswald | Material science of polymers for engineers[END_REF][START_REF] Kennedy | Flow analysis of injection molds[END_REF][START_REF] Zheng | Injection Molding[END_REF]:

1. Mold closing: the mold is closed to obtain the cavity into which the molten material is injected.

2. Injection/Filling: the screw moves forward and forces the melt into the mold cavity.

3. Packing/Holding: the screw is held in its forward position after the mold is filled in order to maintain a holding pressure that compensates for material shrinkage.

4. Cooling: the gate freezes and the cavity is isolated from the pressure applied by the melt, the part cools down and solidifies until ejected.

5. Plastication and screw back: the screw rotates back while plasticating the feed for the next shot which is now in a molten state in front of the screw. 6. Ejection: the mold opens when the part is sufficiently cooled and it is then ejected.

The total cycle time t cycle of an injection molding process is given by [START_REF] Osswald | Material science of polymers for engineers[END_REF]:

t cycle = t closing + t cooling + t ejection (1.1)
where the closing time t closing and ejection times t ejection usually last less than a few seconds whereas the cooling time t cooling dominates the process as the polymer melt is being cooled as soon as it enters the cavity after the closing of the mold till it is ejected from it. Therefore, t cooling includes the cooling taking place during the filling and packing phases as the mold is closed and surrounds the polymer melt during these stages.

During the operation of an injection molding machine, plastic granules are fed through the hopper shown in the simplified sketch in Figure 1.1. The fed material is then transported forward by the rotating screw. The granules melt while being transported due to both the friction created by the screw and the conduction from heating units along the barrel. The molten material is injected into the mold by the screw typically with high flow rates and therefore high shear rates. The main heat transfer mechanism during the filling stage is through convection of the melt in addition to some viscous heating primarily in the runner system and the gates [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. Heat is constantly removed by conduction through the mold wall to the cooling system. This leads to the formation of a frozen layer or skin layer near the wall. During the packing stage, the flow rates are relatively low and conduction becomes the main heat transfer mechanism leading to the increase of this frozen layer until the gate freezes followed by the solidification of the part due to conduction until it is ejected from the mold.

Injection Molding Simulation

To ensure high product quality and minimize design/production costs, injection molding simulation is broadly used in order to efficiently design molds and enable the identification of optimal process settings that mitigate common defects such as warpage, shrinkage, weldline and short shots. The simulation of the injection molding process constitutes using a numerical method with well-defined material models to solve a set of conservation equations [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. Many commercial and academic software have been developed to tackle this problem which involves multiple heating mechanisms as described in Section 1.1 in addition to a phase change and time-dependent boundary conditions at the frozen layer during the various phases of the process. Taking into account all these complex physical phenomena in the simulation is not straightforward and necessitates some simplifying assumptions and approximations.

Main Approximations

These approximations are mainly needed due to the difficulties in determining some material properties experimentally under similar conditions as the ones experienced during injection molding.

INJECTION MOLDING SIMULATION

An overview of these assumptions is summarized in the following points with a focus on the injection molding of an unfilled semi-crystalline thermoplastic in the commercial software Moldflow R ⃝ during the filling, packing and cooling analyses [START_REF] Kennedy | Flow analysis of injection molds[END_REF]:

• Heat sources: the latent heat generated during crystallization is ignored and no heat source term is considered in the energy equation.

• Specific heat capacity: not coupled with crystallization and determined using differential scanning calorimetry for much lower cooling rates than those present during processing.

• Thermal conductivity: defined as a function of temperature without taking into account any pressure dependency.

• Solidification criterion: the use of a no-flow temperature (transition temperature) below which the material exhibits a high viscosity and solidifies with no link between solidification and crystallization.

• Pressure-Volume-Temperature (PVT) data: generated without taking into account the cooling rate or the deformation history of the sample which affects the transition temperature.

Governing Equations

By taking into account the simplifications listed in Section 1.2.1, the governing equations used to simulate the flow in the mold cavity are the three conservation equations. First, the conservation of mass for a fluid is represented as:

∂ρ ∂t + ∇ • ρv = 0 (1.2)
where ρ is the polymer density, t is the time and v is the velocity vector. Second, the conservation of momentum is given by:

ρ ∂v ∂t + ρv • ∇v = -∇P + ∇ • τ + ρg (1.3)
where P is the pressure, τ is the viscous stress tensor and g is the gravitational acceleration vector.

And, third, the conservation of energy is described by:

ρc p ( ∂T ∂t + v • ∇T ) = ∇ • (k∇T ) + τ : ∇v + ζT ( ∂P ∂t + v • ∇P ) (1.4)
where T is the temperature, k is the polymer thermal conductivity, c p is the specific heat capacity of the melt and ζ = -1 ρ ∂ρ ∂T is the polymer expansion.

POLYMER CRYSTALLIZATION

In the case of a 3D analysis, the above-presented equations are typically solved using the finite element method (FEM) with the CAD of the injection molded part being meshed in tetrahedral or hexahedral elements. For a more detailed overview of the injection molding process along with the equations needed to define its simulation refer to Osswald and Menges [START_REF] Osswald | Material science of polymers for engineers[END_REF] as well as Kennedy and

Zheng [START_REF] Kennedy | Flow analysis of injection molds[END_REF].

Polymer Crystallization 1.3.1 Semi-crystalline Thermoplastics

Thermoplastic polymers are divided into amorphous and semi-crystalline polymers. They can mainly be distinguished by the molecular arrangement of their chains as the polymer melt cools down into the solid state. For amorphous polymers, the molecules are randomly oriented and intertwined after cooling. Whereas for semi-crystalline polymers, the molecular structures form ordered regions called crystallites that are surrounded by amorphous regions [START_REF] Osswald | Material science of polymers for engineers[END_REF][START_REF] Kennedy | Flow analysis of injection molds[END_REF]. The coexistence of the two phases provides the polymer with a good balance of stiffness and toughness [START_REF] Michler | Electron microscopy of polymers[END_REF]. The properties of semicrystalline thermoplastics are thus dependent on the content and orientation of both the crystalline and amorphous phases [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. Semi-crystalline polymers are widely used in the plastic industry as both commodity plastics (e.g. PE, PP) as well as technical and functional polymers (e.g. PEEK, PVDF, POM) [START_REF] Michler | Electron microscopy of polymers[END_REF][START_REF] Lamberti | Flow induced crystallisation of polymers[END_REF]. This is mainly due to their versatile nature; where changes in the processing conditions, and thus the thermo-mechanical history, can lead to modifications in the final product's properties e.g. strength, hardness, permeability [START_REF] Santis | Effect of shear flow on spherulitic growth and nucleation rates of polypropylene[END_REF].

The schematic in Figure 1.2 presents the hierarchical arrangement in a semi-crystalline material and the different crystal morphologies observed using polarized optical microscopy in injection molded components. The randomly oriented polymer chain folds into crystal lamella which in turn stacks together to form a lamellar stacking morphology leading to spherulites in the case of quiescent crystallization under no-flow conditions. Whereas in the case of high shear flows, fully extended chains result in extended chain crystals forming parallel to the flow direction along with lamellar crystals forming perpendicular to the flow direction. This leads to the shish-kebab structures characteristic of flow-induced crystallization (FIC) [START_REF] Zuidema | Flow induced crystallization of polymers: application to injection moulding[END_REF].

Crystallization from the Melt

There exist multiple types of crystallization such as crystallization from solution, crystallization by stretching or crystallization from the melt. In this work, the main interest is studying polymer crystallization during the injection molding process and therefore only the crystallization from the melt is considered. This process involves two stages:

1. Nucleation: the formation of active nuclei in the liquid phase acting as starting points for the appearance of crystals [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. Two types of nucleation can be distinguished [START_REF] Hu | Polymer Physics[END_REF][START_REF] Liang | Theoretical Aspects of Polymer Crystallization in Multiphase Systems[END_REF]:

• Homogeneous nucleation: caused by heat motion and starts within a few polymer chains or segments (e.g. from the bulk polymer phase).

• Heterogeneous nucleation: appears on foreign substrates (e.g. nucleating agents, impurities, fillers) or interfaces in multiphase systems.

2. Growth: the growth of the formed nuclei into semi-crystalline morphological structures. These morphologies depend on the thermo-mechanical history experienced by the polymer melt leading to the formation of two distinct nuclei having different growing mechanisms [START_REF] Kennedy | Flow analysis of injection molds[END_REF][START_REF] Zuidema | Flow induced crystallization of polymers: application to injection moulding[END_REF]:

• Spherical nuclei: grow radially in space and form spherical structures known as spherulites, typically seen under no-flow conditions (quiescent crystallization).

POLYMER CRYSTALLIZATION

• Thread-like nuclei: grow perpendicular to the thread and form the shish-kebab structures, commonly seen under flow conditions with the presence of high strains or shears (FIC crystallization).

The presence of flow, as commonly seen during polymer processing, affects the crystallization process as follows [START_REF] Boutaous | Thermally and flow induced crystallization of polymers at low shear rate[END_REF]:

• increase of the nucleation density,

• acceleration of the crystallization kinetics,

• changes in semi-crystalline morphological structures,

• increase of the crystallization and melting temperature.

A lot of work has been done in the last few decades to better understand the physics behind polymer crystallization whether under quiescent or flow conditions. An extensive review of the various theories postulated is given by Zhang et al. [START_REF] Zhang | A Review on Polymer Crystallization Theories[END_REF]. The current understanding of this process under flow conditions is that there exists one threshold below which no changes in crystal morphology are observed and above which the nucleation density and growth rate are altered by the presence of flow in addition to a second threshold where thread-like (fibrillar) morphologies start appearing [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF][START_REF] Pantani | Analysis of flow induced crystallization through molecular stretch[END_REF]. In the following sections, a brief overview is given of the different experimental and modeling techniques used to study both types of crystallization.

Experimental Methods

In order to obtain a better understanding of the crystallization process of semi-crystalline polymers, various experimental techniques are utilized in the literature. Microscopy is one of the main methods used to gather morphological information on developed crystalline structures in a given sample. Polarized-light optical microscopy is most commonly used for these purposes [START_REF] Coccorullo | Crystallization kinetics and solidified structure in iPP under high cooling rates[END_REF][START_REF] Zuidema | Development and validation of a recoverable strain-based model for flow-induced crystallization of polymers[END_REF][START_REF] Zheng | Modeling of flow-induced crystallization of colored polypropylene in injection molding[END_REF][START_REF] Pantani | Thirty years of modeling of injection molding. a brief review of the contribution of UNISA code to the field[END_REF][START_REF] Troisi | Application of a multi-phase multi-morphology crystallization model to isotactic polypropylenes with different molecular weight distributions[END_REF] in addition to scanning electron microscopy [START_REF] Coccorullo | Crystallization kinetics and solidified structure in iPP under high cooling rates[END_REF][START_REF] Parenteau | Structure, mechanical properties and modelling of polypropylene for different degrees of crystallinity[END_REF] and atomic force microscopy [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF][START_REF] Pantani | Effect of flow-induced crystallization on the distribution of spherulite dimensions along cross section of injection molded parts[END_REF]. Light diffraction experiments are utilized to gain crucial knowledge about the crystallinity and the polymer chains' orientation.

The main analysis employed is the wide-angle X-ray diffraction (WAXD) [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF][START_REF] Troisi | Application of a multi-phase multi-morphology crystallization model to isotactic polypropylenes with different molecular weight distributions[END_REF][START_REF] Roozemond | Modeling flow-induced crystallization in isotactic polypropylene at high shear rates[END_REF]. Additionally, IR spectroscopy is occasionally applied for the same purposes but is not able to discern between different crystalline phases such as WAXD [START_REF] Lamberti | Flow induced crystallisation of polymers[END_REF].

POLYMER CRYSTALLIZATION

Differential scanning calorimetry (DSC) measurements are widely performed in the literature to provide quantitative information about the crystallization kinetics which are then used to develop models to describe the process. In addition to DSC, a shearing hot-stage Linkam apparatus is used in combination with optical microscopy to quantify the nucleation and growth steps of crystallization [START_REF] Pantani | Analysis of flow induced crystallization through molecular stretch[END_REF][START_REF] Zheng | Modeling of flow-induced crystallization of colored polypropylene in injection molding[END_REF][START_REF] Koscher | Influence of shear on polypropylene crystallization: morphology development and kinetics[END_REF][START_REF] Mu | Numerical investigation of the crystallization and orientation behavior in polymer processing with a two-phase model[END_REF]. One major shortcoming of the before-mentioned methods is that they mainly characterize and quantify the crystallization kinetics under quiescent conditions and are not yet fully capable of recreating process-relevant conditions in order to fully model flow-induced crystallization.

Nevertheless, there have been some advances in that regard such as the development of fast scanning calorimetry [START_REF] Toda | Insights into polymer crystallization and melting from fast scanning chip calorimetry[END_REF] and the use of an extrusion die in combination with optical microscopy to calculate the critical shear rate of activation [START_REF] Isayev | Injection molding of semicrystalline polymers. II. Modeling and experiments[END_REF][START_REF] Guo | Computer simulation of stress-induced crystallization in injection molded thermoplastics[END_REF].

Modeling Approaches

To study the crystallization kinetics, whether quiescent or flow-induced, various theories and models have been developed in the literature. These models mainly describe the evolution of the relative crystallinity α(t) as a function of time t and temperature T as a state model expressing α as a function of α and T. By definition, the relative crystallinity is the ratio of the crystallized volume X c to the ultimate crystallizable volume X ∞ , such as [START_REF] Boutaous | Thermally and flow induced crystallization of polymers at low shear rate[END_REF]:

α = X c X ∞ .
(1.5)

Quiescent Crystallization

The first modeling approach used to describe the kinetics of quiescent crystallization is based on the Kolmogorov-Avrami-Evans (KAE) theory [START_REF] Avrami | Kinetics of Phase Change. I General Theory[END_REF][START_REF] Kolmogotov | A Statistical Theory for the Recrystallisation of Metals[END_REF][START_REF] Evans | The laws of expanding circles and spheres in relation to the lateral growth of surface films and the grain-size of metals[END_REF], developed for isothermal conditions.

Kolmogorov [START_REF] Kolmogotov | A Statistical Theory for the Recrystallisation of Metals[END_REF] originally used a probabilistic approach to describe the space filling by the formation and growth of spherulites. Independently, Avrami [START_REF] Avrami | Kinetics of Phase Change. I General Theory[END_REF] used a geometrical approach to express this growth while taking into account the impingement of the spherulites as they grow and block the growth of the neighboring spherulites. Later on, Evans [START_REF] Evans | The laws of expanding circles and spheres in relation to the lateral growth of surface films and the grain-size of metals[END_REF] also developed the same model by analogy to expanding circular waves in water. All three models proved to be identical and are referred to as 1.3. POLYMER CRYSTALLIZATION the KAE theory and is given by:

dα = (1 -α)dΦ (1.6a) α(t) = 1 -exp(-Φ(t)) (1.6b)
where Φ is a fictitious crystallinity obtained by assuming an unrestricted growth G(t) of the spherulites which are represented as spheres. This fictitious crystallinity is defined as follows for a time interval between s and t:

Φ = 4π 3 ∫ t 0 Ṅq (s) [∫ t s G(u)du ] 3 ds (1.7)
where Ṅq (t) is the rate of nuclei number per unit volume. Under quiescent conditions and in the case of instantaneous nucleation, the nuclei number density N q (t) = N 0 H(t) where N 0 is the constant nuclei density also known as the number of activated nuclei, H(t) is the Heaviside unit step function such as H(t) = 0 for t < 0 and H(t) = 1 for t ≥ 0. In turn, Ṅq (t) = N 0 δ(t) where δ(t) is the Dirac delta function centered at t = 0. By integrating equation 1.7, equation 1.6b is expressed for this considered case as:

α(t) = 1 -exp ( -4π 3 N 0 G 3 t 3 ) (1.8)
where the time exponent represents the Avrami index n which is 3 in the case of spherulitic growth and instantaneous nucleation as shown in equation 1.8. In general, the KAE model can be represented using the following expression:

α(t) = 1 -exp(-k Avrami t n ) (1.9)
where k Avrami is the Avrami kinetic constant.

According to experimental findings by Koscher and Fulchiron [START_REF] Koscher | Influence of shear on polypropylene crystallization: morphology development and kinetics[END_REF], for instantaneous nucleation, the nuclei number density increases when the degree of supercooling ∆T = T 0 m -T increases where T 0 m is the equilibrium melting temperature, such as:

ln N 0 (T ) = a N ∆T + b N (1.10)
where a N and b N are experimentally-fitted parameters. More complex temperature dependencies can be found in the literature given by Coppola et al. [START_REF] Coppola | Effects of the degree of undercooling on flow induced crystallization in polymer melts[END_REF].

As for the temperature dependency of the spherulite growth, the Hoffman-Lauritzen theory is widely used to describe it [START_REF] Lauritzen | Theory of formation of polymer crystals with folded chains in dilute solution[END_REF]:

G(T ) = G 0 exp ( - U ⋆ R(T -T ∞ ) ) exp ( - K g T (T 0 m -T ) ) (1.11) 1.3. POLYMER CRYSTALLIZATION
where G 0 and K g are experimentally-determined parameters, U ⋆ is the activation energy for segmental jump of polymer molecules with a universal value of 6270 J/mol, R is the gas constant and T ∞ is the temperature at which no molecular displacement happens defined using the glass transition temperature T g as T ∞ = T g -30.

Ozawa [START_REF] Ozawa | Kinetics of non-isothermal crystallization[END_REF] extended the phenomenological KAE model to adhere to the non-isothermal case such as:

α(t) = 1 -exp ( - k Ozawa Ṫ n ) (1.12)
where k Ozawa is the Ozawa kinetic constant and Ṫ is the cooling rate. However, this relation showed some limitations as it applies to a limited range of cooling rates.

Nakamura et al. [START_REF] Nakamura | Some aspects of nonisothermal crystallization of polymers. I. Relationship between crystallization temperature, crystallinity, and cooling conditions[END_REF] introduced an isokinetic approach to describe the kinetic constant in the KAE model under non-isothermal conditions:

α = 1 -exp [ - (∫ t 0 K(T )dt ) n ] (1.13)
where K(T ) is the Nakamura kinetics crystallization function and can be related to the Avrami isothermal kinetic constant by K(T ) = k Avrami (T ) 1/n . The differential form of the Nakamura model is commonly used in literature such as:

α(t) = K(T )n(1 -α(t))[-ln(1 -α) 1-1 n ]. (1.14)
An alternative approach for the modeling of non-isothermal crystallization is proposed by Schneider et al. [START_REF] Schneider | Non -Isothermal Crystallization Crystallization of Polymers: System of Rate Equations[END_REF]. This approach uses a set of four first-order differential equations derived from the KAE model to solve for the relative crystallinity. This is done by differentiating equation 1.7 four times with respect to time as follows:

Φ 0 = Φ = 4π 3 [∫ t s G(u)du ] 3 dN q (s), (1.15a
)

Φ 1 = 1 G Φ0 = 4π [∫ t s G(u)du ] 2 dN q (s), (1.15b)
Φ 2 = 1 G Φ1 = 8π [∫ t s G(u)du ] dN q (s), (1.15c)
Φ 3 = 1 G Φ2 = 8πN q (t). (1.15d) 1.3. POLYMER CRYSTALLIZATION
Therefore, the Schneider rate equations are defined as:

Φ3 = 8π Ṅq (t), (1.16a)
Φ2 = GΦ 3 , (1.16b) Φ1 = GΦ 2 , (1.16c) Φ0 = GΦ 1 , (1.16d)
The Schneider functions Φ i provide morphological information concerning the formed crystalline structures such as:

• Φ 0 : the total volume of the undisturbed spherulites (no impingement) per unit volume,

• Φ 1 : the total surface area of the undisturbed spherulites per unit volume,

• Φ 2 : 8π times the sum of the radii of the undisturbed spherulites per unit volume,

• Φ 3 : 8π times the number of the undisturbed spherulites per unit volume.

The solution of the system of equations 1.16 is Φ 0 which can be written in the form of equation 1.6b to take into account the impingement of the spherulites, such as:

α(t) = 1 -exp(-Φ 0 (t)) (1.17a) Φ 0 (t) = -ln(1 -α(t)).
(1.17b)

Flow-induced Crystallization

As already mentioned in Section 1.3.2, the presence of flow greatly affects the crystallization kinetics. This has been observed experimentally, however, the exact physics behind this phenomenon is still not fully understood [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF]. Therefore, literature works proposed different approaches to model flow-induced crystallization. These can be clustered in two main categories according to which type of model the work is based on:

• Nakamura's model: the increase of crystallization kinetics is taken into account by multiplying the kinetic function K(T ) in equation 1.14 by an enhancement factor depending on the postulated driving force for flow-induced crystallization (stress, strain, shear rate, melting temperature increase). This approach neglects the changes in morphology.

POLYMER CRYSTALLIZATION

• KAE model and/or Schneider's rate equations: the increase of crystallization kinetics is connected to enhancements in the nucleation rate and/or growth rate of the crystal morphologies as well as the volumetric free energy difference between the molten and crystalline phases used while describing the nucleation rate.

Table 1.1 presents an overview of the literature in the field of modeling the FIC categorized according to the two approaches presented previously in addition to the proposed driving force used to describe the effect of the flow on the crystallization kinetics. These two approaches are not strictly separate as some of the functions used to describe the driving force proposed by the authors using the Nakamura model can be easily used to obtain morphological information using a modified form of the Schneider rate equations. In the following, an example of that is shown, starting first by using an enhancement function derived by Eder et al. [START_REF] Eder | Crystallization processes in quiescent and moving polymer melts under heat transfer conditions[END_REF] followed by one proposed by Tanner et al. [START_REF] Tanner | A comparison of some models for describing polymer crystallization at low deformation rates[END_REF].

In the presence of flow, the total nuclei number density can be expressed as:

N (t) = N q (t) + N f (t) (1.18)
where N q is the number of existing nuclei at the beginning of the process equivalent to those present under quiescent conditions with a temperature dependence given in equation 1.10 and N f is the flowinduced nuclei number density. Eder et al. [START_REF] Eder | Crystallization processes in quiescent and moving polymer melts under heat transfer conditions[END_REF] proposed the following differential equation to represent the flow-induced nucleation:

Ṅf + 1 λ N N f = Υ (1.19)
where λ N is a temperature-dependent relaxation time and Υ is a function dependent on flow variables and temperature. This function represents the main driving force that enhances nucleation in the presence of flow. Eder et al. [START_REF] Eder | Crystallization processes in quiescent and moving polymer melts under heat transfer conditions[END_REF] postulated that the shear rate is that driving force and defined this function as:

Υ = g n ( γ γn ) 2 (1.20)
where γ is the shear rate, γn is the critical shear rate of activation and g n is a factor. By assuming that each shish-kebab is represented as a cylinder, it is possible to modify the Schneider rate equations to describe the evolution of these crystals in the presence of flow using the enhanced nucleation expression 1.3. POLYMER CRYSTALLIZATION Table 1.1 -Overview of the literature in the field of modeling the flow-induced crystallization during polymer processing, based on [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF].

Crystallization model based on

Effect of flow on crystallization kinetics using Author(s)

Nakamura's model [START_REF] Nakamura | Some aspects of nonisothermal crystallization of polymers. I. Relationship between crystallization temperature, crystallinity, and cooling conditions[END_REF] Multiplying factor function of the extra stress tensor Doufas et al. [START_REF] Doufas | Simulation of melt spinning including flow-induced crystallization: Part I. Model development and predictions[END_REF], Zinet [START_REF] Zinet | Numerical modeling of nonisothermal polymer crystallization kinetics: Flow and thermal effects[END_REF] Multiplying factor function of the strain and the shear rate

Φ3,f + Φ 3,f λ N = 8π ( γ γn ) 2 g n , (1.21a) Φ2,f + Φ 2,f λ l = Φ 3,f ( γ γn ) 2 g l , (1.21b) Φ1,f = GΦ 2,f , (1.21c) Φ0,f = GΦ 1,f , (1.21d)
where λ l is the temperature and shish length dependent relaxation time during axial growth and g l is a factor for the axial growth of the shish. As for the crystal growth G, Eder assumed that it still follows the Hoffman-Lauritzen theory presented in equation 1.11. Additionally, the flow-induced Schneider functions represent the following:

• Φ 0,f : the total volume of the undisturbed shish-kebabs (no impingement) per unit volume,

• Φ 1,f : the total surface area of the undisturbed shish-kebabs per unit volume,

• Φ 2,f : 8π times the total length of the shish per unit volume,

• Φ 3,f : 8π times the number of the undisturbed flow-induced nuclei N f per unit volume.

To take into account impingement and the fact that typically both spherulite and shish-kebab morphologies are present during polymer processing, the relative crystallinity evolution can be represented by including both the quiescent undisturbed volume Φ 0 and that of the flow-induced Φ 0,f :

α(t) = 1 -exp(-Φ 0 (t) -Φ 0,f (t)) (1.22a) Φ 0 (t) + Φ 0,f (t) = -ln(1 -α(t)).
(1.22b)

The equation system given by equations 1.21 and 1.22 is presented using the function Υ as defined by Eder et al. [START_REF] Eder | Crystallization processes in quiescent and moving polymer melts under heat transfer conditions[END_REF]. However, this can be easily generalized by using: A different approach based on thermodynamics is used by Poitou et al. [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF] to describe the flowinduced crystallization in the framework of generalized standard materials. The standard material formalism is commonly used in solid mechanics to describe various coupled phenomena. It necessitates two potentials, a thermodynamic potential and a pseudo-potential, in order to describe the behavior of a material. The first potential helps in quantifying the capability of the material to store energy whereas the second potential quantifies the capability of the material to dissipate energy [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF]. In this case, Poitou et al. [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF] used the free energy Ψ as the thermodynamic potential, defined as:

Φ3,f =8π Ṅf (1.23a) Φ 3,f =8πN f . ( 1
Ψ(T, α) = αX ∞ Ψ c (T ) + (1 -αX ∞ )Ψ a (T ) (1.25)
where Ψ c and Ψ a correspond to the free energy of the crystalline and amorphous phases, respectively.

By using the standard material formalism, it is therefore possible to fully describe a coupled phenomena such as the flow-induced crystallization since the mechanical parameters are dependent on the degree of crystallinity. The coupling is taken into account by adding up the potential representing the quiescent kinetics given by the Nakamura model [START_REF] Nakamura | Some aspects of nonisothermal crystallization of polymers. I. Relationship between crystallization temperature, crystallinity, and cooling conditions[END_REF] and the potential referring to the mechanical constitutive behavior [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF]. By assuming that the material is a Newtonian fluid, this mechanical dissipation potential is thus approximated using a simple relation between the strain rate tensor and the viscosity. For more detailed information concerning the derivation of this model, refer to [START_REF] Poitou | Cristallisation induite par écoulement ou par déformation d'un polymère -une approche thermodynamique[END_REF][START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF]. The model is given by:

α(t) = (1 -α) [ - 1 β(T ) ∂η ∂α trD 2 + n χ(T ) (-ln(1 -α)) 1-1 n ] , (1.26)
where D is the strain rate tensor, β(T ) is a function taking into account the temperature-dependency of the flow-induced crystallization and ∂η ∂α is the viscosity derivative in respect to the relative crystallinity α. χ(T ) is the temperature-dependent kinetic function as defined by Hieber [START_REF] Hieber | Correlations for the quiescent crystallization kinetics of isotactic polypropylene and poly (ethylene terephthalate)[END_REF]. The first term in equation 1.26 is the contribution of the flow to the relative crystallinity evolution whereas the second term is the original Nakamura model given in equation 1.14 with χ(T ) = 1/K(T ).

POLYMER CRYSTALLIZATION

Γ = 1 + b 1 exp ( -b 2 /α b 3 ) Empirical Tanner [70] Γ = 1 + b 1 α + b 2 α 2
Empirical, based on suspension Zheng et al. [START_REF] Zheng | Modeling of flow-induced crystallization of colored polypropylene in injection molding[END_REF] Γ

= 1 + (α/b 1 ) b 2 / (1 -α/b 1 ) b 3 Empirical

Effect on Other Material Properties

During processing, the development and evolution of the different crystal structures influence other material properties such as viscosity, specific volume, thermal capacity and thermal conductivity.

Therefore, in the last two decades, various efforts have been made to develop models that describe the effect of crystallization on these material properties. In the following, some of these models are presented for describing viscosity and specific volume. For a more detailed look into the effect of crystallization on the thermal capacity and the thermal conductivity refer to [START_REF] Kennedy | Flow analysis of injection molds[END_REF].

Viscosity

For semi-crystalline materials, an accurate prediction of the solidification of the material during processing is imperative to achieving precise simulation results. However, as the solidification is mainly due to crystallization in such materials, specifying a single no-flow temperature is highly inaccurate especially since the crystallization is affected by both the thermal and flow history [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. An alternative way to determine when solidification happens is by using a solidification criterion based on the viscosity.

For that, it is then important to couple the viscosity and the relative crystallinity of the material to obtain more realistic results.

As it is difficult to measure the viscosity and the crystallinity simultaneously, separate testing methodologies are typically performed to couple these two phenomena. Rheological measurements are done to describe the viscosity evolution as a function of temperature and shear rate. Whereas, DSC 1.3. POLYMER CRYSTALLIZATION measurements are carried out to obtain the crystallinity evolution with temperature. Special care is taken to assure that these measurements experience the same thermal histories [START_REF] Lamberti | Crystallinity and Linear Rheological Properties of Polymers[END_REF]. The experimental data obtained from the previously mentioned methods are used to obtain models describing the viscosity enhancement due to crystallization. These models or coupling functions are mainly based on two different approaches:

• derived from models developed in the field of suspension rheology,

• empirical equations recreating the abrupt increase in viscosity due to crystallization.

Table 1.2 presents some of the commonly used models to couple the viscosity and the crystallization using the following normalized rheological function:

Γ(α) = η η 0 (1.27)
where η is the melt viscosity and η 0 is the zero-shear-rate viscosity. The list is not comprehensive, refer to [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF] and [START_REF] Lamberti | Crystallinity and Linear Rheological Properties of Polymers[END_REF] for more details concerning this subject.

Specific volume

During polymer processing, the specific volume is mainly affected by the temperature, pressure and flow history of the material which in turn influences the produced part's shrinkage behavior [START_REF] Zuidema | Flow induced crystallization of polymers: application to injection moulding[END_REF].

For semi-crystalline polymers, the specific volume is additionally dependent on the cooling rate and the crystallinity evolution as they affect the shifting of the transition zone. This is represented in Therefore, it is essential to model accurately the specific volume using PVT models that take these dependencies into account. Some efforts have been made in this regard by Luyé et al. [START_REF] Luye | PVT measurement methodology for semicrystalline polymers to simulate injection-molding process[END_REF], Fulchiron et al. [START_REF] Fulchiron | Analysis of the pressure effect on the crystallization kinetics of polypropylene: dilatometric measurements and thermal gradient modeling[END_REF], Zheng et al. [START_REF] Zheng | Modeling of flow-induced crystallization of colored polypropylene in injection molding[END_REF] and Zhao et al. [START_REF] Zhao | A novel method for predicting degrees of crystallinity in injection molding during packing stage[END_REF]. They included the relative crystallinity into the calculation of the specific volume v by assuming a simple mixing rule of the molten and solidified phases' specific volumes, represented respectively as v m and v s . This law is written as: The specific volumes v m and v s are described using the empirical Tait equation, such as:

v = αv s + (1 -α)v m . ( 1
V x (T, P ) = v 0 (T ) [ 1 -C ln ( 1 + P B(T ) )]
(1.29) with x ∈ m, s. In equation 1.29, C is a universal constant equal to 0.0894, v 0 (T ) is the specific volume at zero gauge pressure and B(T ) describes the pressure sensitivity of the studied material. These temperature-dependent functions are defined as follows:

v 0 (T ) = b 1x + b 2x (T -b 5 ), (1.30)
B(T ) = b 3x exp (-b 4x (T -b 5 )) (1.31)
where b 1x , b 2x , b 3x , b 4x and b 5 are data-fitted coefficients usually obtained using dilatometry measurements [START_REF] Troisi | Application of a multi-phase multi-morphology crystallization model to isotactic polypropylenes with different molecular weight distributions[END_REF][START_REF] Zinet | Thermophysical Characterization and Crystallization Kinetics of Semi-Crystalline Polymers[END_REF][START_REF] Roozemond | Flow-induced crystallization of isotactic polypropylene: Modeling formation of multiple crystal phases and morphologies[END_REF].

Another approach to improve the specific volume predictions is to include the cooling rate Ṫ into the PVT model as well as the transition temperature calculations, as developed by Cook et al. [START_REF] Cook | A research framework for cooling rate-dependent PVT models[END_REF],

such as:

T t (P, Ṫ ) = b 5 + b 6 P + c log ( Ṫ Ṫ0 ) (1.32)
where b 6 and c are material parameters and Ṫ0 is a low cooling rate for which the polymer maintains an equilibrium state. Other authors such as Wang et al. [START_REF] Wang | Continuous Two-Domain Equations of State for the Description of the Pressure-Specific Volume-Temperature Behavior of Polymers[END_REF] and Hopmann et al. [START_REF] Hopmann | Prediction and validation of the specific volume for inline warpage control in injection molding[END_REF] proposed more complex cooling rate dependencies for the PVT modeling, for more details refer to [START_REF] Wang | Continuous Two-Domain Equations of State for the Description of the Pressure-Specific Volume-Temperature Behavior of Polymers[END_REF][START_REF] Hopmann | Prediction and validation of the specific volume for inline warpage control in injection molding[END_REF].

SURROGATE MODELING

Surrogate Modeling

Currently, in science and engineering, there still exists a large number of physical problems that high performance modern computing is not able to accurately solve despite the advances made in modeling, numerical and discretization techniques [START_REF] Hijazi | Data-driven POD-Galerkin reduced order model for turbulent flows[END_REF]. This is mainly due to the high complexity and multi-dimensionality of such physical models [START_REF] Chinesta | PGD-Based Computational Vademecum for Efficient Design, Optimization and Control[END_REF]. Therefore, in the last couple of decades, there has been a growing interest in surrogate modeling techniques in order to address these challenges.

Such methods approximate the response of complex models using a surrogate model also known as a metamodel after being trained with a limited number of input and output results. This makes surrogate models cheaper to run and are thus used instead of the complex model in various fields such as engineering design optimization, uncertainty quantification, sensitivity analysis ...

General Overview

The construction of a surrogate model consists of multiple steps. The basic process can be summarized as follows [START_REF] Wang | Review of Metamodeling Techniques in Support of Engineering Design Optimization[END_REF][START_REF] Forrester | Recent advances in surrogate-based optimization[END_REF][START_REF] Iuliano | Proper Orthogonal Decomposition, surrogate modelling and evolutionary optimization in aerodynamic design[END_REF]:

1. Design Variables Choice: Selection of variables, which presumably have a non-negligible impact on the model output, this choice is usually supported by preliminary experiments, whether physical or numerical experiments. needing to be approximated by the surrogate model. Following this, the surrogate model is chosen along with the appropriate fitting method. In the following sections, a closer look is taken into each of the four steps shown in Figure 1.4.

Design of Experiments

Design of experiments (DoE) or sampling plan are terms used to represent the physical or computational experiments needed to be run in order to capture the behavior of an underlying system over a limited number of variables [START_REF] Forrester | Recent advances in surrogate-based optimization[END_REF][START_REF] Razavi | Review of surrogate modeling in water resources: REVIEW[END_REF]. To build a surrogate model of a system, it is crucial to cleverly select the design points in order to cover the complete design space using the least possible number of samples. There exist two main categories of DoEs:

• Classical DoE methods: used typically for physical experiments since these experiments are nonrepetitive due to the presence of random error sources (full/fractional factorial, Box-Behnken, 1.4. SURROGATE MODELING central composite...)

• Modern DoE methods: applied to deterministic computer experiments where systematic errors are mostly involved (Taguchi method, Latin hypercube, uniform designs, Monte Carlo...).

There exists a high number of DoE methods present in literature, however, in this thesis, Latin

Hypercube sampling and Quasi-Monte Carlo sampling are used and are presented in detail in the following sections. For more information concerning the other DoE techniques, refer to [START_REF] Gao | Intelligent methods for the process parameter determination of plastic injection molding[END_REF][START_REF] Simpson | Metamodels for Computer-based Engineering Design: Survey and recommendations[END_REF][START_REF] Giunta | Overview of Modern Design of Experiments Methods for Computational Simulations[END_REF].

Latin Hypercube Sampling

The Latin hypercube sampling (LHS) technique is widely used in computational applications as it can work with any sample size which gives the user the freedom to choose the number of samples according to the available computational resources [START_REF] Giunta | Overview of Modern Design of Experiments Methods for Computational Simulations[END_REF][START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF]. LHS is a space-filling method, as it tries to find a design that fills the design space given a specific number of samples. To demonstrate the technique, let us consider a 2-dimensional design space Ω = [0, 1] 2 with variables having uniform probability distribution functions. If the desired number of sample points is, for example, N = 4, then the design space is divided into 4 2 equally sized cells. Next, values from 1 through N are placed in each row so that no two columns have an integer repeating, similar to Sudoku [START_REF] Giunta | Overview of Modern Design of Experiments Methods for Computational Simulations[END_REF][START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF]. Figure 1.5

shows two of the many possible arrangements. Next, a random integer between 1 and N is selected to specify the N cells in which a sample point is randomly picked. In the examples shown in Figure 1.5, the chosen integer is 2 such as the shaded cells representing the sampling sites.

The Latin square presented in Figure 1.5 generalizes to a Latin hypercube for higher dimensional (n > 2) design spaces. For instance, a design space with n = p design variables requiring N sample points will form N p hypercubes [START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF]. It is also worth noting that there can exist some cases, such as diagonal arrangements, for which the chosen sample sites are not optimally positioned and do not fill the design space. One way to eliminate such an arrangement is by introducing additional conditions to check the minimum distance between the design points and selecting the largest one [START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF].

Quasi-Monte Carlo Sampling

Quasi-Monte Carlo (QMC) sampling is the deterministic counterpart of the classical Monte Carlo (MC) sampling method [START_REF] Jank | Quasi-Monte Carlo sampling to improve the efficiency of Monte Carlo EM[END_REF]. In classical MC, sample points are randomly selected in the design space for a given interval which can lead to having some regions of the design space unexplored [START_REF] Metropolis | The Monte Carlo Method[END_REF]. Therefore, to overcome this shortcoming, QMC methods were developed to provide deterministic sample points with the optimal spread over the design space. These methods are also known as lowdiscrepancy sequences since they fill the space with some uniformity in order not to leave big gaps [START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF].

There exist various numbers of such sequences such as the Halton sequence [START_REF] Halton | On the efficiency of certain quasi-random sequences of points in evaluating multi-dimensional integrals[END_REF], the Niederreiter sequence [START_REF] Niederreiter | Society for Industrial and Applied Mathematics, Random number generation and quasi-Monte Carlo methods: Author[END_REF], the Sobol sequence [START_REF] Sobol | On the distribution of points in a cube and the approximate evaluation of integrals[END_REF]... Many of the mentioned sequences are based on the van der Corput one, which is the simplest one-dimensional low-discrepancy sequence. One can refer to [START_REF] Faure | From van der Corput to modern constructions of sequences for quasi-Monte Carlo rules[END_REF] for a detailed explanation concerning this sequence and its generalizations. For this work, the Sobol sequence is mainly used to obtain quasi-random DoE samples. In this type of sequence, the prime number 2 is used as the base for all the dimensions of the sequence. The first dimension is the van der Corput sequence with a base of 2 and the higher dimensions are permutations of this first dimension [START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF][START_REF] Sobol | On the distribution of points in a cube and the approximate evaluation of integrals[END_REF]. Figure 1.6a shows an example of 1000 sample points obtained using QMC with the Sobol sequencing method in comparison to those gotten using LHS in Figure 1.6b for a 2-dimensional design space Ω = [0, 1] 2 with variables having uniform probability distribution functions.

Model Order Reduction

Model order reduction (MOR) techniques provide a way to reduce the complexity of high-fidelity models by representing them in a reduced form. This allows for a more efficient model evaluation for offline and online applications. MOR methods work on the discretization of a state equation's dimensionality instead of on its design space which differentiates them from data fitting methods [START_REF] Iuliano | Proper Orthogonal Decomposition, surrogate modelling and evolutionary optimization in aerodynamic design[END_REF].

There exist multiple MOR techniques in the literature such as Proper Orthogonal Decomposition (POD), Proper Generalized Decomposition (PGD), Krylov subspace methods... In the following, the discussion is limited to describing the POD method for the case of a time-dependent problem. For this discussion, let us consider Y to be one model output and X = (X 1 , • • • , X p ) the set of input parameters defined in the design space D X i such as:

Y = f (X) = f (X 1 , • • • , X p ) ≈ F (X 1 , • • • , X p ) (1.33)
where the function f (X) represents the system or simulation output that is being approximated using a surrogate model F (X) and p the parameter space design size. These notations are used to describe the methods presented in the following sections for the case of one output result which can be extended to the multi-output case.

SURROGATE MODELING

Proper Orthogonal Decomposition

The POD technique is one of the most commonly used model order reduction method. It is known by different names such as the Karhumen-Loève (KL) expansion as well as Principle Component Analysis (PCA) and an extension of the Singular Value Decomposition [START_REF] Nouy | A priori model reduction through Proper Generalized Decomposition for solving time-dependent partial differential equations[END_REF]. POD is a mathematical procedure that provides orthonormal basis functions, known as empirical eigenvectors, in order to obtain a simplified representation of a set of data or a state's evolution [START_REF] Iuliano | Proper Orthogonal Decomposition, surrogate modelling and evolutionary optimization in aerodynamic design[END_REF][START_REF] Ryckelynck | On the a priori model reduction: Overview and recent developments[END_REF]. These obtained eigenvectors correspond to the highest eigenvalues and they represent the basis functions that are able to describe the main modes or events involved in a certain state evolution [START_REF] Ryckelynck | On the a priori model reduction: Overview and recent developments[END_REF]. Let M be a model state variable considered in a specific system. The variable can be represented at specific time steps and for the different nodes present using a matrix A M such as:

A M = ⎛ ⎜ ⎜ ⎜ ⎜ ⎝ M 1 1 M 2 1 • • • M κ 1 M 1 2 M 2 2 • • • M κ 2 . . . . . . . . . . . . M 1 ν M 2 ν • • • M κ ν ⎞ ⎟ ⎟ ⎟ ⎟ ⎠ (1.34)
where κ is the number of time steps and ν is the number of nodes. The main objective of the POD technique is to find a set of orthogonal basis functions φ i (i = 1, • • • , κ) able to describe the deviation Mi of the model variables from their mean Mi .

These basis functions can be represented through a linear relation with the deviation of the model variables as follows:

φ = κ ∑ i=1 a i Mi . (1.35)
To obtain the POD basis vectors, also known as POD modes, the following eigenproblem needs be solved to obtain a subspace with a low dimension m able to provide a good approximation of the true data:

AA T φ = λ i φ (1.36)
with λ i as the eigenvalues. The state variable M can be then represented using a linear combination of the calculated POD basis functions with the following relation:

M (x, t) = M + m ∑ i=1 α i (t)φ i (x) (1.37)
where α i (t) are the POD coefficients [START_REF] Du | POD reducedorder unstructured mesh modeling applied to 2d and 3d fluid flow[END_REF]. More in-depth information concerning POD and model order reduction can be found in [START_REF] Iuliano | Proper Orthogonal Decomposition, surrogate modelling and evolutionary optimization in aerodynamic design[END_REF][START_REF] Nouy | A priori model reduction through Proper Generalized Decomposition for solving time-dependent partial differential equations[END_REF][START_REF] Du | POD reducedorder unstructured mesh modeling applied to 2d and 3d fluid flow[END_REF].

SURROGATE MODELING

The POD model by itself is unable to approximate the state variable at sites not included in the original data set [START_REF] Iuliano | Proper Orthogonal Decomposition, surrogate modelling and evolutionary optimization in aerodynamic design[END_REF]. Therefore, the next step is to generate a surrogate model to predict the POD coefficients. For this, the calculated coefficients are used to train a surrogate model. This is done by first choosing a model to fit the POD parameters such as Kriging, radial basis functions, polynomial functions then use an appropriate fitting method to determine the surrogate model's coefficients. As such, one can obtain a POD-based surrogate model.

Surrogate Modeling Methods

The first approximation methods such as polynomial response surface models and neural networks appeared prior to 1990 and were mainly applied to structural design problems [START_REF] Simpson | Design and analysis of computer experiments in multidisciplinary design optimization: a review of how far we have come-or not[END_REF]. Later on in the 1990's, interest in such surrogate modeling and approximation techniques increasingly grew.

Researchers started investigating additional methods such as higher order response surface models, polynomial models [START_REF] Simpson | Design and analysis of computer experiments in multidisciplinary design optimization: a review of how far we have come-or not[END_REF]. At the beginning of the 21 st century, alternative surrogate modeling techniques emerged such as Kriging, splines, radial basis functions, support vector machine, polynomial chaos expansions...

Response surface methodology (RSM) represents the high-fidelity model response as:

Y = g(X) + ϵ (1.38)
where ϵ is the normally distributed random error with zero mean and standard deviation [START_REF] Alizadeh | Managing computational complexity using surrogate models: a critical review[END_REF] and g(X) is a polynomial function used to approximate the response. Typically low-order polynomials are used to express g(X) [START_REF] Simpson | Metamodels for Computer-based Engineering Design: Survey and recommendations[END_REF].

Kriging is an interpolation-based method derived from statistical theory [START_REF] Liu | Integration optimization of molding and service for injection-molded product[END_REF]. A Kriging surrogate model is defined as a combination of a polynomial part and a random error part such as:

Y = g(X) + Z(X) (1.39)
where Z(X) is a normally distributed Gaussian random process with a zero mean, variance σ 2 and a nonzero covariance [START_REF] Chen | Multi-objective optimal approach for injection molding based on surrogate model and particle swarm optimization algorithm[END_REF][START_REF] Simpson | Metamodels for Computer-based Engineering Design: Survey and recommendations[END_REF]. The function g(X) is similarly to RSM defined using a polynomial expression of the design variables X.

Radial basis functions (RBF) provide an approximation F (X) for an output response such as:

F (X) = p ∑ i=1 w i ϱ(∥ X -X i ∥) (1.40)
where w i is a weight coefficient, ϱ is a nonlinear basis function and ∥ • ∥ represents the Euclidean distance [START_REF] Li | A modified global optimization method based on surrogate model and its application in packing profile optimization of injection molding process[END_REF].

Artificial neural networks (ANN) constitute a special architecture of neurons that are each defined by linear regression models along with a nonlinear transform to approximate Y [START_REF] Simpson | Metamodels for Computer-based Engineering Design: Survey and recommendations[END_REF]. They learn the output response through the training process while adjusting their weights until an error threshold is reached [START_REF] Zhao | Intelligent Injection Molding on Sensing, Optimization, and Control[END_REF].

Simpson et al. [START_REF] Simpson | Metamodels for Computer-based Engineering Design: Survey and recommendations[END_REF] as well as Wang and Shan [START_REF] Wang | Review of Metamodeling Techniques in Support of Engineering Design Optimization[END_REF] provide additional details on the various surrogate modeling methods available in the literature. In the following section, the polynomial chaos expansions are presented more in detail.

Polynomial Chaos Expansion

Polynomial chaos expansions (PCE) are surrogate modeling techniques that expand Y = f (X) using a series of multivariate basis functions orthogonal to the probability density function g X i of an input variable X i [START_REF] Teixeira | Adaptive approaches in metamodel-based reliability analysis: A review[END_REF]. In general, polynomial chaos expansions can be considered as a special case of KL or PCA since orthogonal polynomial families are used to describe the basis functions instead of eigenvectors. The resulting random model response Y is assumed to have a finite variance and thus belonging to the so-called Hilbert space, allowing for the following spectral representation [101]:

Y = ∞ ∑ j=0 c j ψ j . (1.41)
The random variable Y is an infinite series, where {ψ j } ∞ j=0 are a set of countable random variables forming the basis of the Hilbert space and {c j } ∞ j=0 are the coefficients representing the coordinates of Y in this basis [101]. With the assumption that the input variables are independent, an inner product can be defined for each variable X i with any two functions ϕ 1 , ϕ 2 such as:

⟨ϕ 1 , ϕ 2 ⟩ i = ∫ D X i ϕ 1 (x)ϕ 2 (x)g X i (x)dx.
(1.42)

When replacing the arbitrary functions in equation 1.42 with orthogonal polynomials

P (i) k : ⟨P (i) j , P (i) k ⟩ i = ∫ D X i P (i) j (x)P (i) k (x)g X i (x)dx = a (i) j δ jk (1.43)
where k is the polynomial degree and δ jk is the Kronecker symbol equal to 1 for j = k and 0 otherwise.

As for the term a

(i) j , it is the squared norm of P (i)
j and is equal to 1 for orthonormal polynomials:

a (i) j = ||P (i) j || 2 i = ⟨P (i) j , P (i) j ⟩ i .
(1.44)

SURROGATE MODELING

The family of orthonormal polynomials {ψ (i) j } is obtained by normalizing the P (i) j functions :

ψ (i) j = P (i) j / √ a (i) j , i = 1, 2, • • • , p. (1.45)
Depending on the distribution of the input variables, there exist well-known orthogonal polynomial families. For example, if X i has a uniform distribution, the corresponding family is the Legendre polynomials or if X i has a Beta distribution then the Jacobi polynomials constitute the basis functions for PCE [START_REF] Crestaux | Polynomial chaos expansion for sensitivity analysis[END_REF].

In order to estimate the polynomial chaos coefficients, there exists various intrusive (Galerkin projection) and non-intrusive methods. Popular non-intrusive techniques are error minimization ones which solve a least squares (LS) or least angle regression (LAR) problem [START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF]. For a more detailed description concerning PCE, the reader is advised to refer to [START_REF] Teixeira | Adaptive approaches in metamodel-based reliability analysis: A review[END_REF][101][START_REF] Crestaux | Polynomial chaos expansion for sensitivity analysis[END_REF].

Fitting Techniques

The different methods presented in Sections 1.4.3 and 1.4.4 provide a way to approximate a model's response through various mathematical relations. To be able to use these equations to calculate this approximated response, the parameters or bases need to be determined. This is done by a variety of fitting techniques some of which are listed in Figure 1.4.

A commonly used fitting technique for polynomial functions is the least-squares method with which the needed parameters are obtained by minimizing the difference between the approximated response and the high-fidelity one [START_REF] Alizadeh | Managing computational complexity using surrogate models: a critical review[END_REF]. Another highly utilized method is back-propagation usually used to fit data for an ANN. Whereas, when generating a model using Kriging, the best linear unbiased predictor is typically used [START_REF] Simpson | Metamodels for Computer-based Engineering Design: Survey and recommendations[END_REF].

1.4.6 Some Applications

Variance-Based Sensitivity Analysis

A surrogate model can be easily utilized to efficiently perform sensitivity analysis on the chosen design variables. One type of such analyses is variance-based ones which are used to quantify the variance contribution of an input parameter to the unconditional variance of the model output [START_REF] Nossent | Sobol sensitivity analysis of a complex environmental model[END_REF].

A commonly utilized method in this field is the Sobol method [START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates[END_REF]. An attractive feature of this technique is its ability to not only measure the amount of variance caused by one input but also the 1.4. SURROGATE MODELING interaction of two or more inputs and their contribution to the output. These are known as the Sobol sensitivity indices. The method utilizes approximate Monte Carlo integration to calculate the different indices [START_REF] Rosolem | A fully multiple-criteria implementation of the Sobol method for parameter sensitivity analysis: MULTICRITERIA SOBOL SENSITIVITY ANALYSIS[END_REF].

The Sobol approach decomposes the function f into terms of increasing dimensionality [START_REF] Nossent | Sobol sensitivity analysis of a complex environmental model[END_REF] f

(X) = f 0 + p ∑ i=1 f i (X i ) + p ∑ i=1 p ∑ j=i+1 f ij (X i , X j ) + • • • + f 1,••• ,p (X 1 , • • • , X p ) (1.46)
such as each successive term represents the increasing degrees of interactions between the various parameters. The total variance V (Y ) can be then defined using the sum of the partial variances through a similar decomposition to equation 1.46 with the assumption that the parameters are mutually orthogonal:

V (Y ) = p ∑ i=1 V i + p-1 ∑ i=1 p ∑ j=i+1 V ij + • • • + V 1,••• ,p . (1.47)
As such, the Sobol sensitivity indices can be formulated as follows:

S i = V i V (1.48)
for the first order sensitivity indices and

S ij = V ij V (1.49)
for the second order sensitivity indices.

Optimization

Many engineering problems revolve around optimization routines whether it is to obtain a better geometrical design for a part or to improve the production of that part by finding the optimal process parameters or to better predict the production process using a simulation with accurate model parameters... Surrogate models provide a means to perform these optimizations in the order of seconds instead of hours or even days. A typical surrogate model based optimization problem can be represented as:

min f (X) s.t. mi (X) ≤ 0 (i = 1, ..., p) X ∈ [X L , X U ]

SURROGATE MODELING IN INJECTION MOLDING SIMULATION

where the tilde symbol signifies the surrogate model of a function, X L and X U are the lower and upper bounds of the design variables, respectively [START_REF] Wang | Review of Metamodeling Techniques in Support of Engineering Design Optimization[END_REF].

Uncertainty Quantification

In every numerical simulation, there exist uncertainties in the input. These can be anything from dimensional errors, constitutive material properties to parameter uncertainties [START_REF] Mcclarren | Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers[END_REF]. Uncertainty quantification (UQ) is the study of how these uncertainties affect the simulation output. To perform such a study, it is necessary to run a high number of simulations while varying the selected quantities of interest (QoI). This step is highly challenging since it is computationally very expensive. To overcome this challenge, surrogate modeling is used to represent the simulation output using the chosen QoI.

The generated surrogate model is then used to study the propagation of the input uncertainties to the simulation results using various statistical methods.

Surrogate Modeling in Injection Molding Simulation

In the last two decades, the use of surrogate models to approximate outputs from the injection molding simulation has been growing steadily especially in the field of process parameter optimization to enhance product quality and molding efficiency. Gao and Wang [5] employed a Kriging approximation model along with an adaptive optimization technique to minimize the warpage in produced parts by varying process parameters such as the mold and melt temperature, injection time as well as the holding pressure and time. Similar works were performed by Chen et al. [START_REF] Chen | Multi-objective optimal approach for injection molding based on surrogate model and particle swarm optimization algorithm[END_REF], Wang et al. [START_REF] Wang | Warpage optimization with dynamic injection molding technology and sequential optimization method[END_REF] and Kang et al. [START_REF] Kang | Metamodel-based design optimization of injection molding process variables and gates of an automotive glove box for enhancing its quality[END_REF]. Others used radial basis function [2,[START_REF] Li | A modified global optimization method based on surrogate model and its application in packing profile optimization of injection molding process[END_REF][START_REF] Kitayama | Multi-objective optimization of injection molding process parameters for short cycle time and warpage reduction using conformal cooling channel[END_REF], artificial neural networks [START_REF] Cheng | Multiobjective optimization of injection molding parameters based on soft computing and variable complexity method[END_REF][START_REF] Shi | A warpage optimization method for injection molding using artificial neural network with parametric sampling evaluation strategy[END_REF],

Gaussian process [START_REF] Xia | An enhanced optimization approach based on Gaussian process surrogate model for process control in injection molding[END_REF] as surrogate modeling technique to optimize process parameters for controlling shrinkage and warpage in the final part. Additional applications for surrogate models seen in the literature were used for the optimization of cycle time [START_REF] Kitayama | Multi-objective optimization of injection molding process parameters for short cycle time and warpage reduction using conformal cooling channel[END_REF] and part weight [START_REF] Villarreal-Marroquín | A comparison of two metamodel-based methodologies for multiple criteria simulation optimization using an injection molding case study[END_REF].

All of the above-mentioned publications use a surrogate model to perform an optimization to determine optimal process parameters. However, another interesting utilization of surrogate modeling in injection molding was recently published by Ivan et al. [START_REF] Ivan | Improving Numerical Modeling Accuracy for Fiber Orientation and Mechanical Properties of Injection Molded Glass Fiber Reinforced Thermoplastics[END_REF] where the surrogate model is used to identify two fiber orientation model parameters. The authors used experimental fiber orientation data obtained by micro-computed tomography to calibrate the fiber orientation model defined in

Summary and Thesis Structure

As already stated in Section 1.2, injection molding simulation provides a cost-efficient means to design and improve part quality before producing a mold as well as to identify optimal process settings for the production. However, to be able to achieve that, simulations have to be highly accurate and include less and less approximations. In the case of semi-crystalline thermoplastics, this means taking into account the complex crystallization process taking place during the processing of the material along with its effect on other properties such as viscosity and specific volume. Additionally, in order to improve the overall fill predictions as well as the warpage and shrinkage estimations, the solidification criterion needs to be coupled to the crystallization instead of using a constant no-flow temperature value to determine solidification. Nevertheless, as seen in Section 1.3, the crystallization process under flow conditions is complex and the physics behind it is still not yet fully understood making this problem not straightforward to solve. Also, most of the available models are phenomenological ones with a high number of parameters that are difficult to determine experimentally due to the inability of recreating the extreme conditions experienced during the injection molding process using current experimental devices. One way to overcome this drawback is by using a surrogate model to represent a simulation output which can be measured experimentally on the injection molding machine and using it to identify the model parameters such as Ivan et al. [START_REF] Ivan | Improving Numerical Modeling Accuracy for Fiber Orientation and Mechanical Properties of Injection Molded Glass Fiber Reinforced Thermoplastics[END_REF] did for the fiber orientation model parameters.

Therefore, in this thesis, the main aim is to implement a crystallization model into a 3D injection simulation routine while taking into consideration the effect of crystallization on the viscosity, PVT and solidification models. In order to identify the experimentally-challenging-to-determine parameters, surrogate models are generated to represent the pressure results obtained by the simulation and to perform a calibration using experimental pressure signal measurements obtained during injection molding trials.

To reach this goal, the theory and the current state of the art in the field of injection molding, polymer crystallization and surrogate modeling were first presented in this chapter. Chapter 1 finished off by introducing some works from the literature that utilized surrogate modeling to represent injection molding simulation results similar to what will be presented in this work.

Following this, the material characterization along with the performed injection molding trials are 1.6. SUMMARY AND THESIS STRUCTURE described in Chapter 2. Additionally, in this chapter, some of the experimental pressure signals are analyzed for the different geometries considered in this thesis as well as some micrographs produced to observe the crytalline morphology along the thickness of the produced parts.

In Chapter 3, a feasibility study is performed to test the use of surrogate modeling for the purpose of calibrating some modeling parameters in the injection molding simulation. This study compares three different surrogate modeling approaches and tests the effect of the number of training simulations on the accuracy of the predicted pressure results.

And, finally, in Chapter 4, the implementation of a thermo-mechanical based crystallization model is described along with crystallization-dependent viscosity, PVT and solidification models in Moldflow R

⃝ to simulate the injection molding of an unreinforced thermoplastic material. Five modeling parameters are identified using three generated surrogate models by a multi-objective optimization routine. The simulation results using the calibrated models are then analyzed for different processing conditions and mold geometries.

Chapter 2

Material Characterization and Injection Molding Trials

Material Description

The polymer material used in this work is an industrial grade of an unreinforced polyoxymethylene (POM) homopolymer (-CH 2 -O-) n in its granular form with a density of 1.42 g/cm 3 . The polymer is also known as polyacetal and is an engineering semi-crytalline thermoplastic material widely used in the automotive industry due to its high stiffness, low friction and excellent anti-wear properties [START_REF] Santis | Effect of shear flow on spherulitic growth and nucleation rates of polypropylene[END_REF][START_REF] Li | Development of Polyoxymethylene/Polylactide Blends for a Potentially Biodegradable Material: Crystallization Kinetics, Lifespan Prediction, and Enzymatic Degradation Behavior[END_REF]. Additionally, in comparison to other semi-crystalline polymers, POM presents with a high degree of crystallinity around 75 to 85 % [START_REF] Santis | Effect of shear flow on spherulitic growth and nucleation rates of polypropylene[END_REF][START_REF] Durmus | Effect of polyhedral oligomeric silsesquioxane (POSS) reinforced polypropylene (PP) nanocomposite on the microstructure and isothermal crystallization kinetics of polyoxymethylene (POM)[END_REF].

Differential Scanning Calorimetry

Technique Overview

Differential scanning calorimetry (DSC) is a thermal analysis technique used to measure a thermophysical property of a material as a function of temperature [START_REF] Sandler | Differential Scanning Calorimetry[END_REF]. A studied sample along with a reference are subjected to a predetermined temperature profile during which the difference in heat input between the samples is measured [START_REF] Cheremisinoff | Thermal Analysis[END_REF]. From these measurements, various physical characteristics of the polymer can be quantified, such as the melting enthalpy and the melting temperature... Additionally, the measured enthalpy changes can be utilized to evaluate the relative crystallinity evolution as a function of temperature and time then crystallization kinetics.

Experimental Set-up

The DSC measurements are performed using a TA Instruments, Inc. DSC Q1000 machine with nitrogen as purge gas with a flow rate of 50 mL/min. The granular form of POM is used for the experiments. The prepared samples are encapsulated in Aluminum pans and lids after being weighed.

Ten runs are performed with five different applied cooling rates (1, 3, 5, 10, 40 • C/min) as summarized in Table 2.1. As can be seen in Table 2.1, the sample weights vary for the different cooling rates. This is done as a compromise between the heat flux measurement sensitivity and the low temperature gradient in the samples. Before starting the measurement runs, the DSC machine is calibrated using an Indium sample having a weight of 10.36 mg. In total, three calibration runs are performed:

1. at 1 • C/min for the measurement run with a cooling rate of 1 To set-up the different runs, the following steps are followed:

1. Equilibrate at 25 • C.

2. Isothermal for 2 min.

Results

The measured heat flows of the DSC runs summarized in Table 2.1 are analyzed in order to quantify important physical and kinetic characteristics of the crystallization of the studied POM material. from the heating cycles of the DSC measurements, the melting temperature (T m ) of the material is determined. Table 2.2 summarizes these quantities obtained from the DSC data. In practice, the crystallization onset temperature depends linearly to the logarithm of its corresponding cooling rate. Therefore, in Figure 2.3, the mean onset temperature of each studied cooling rate ( Ṫ ) is plotted as a function of ln( Ṫ ). However, as can be seen in this figure, the Linear fit 1 used to fit all five data points is not optimal as that of Linear fit 2 which fits only the first four points. This shows the possibility of having inaccurate DSC results especially with the 40 • C/min cooling rate run. Therefore, to determine the crystallization temperature and enthalpy of the studied POM material, Runs 9 and 10 are not taken into consideration as they showed the previously described deviation. Whereas for the determination of the melting temperature, all runs are considered as the heating rate was constant at 10 • C/min for all of the performed runs.

To get the relative crystallinity evolution from the DSC runs, it is imperative to determine an accurate baseline from which the area under the crystallization peak is calculated. This is done using a MATLAB R2019b code with an example illustration presented in Figure 2.4. After running the MATLAB script, the user is asked to specify four points; two before and two after the crystallization peak (turquoise crosses in Figures 2.4a and 2.4b). These points are used to fit two different lines before and after the peak (pink dotted lines). The user is then asked to specify the start and end of the peak (purple circles). In between these two points, a first approximation of the relative crystallinity is calculated and used to correct the purple baseline by the turquoise one in Figure 2.4 below the 2.2. DIFFERENTIAL SCANNING CALORIMETRY Table 2.2 -Overview of the material characteristics of the studied POM determined using the DSC measurements. crystallization peak.

Run Ṫ ( • C/min) T c,onset ( • C) T c ( • C) T m ( • C) ∆H c (J/g)
To compute the relative crystallinity evolution for the different cooling rates, one run is chosen to represent each cooling rate since the error between the different runs at the same cooling rate was considered negligible. The relative crystallinity α is computed using the corrected baseline by:

α(t or T ) = Area below

Experimental Set-up

The dynamic temperature ramp tests are performed using an ARES rheometer from TA Instruments, Inc. with a plate-plate geometry. Two ramp tests are carried out one with a constant frequency of 1 rad/s and a 0.3 % strain and another with a frequency of 3 rad/s and a 0.03 % strain. The heating/cooling rate used is 3 • C/min. To set-up the test, the following steps are followed:

1. Heat up the apparatus to 210 • C.

2. Set the zero gap between the plates.

3. Open the furnace, increase the gap and place the POM granules on the lower plate.

4. Reheat the apparatus.

5. Decrease the gap to 0.7 mm and clean-up in case of polymer overflow.

6. Run the test using the software interface. • C and 147 • C. Since these measurements were performed with a cooling rate of 3

Results

Injection Molding Experiments 2.4.1 Mold Description

A simple mold geometry is used to perform injection molding trials to obtain experimental data that are needed to calibrate model parameters in the injection molding simulation presented in the next chapters. The mold is designed to fulfill the following requirements:

• a long curved channel with the possibility of producing short shots,

• two different mold thicknesses (3 mm & 1.5 mm) and lengths.

To obtain the final design, numerical DoEs are carried out in Moldflow R ⃝ where results such as the freezing time of the gate and the pressure at multiple locations are compared and deemed acceptable for the presented design. The produced mold is shown in Figure 2.8 whereas Figure 2.9 shows a sketch to note is that for the full length part with a 1.5 mm thickness, it is not possible to fill the cavity completely with the studied POM material under the recommended processing conditions. Therefore, for this geometry a shot volume of 32 cm 3 is used for the trials but it does not represent the cavity's volume.

Design of Experiments

The same electrical injection molding machine is used to run the experiments used in this work.

For all the performed runs, a melt temperature of 220 • C is chosen and the switch-over condition is set by a volume of 10 cm 3 . The mold temperature, the injection velocity and the holding pressure are varied in these runs. The upper and lower bounds of these process settings are different for the four available mold geometries. This is due to the geometrical changes affecting the injection pressure and solidification of the material. Table 2.4 summarizes the bounds of these variables. A design of experiments (DoE) is defined for each geometry where the process settings are varied between their upper and lower bounds as well as the middle values. These DoE are presented in Tables 2.5 to 2.8.

Fewer runs are performed for the 3 mm full length part due to the small injection velocity range (5 to 15 cm 3 /s) for which the injection process can happen without crossing the machine's maximum pressure limit of 240 MPa.

Pressure Signals

As mentioned in Section 2.4.1, pressure and temperature sensors are located at four different positions in the mold. In this work, the main focus is on the pressure signal results measuring the pressure in the cavity before the gate (P1), after the gate (P2), before half the length of the channel (P3) and near the end of the channel (P4). Figures 2.12, 2.13 and 2.14 show some example pressure 1. Velocity-controlled filling phase: the pressure increases linearly until the switch-over condition is reached.

2. Pressure-controlled filling phase: a sudden pressure drop due to the usually smaller specified holding pressure than the one experienced during the velocity-controlled filling phase.

3. Packing phase: a gradual increase of pressure as the packing starts, followed by a constant pressure observed before the gate whereas after the gate a gradual decrease is observed as the material solidifies.

4. End of packing phase: a sudden decrease in pressure before the gate due to a partially solidified material followed by a gradual decrease to reach zero pressure as seen also after the gate. 2.5), respectively. It is apparent that the low injection velocity of 5 cm 3 /s for the S3 run presented in Figure 2.12a led to a long part filling time of around 9 s. Additionally, for the S3 run, the pressure plateau seen after the switch-over highlights that the part was not fully filled before switching to the pressure-controlled phase and therefore the holding pressure is first used to complete the filling and afterward for the packing. This behavior is not observed in the S6 run presented in Figure 2.12b where the injection velocity is higher at 15 cm 3 /s and the filling is nearly finished before switching over. In both examples, there is an expected 20 MPa pressure drop at P1 from the set holding pressure. of a partial remelting of the gate leading more material to enter the cavity. Additionally, near the end of the packing phase, a slight increase in the P1 signal is observed in all example plots in Figure 2.13 affecting in turn the pressure decay slope in P2 and P3. This behavior can be explained by an over-compacting before the end of the packing phase. Similarly as seen in Figure 2.12, since we use a fixed switch-over condition for all experimental runs, a low injection velocity leads to having a less than fully filled part before switching and therefore the pressure plateau seen in Figure 2.13d for V22 in comparison to the smooth switch-over seen in Figure 2.13c for V15.

As already mentioned in Section 2.4.2, it is not possible to produce fully filled parts using the full length mold cavity with a 1.5 mm thickness. Therefore, all runs performed with this geometry under the conditions presented in Table 2.7 created short shoted U-shaped channels. One example pressure signal result is given in Figure 2.14a obtained during the W3 run. In this plot, P0 represents the specific injection pressure and as can be seen it reaches the machine's maximum injection pressure of 240 MPa without being able to fill the part. However, when using the insert to cut the channel length by half, it became possible to fill the cavity for most of the process settings summarized in Table 2.8.

One representative result is presented in Figure 2.14b for the T20 run. One can observe the immense loss of pressure in the cavity by observing the P3 signal. This is due to the thinness of the cavity leading to faster polymer solidification.

Polarized Optical Microscopy

In order to visualize the different crystal morphologies that developed during the injection molding process, microtome cuts are made at the sensor locations from four produced samples and analyzed by polarized optical microscopy. The samples are taken from the injection of the half length 3 mm thick cavity having V1, V3, V19 and V21 as process settings (refer to Table 2.6). Figure 2.16a shows the location of the microtome cuts in the case of a fully filled part. Whereas, in the case of a short shot such as for V1 and V19, the third cut is made near the end of the part as shown in Figure 2.16b and annotated as Cut 3.

Figure 2.17 shows the results for V1 at the second pressure sensor located after the gate (P2) for three different magnifications (Figure 2.17a: 25x, Figure 2.17b: 50x, Figure 2.17c: 100x). As a first observation, it is apparent that as we move away from the wall to the core of the part the spherulites' size increases since the crystals have time to form and grow under no shear conditions. Another observation is the presence of a translucent skin layer indicating that the crystallites in this region are smaller than the wavelength of the utilized light. These small crystallites in this layer are postulated to have been formed under high shear conditions during the filling stage. In Figure 2.17, this layer's thickness is around 0.2 to 0.3 mm. Table 2.9 presents a summary of the approximately measured thicknesses for the different samples and cuts.

In addition to having

Conclusion

In this chapter, the material used in this work was presented along with two types of material characterization experiments. The first type was the DSC measurements which provided information about the material's crystallization behavior under different cooling rates. This made it possible to determine the melting temperature of the material at 178.88 • C and the crystallization enthalpy of 155.35 J/g. The relative crystallinity evolution with temperature and time under the various cooling rates was also determined and will be used to quantify the crystallization kinetics of POM in Chapter 4.

The second type of experiment is a dynamic temperature ramp test used to identify the temperaturedependency of the viscosity and observe its behavior after the onset of crystallization. This data will similarly be used to determine the viscosity model parameters later on in the work.

A mold was designed to provide the possibility to produce four different geometries varying by length and/or thickness in addition to measuring the pressure signals during the injection molding process at four locations in the cavity. To set-up the injection molding trials, filling and packing studies were performed followed by the official trials where for each geometry three process settings are varied. The variable settings are the mold temperature, the injection velocity and the holding pressure. Some characteristic pressure signals were analyzed in this chapter to give an overview of the obtained results that are used in this work for the calibration of model parameters.

Some of the injection molded samples were studied using a polarized optical microscope to get a deeper insight into the crystal morphologies and the different layers that form during the process.

Microtome cuts were made at the sensor locations for different process settings to study the effect of mold temperature and injection velocity on the morphology and skin layer thickness. By observing

INTRODUCTION

Introduction

In this chapter, the feasibility of using surrogate modeling to calibrate parameters in an injection molding simulation is studied. Three modeling parameters are chosen for the calibration in this study.

These are the heat transfer coefficients during filling and packing as well as the pressure-dependency coefficient in the Cross-WLF viscosity model [START_REF] Cross | Rheology of non-Newtonian fluids: A new flow equation for pseudoplastic systems[END_REF][START_REF] Williams | The Temperature Dependence of Relaxation Mechanisms in Amorphous Polymers and Other Glass-forming Liquids[END_REF]. The before-mentioned parameters are known to affect the pressure results in the simulation and their experimental determination is challenging and time-consuming [START_REF] Nguyen-Chung | Determination of the heat transfer coefficient from short-shots studies and precise simulation of microinjection molding[END_REF][START_REF] Stricker | Determination of heat transfer coefficients at the polymermold-interface for injection molding simulation by means of calorimetry[END_REF][START_REF] Kleindel | The Elastic Mold Deformation During the Filling and Packing Stage of the Injection Molding Process[END_REF]. For this calibration, the half length 3 mm thick geometry presented in Figure 2.11 is utilized and the three modeling parameters are calibrated using the experimental P2 pressure sensor data (refer to Figure 2.9).

In this feasibility study, the accuracy of three regression-based surrogate modeling techniques and two design of experiment methods in their ability to approximate the pressure signal at a single node in the simulation is studied. The proposed modeling approaches are the following: The models are trained using 30, 60 and 120 simulation runs and their predictions are assessed accordingly. The three models are built by varying a total of six input parameters, including process settings and modeling variables. The process parameters are the injection velocity, coolant inlet temperature, holding pressure and the modeling parameters, as already mentioned, are the heat transfer coefficients during filling and packing as well as the pressure-dependent coefficient in the Cross-WLF viscosity model. After obtaining the surrogate models, a sensitivity analysis is performed to study the contribution of each varied parameter to the pressure signal as well as their interaction with each other. Finally, using the measured experimental pressure data and one of the generated surrogate models, the modeling parameters are calibrated and the error between the optimized simulation and the default one is compared.

SOFTWARE TOOLS

Software Tools

Injection Molding Simulation

The injection molding simulations are performed using Autodesk R ⃝ Moldflow R ⃝ Insight 2021.1 (AMI2021.1). The simulations are automated with the help of the three command line control functions. This set of utilities allow the use of third-party software in tandem with the simulation runs in order to automatize extensive studies such as optimizations and sensitivity analyses... The main functions are editing, running, and retrieving results from a user-customized Moldflow R ⃝ analysis automatically through a command line interface. The following are the available utilities:

1. studymod : This command generates a modified simulation study from a base study using an XML format modifier file (modified geometry, boundary condition, process/model parameters).

runstudy:

This command launches a new simulation analysis.

studyrlt:

This command extracts simulation results from a finalized analysis in simple text or XML file format.

Python Uncertainty Quantification Library

The Python Uncertainty Quantification (pyUQ) library is a Bosch proprietary Python tool providing the user the ability to perform uncertainty quantification (UQ) on simulation models using state-of-the-art UQ models and methods. The methods cover the design of experiments, surrogate modeling, sensitivity analysis and statistical analysis. The pyUQ module provides functions that can be modified or extended to act as an interface between simulation tools and UQ methods. Parameter uncertainties are defined in this tool using standard distributions (uniform, normal, Gamma, Beta) by stating the parameter's upper and lower limits. The simulation results are stored in a Hierarchical Data Format (HDF5) to be able to efficiently perform statistical analysis or surrogate model generation using this tool. HDF5 is a folder-like storage structure, which saves data in compressed format and allows data slicing leading to efficient memory usage. Apart from these features, the main highlight of this tool is its active learning algorithm with an adaptive DoE generation which is one of the methods evaluated in this work.

MATLAB

The POD model order reduction and the non-linear regression of the respective POD-based coefficients are implemented in MATLAB R2019b. The fmincon built-in function of the optimization toolbox is used for the optimization of the POD-based surrogate model.

Simulation Environment

The 3D high-fidelity simulation model is set-up in AMI2021. The various surrogate models are obtained by modifying six parameters in the simulation. Three of which are the process parameters changed during the experimental trials as presented in Table 2.4

for the half length 3 mm thick part (HL3). Since a cooling analysis is performed in the simulations,

PROPOSED DOE AND SURROGATE MODELING TECHNIQUES

Proposed DoE and Surrogate Modeling Techniques

In this study, three different combinations of DoEs and surrogate modeling methods are used in order to compare their performance in predicting the pressure signal results at one surface node corresponding to the location of the second sensor in Figure 2.9. Additionally, the total number of simulations used for the training phase is varied such as 30, 60 and 120 simulation runs are considered and their performance is tested using four additional runs. The testing simulations are summarized in Table 3.2.

POD-NLR

The POD-based surrogate model whose parameters are fitted using non-linear regression (NLR) is referred to as POD-NLR. To generate this model, LHS is used to obtain the DoEs according to the needed total number of training simulations S. The obtaining of the final model requires the following steps after running all the simulations according to the DoE:

1. Pre-processing: the various nodal pressure results P i (i = 1, • • • , S) are pre-processed by shifting all corresponding time series by t start so that all non-zero signals start at t 0 = 0 s and followed by the resampling of the pressure data using a common time vector from 0 to 50 s with 0.02 s time steps. This step is important given that the different simulation settings lead to different times at which the flow front reaches the sensor node and produces a non-zero pressure signal.

2. Model order reduction: the POD basis functions φ are calculated by solving the eigenvalue problem P P T V = λV where the eigenvectors V corresponding to the most influential eigenvalues

PROPOSED DOE AND SURROGATE MODELING TECHNIQUES

λ constitute the basis functions able to reconstruct the pressure signal by:

P (s) reconstructed = n ∑ i=1 Γ (s) n φ n (3.1)
where s is the simulation number, n is the number of modes or basis obtained according to a specified error value and Γ (s) n is the POD parameter for a specific mode n and simulation s. The truncation criterion is done according to an error value of ϵ = 10 -3 which only selects the modes corresponding to eigenvalues that fulfill the condition λ λmax ≥ ϵ.

3. Model fitting: a least-squares regression of a second-order polynomial is used to train a model to predict the POD parameters Γ and the time shift value t start as follows:

Γ (s) n = a (n) + b (n) i X (s) i + c (n) ij X (s) i X (s) j (3.2) t (s) start = d + e i X (s) i + f ij X (s) i X (s) j (3.3)
where X i and X j are the surrogate model input variables with i ̸ = j representing the number of variables (1, • • • , 6) and a, b, c, d, e, f are the surrogate parameters. In equations 3.2 and 3.3 the implicit Einstein summation convention is used for indexes i and j. This implementation is restricted to bilinear regression in order to mitigate overfitting the data.

Regular-PCE-LS

The same LHS DoEs used to train the POD-NLR model are employed as inputs to a second-order PCE-based surrogate model. The model parameters are trained using least-squares regression and thus referred to as regular-PCE-LS since the DoE is fully predefined. The PCE model is generated using the pyUQ library using the following steps:

1. Pre-processing: similar to POD-NLR, a resampling step is performed to have a unique time index starting with t 0 = 0 s and spanning to 50 s with 0.02 s time steps. signifies that all interaction effects are considered in the PCE basis. Additionally, the PCE regression order is estimated by pyUQ using heuristics and in this case a second-order regression is used.

Adaptive-PCE-LS

A sequential design of experiment (SDOE) is used to build iteratively a DoE from a big candidate DoE obtained with QMC according to the importance of the input variables. The surrogate model utilizes PCE and is trained using least-squares regression. This model is generated using the pyUQ library and is referred to as adaptive-PCE-LS since its DoE is enriched as the model is being generated making it an adaptive approach. Therefore, given a sample set, the surrogate model generation follows the same steps presented in regular-PCE-LS. The main difference between the two methods lies in the iterative steps needed to build the whole surrogate model by updating the DoE.

The algorithm starts with an initial small LHS DoE of 10 simulations and is enriched after each iteration with 5 additional sample points from a large QMC candidate DoE (1000 samples). The enrichment is done via the bootstrap method [START_REF] Efron | Nonparametric estimates of standard error: The jackknife, the bootstrap and other methods[END_REF] which calculates the local variances of the surrogate model and chooses the new sample points that lead to the maximum variance. In other words, it performs a UQ analysis on the simulation results to estimate which parameters contribute most to the statistical fluctuations. This process is repeated until convergence in terms of maximum simulations S of 120 or a global error threshold E of 10 -2 is achieved.

The surrogate model is assessed by cross-validation using the training data set with the leaveone-out error method such as the model is generated from S -1 simulations and evaluated on the remaining one simulation result. The standard error metric used to evaluate the surrogate model after each iteration is its global error. This error value is an average of a scalar model output predefined by the user. In this work, the selected scalar quantity corresponds to the mean value of the pressure vector at the studied node.

Polynomial chaos expansion

Least-squares regression

The flowchart in Figure 3.2 presents a summary of the steps followed while generating the adaptive-PCE-LS surrogate model using an active learning/adaptive approach.

For the sake of clarity, the characteristics of the three studied surrogate models are summarized in Table 3.3.

Evaluation Strategy

To evaluate and compare the various surrogate models, the root mean squared error (RMSE) is used. Let the surrogate model response and the true values of the high-fidelity simulation model or experiment be represented by Ŷ and Y , respectively. Using these notations, the RMSE is computed by evaluating the surrogate model pointwise on a given test data set, such as:

RM SE(Y, Ŷ ) =    √ 1 n n ∑ i=1 (Y i -Ŷ i ) 2 (3.4)
where n is the total number of data points. In this work, this error is normalized using the min-max normalization method to be able to compare the various surrogate modeling techniques:

N ormalized RM SE(Y, Ŷ ) = RM SE(Y, Ŷ ) Y max -Y min .
(3.5)

Optimization

In this feasibility study, an optimization routine is performed to identify the uncertain simulation model parameters: HT C f illing , HT C packing and D 3 . For these purposes, the sensor data from the 27 experimental runs presented in Table 2.

Results and Discussion

Comparison of Surrogate Models

The performance of the three surrogate model techniques is assessed by their ability to predict the pressure signal at a surface node in simulations unseen by the model. In this work, those test simulations are presented in Table 3.2. Although the total simulated time is around 48 s, the main focus is on the pressure signal between 0 and 22 s. This interval includes the filling phase between 0 and ∼2 s, the packing phase till ∼18 s and the start of the cooling phase until the pressure goes to 0 MPa before 22 s. Therefore, all the error metrics and plots consider only this time interval.

Sensitivity Analysis

A Sobol sensitivity analysis is performed using the generated adaptive-PCE-LS surrogate model trained with 120 simulation runs. pressure development. Interestingly, at the time when pressure strongly drops during the packing phase, probably associated to the freezing of the gate, pressure development is controlled by the interaction between the holding pressure and the coolant inlet temperature. The interaction between HT C packing and D 3 after 20 s is difficult to explain in terms of polymer physics, but eventually reveals some subtleties caught in the complexity of the simulated phenomena.

In general, such sensitivity analyses are helpful in understanding the effect of parameters on specific output results. One way they can play an important role in surrogate model generation is by pinpointing the most influential parameters at the regions in which the model exhibits the large errors.

For example, at the end of the filling phase, the models are not fully capable to predict the correct pressure drop, as can be seen in Figures 3.7c and 3.7d . A possible way to improve these predictions is by adding DoE points where we mainly vary the injection velocity, holding pressure and D 3 as their first and second order sensitivity indices show their high contribution to the output in this region.

Modeling Parameter Calibration

All three models performed similarly after training with 120 high-fidelity simulations as shown in One of the main objectives of this optimization is the identification of material parameters or boundary conditions for enabling a more accurate estimation of the pressure field in a high-fidelity injection molding simulation. To show the impact of this parameter calibration, the experimental DoE is simulated using both the optimized parameters and the default parameters given in Table 3.4.

The simulated pressure signals at the sensor location (taken from a surface node) are compared to the experimental pressure signals at that position. Four example DoEs are shown in Figure 3.12. In all cases, the simulation using the optimized parameters is closer to the experimental data than using the default parameters. The experimental setting showed in Figure 3.12a corresponds to V1 (refer to Table 2.6), for which the cavity was not filled completely (short shot) as seen in Figure 2.15a and both default and optimized simulations were not able to predict this event. This failure in the simulation can be related to the default solidification criterion, which defines the melt-to-solid transition at a constant temperature and oversimplifies the crystallization phenomena occurring in a wide range of cooling rates. In any case, the simulation results using the optimized parameters show more accurate pressure estimations during the packing phase as seen in Figures 3.12b and 3.12d. However, the estimation of the filling time and maximum pressure at the end of filling are less accurate than those obtained with default parameters. This issue is certainly due to the high value of the optimized pressure-dependent viscosity coefficient, which leads to higher viscosities and thus longer filling times. The prediction of the time when the pressure goes to zero is also prone to further improvement for both default and optimized simulations. This shortcoming should be also strongly related to the current simple solidification criterion. The normalized RMSE comparing the pressure predictions using the default and optimized parameters for all 27 experimental DoEs are given in Figure 3.13. The use of the optimized parameters led to a huge decrease in the RMSE over all considered DoE cases. The simulation using default values exhibited for some cases errors around 40 % which are decreased by more than half while using the optimized parameters. However, both simulations are unable to recreate the short shots that were observed experimentally for the molding trials with a low injection velocity of 10 cm 3 /s and a holding pressure of 20 MPa, as seen in Figure 3.12a. Since the optimized value for pressure-dependent viscosity parameter D 3 is found at the upper bound of the surrogate model generation range, it can indicate that the chosen intervals are insufficient or that the implemented models are intrinsically limited to mimic all physical phenomena during the injection molding process.

Conclusion

Surrogate models of high-fidelity simulations, where the variables are material-dependent parameters, offer an alternative to standard experimental identification methods, as shown for the case of injection molding in this feasibility study. The main aim of the method is the generation of an accurate surrogate model, which enables the reverse engineering of the material-dependent parameters using optimization techniques. Considering six independent variables, the POD-NLR and the regular-PCE-LS surrogate models exhibit basically the same performance especially when using a low number of

Introduction

In the previous chapter, the feasibility of using surrogate modeling to identify three modeling parameters was studied. The study showed that it is possible to optimize these parameters using experimental pressure results with the help of a surrogate model representing the simulated pressure results at one surface node. Therefore, in this chapter, a similar methodology is followed in order to calibrate a flow-induced crystallization model along with crystallization-dependent PVT, viscosity and solidification models. These models are implemented in Moldflow R ⃝ with the help of a user code.

In total, three surrogate models are used to identify the five uncertain modeling parameters and these surrogates are generated using an analogous POD-NLR method as the one presented in Chapter 3.

Such as the pressure signals are first reduced using the proper orthogonal decomposition technique followed by a least-squares regression of a second-order polynomial. This methodology is chosen as it exhibited a good performance when using a low number of training data and a low number of input parameters. After obtaining the optimized modeling parameters, the simulation performance is assessed by studying how well it is able to predict the pressure signals, filling status and skin layer thickness using the user-defined models in comparison to the experimental results. Additionally, the calibrated models are tested using three various mold geometries in order to check the generality of the identified parameters.

Description of the Implemented Models

Crystallization Model

The implemented crystallization model is a thermo-mechanical based model developed by Poitou et al. [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF] in the framework of irreversible thermodynamics given by equation 1.26. This model is chosen since it takes into account both the flow-induced and quiescent crystallization in addition to having a limited number of parameters to identify. In Poitou et al.'s model [START_REF] Poitou | Crystallization of polymers under strain: from molecular properties to macroscopic models[END_REF], the crystallization's kinetics are modeled using a temperature-dependent function χ(T ) defined by Hieber [START_REF] Hieber | Correlations for the quiescent crystallization kinetics of isotactic polypropylene and poly (ethylene terephthalate)[END_REF]. However, in this work, this relation is substituted by the Avrami kinetic constant which is proportional to the Avrami growth rate such as k Avrami ∝ G 3 making it possible to be described using the Hoffman-Lauritzen model and related to the previous function by χ(T ) = 1/k Avrami (T ) 1/n . Additionally, the β is considered a constant value to be determined and not dependent on temperature as presented in

DESCRIPTION OF THE IMPLEMENTED MODELS

the original model. The crystallization model only applies for temperatures lower than the material's melting temperature T m (P ) and can be rewritten as follows for T < T m (P ):

α(t) = (1 -α) [ - 1 β ∂η ∂α trD 2 + nk Avrami (T ) 1 n (-ln(1 -α)) 1-1 n ] , (4.1) k Avrami (T ) = K 0 exp ( - 3K g T (T 0 m (P ) -T ) ) exp ( - 3U ⋆ R(T -T ∞ (P )) ) (4.2)
where β, n, K 0 and K g are data-fitted parameters, D is the strain rate tensor, T 0 m (P ) is the pressuredependent equilibrium melting temperature, U ⋆ is the activation energy for segmental jump of polymer molecules with a universal value of 6270 J/mol, R is the gas constant and T ∞ = T g (P )-30 with T g (P )

as the pressure-dependent glass transition temperature. The first term in equation 4.1 represents the contribution of the flow-induced crystallization to the relative crystallinity evolution α by computing the viscosity derivative in respect to the relative crystallinity α represented by ∂η ∂α . Whereas the second term in equation 4.1 is the contribution of the quiescent crystallization to α.

In order to reduce the number of model parameters, a linear relation is deduced between the kinetic constants of the Hoffman-Lauritzen equation for the used POM material such as:

K g = 6732 ln(K 0 ) + 11296. (4.3)
This dependency is obtained by fitting equation 4.2 to the Avrami kinetic constant calculated using the Avrami-Ozawa relation which relates the kinetic constants of both models using:

k Avrami (T ) = [ - d dT ( k Ozawa (T ) 1/n ) ] n (4.4)
where n is the Avrami exponent with a typical value of 3 representing instantaneous nucleation, same parameter as in equation 4.1. In this work, since the DSC measurements presented in Section 2.2 are performed under non-isothermal conditions, k Ozawa is first calculated under these conditions and following that the k Avrami is obtained using equation 4.4.

Additionally, the pressure dependencies of the temperatures are defined by:

T m (P ) = T m (0) + b 6 P, ( 4.5)
T 0 m (P ) = T 0 m (0) + aP, ( T g (P ) = T g (0) + bP (

where b 6 , a and b are usually experimentally-determined parameters. The pressure dependency of the melting temperature b 6 is typically obtained from PVT measurements. However, the two other dependencies are harder to determine experimentally for semi-crystalline materials. Therefore, for this work, a relation between the three parameters is proposed by computing the crystal growth curves as function of temperature for different pressure levels using the Hoffman-Lauritzen model given in equation 4.2. For various combinations of a and b, the maximal growth rate is determined for the different pressures and the corresponding temperatures are used to calculate b 6 . Figure 4.1 shows the surface representing the relation between the dependencies such as:

b 6 = f • a + g • b (4.8)
where f = 0.629 and g = 0.4015 are the obtained fitted parameters. The melting temperature and the equilibrium melting temperature are assumed to have the same pressure dependency such as b 6 = a = 0.175 K/MPa (value obtained from material's supplier). As for the dependence of the glass transition temperature, it can be deduced using equation 4.8 and found to be b = 0.161 K/MPa.

DESCRIPTION OF THE IMPLEMENTED MODELS

Table 4.1 -The model parameters of the implemented crystallization model.

Model Parameter Value (units)

n 3 (-) U ⋆ 6270 (J/mol) T 0 m (0) 199.7 ( • C) [124]
T m (0) 178.9 ( • C) As the polymer experiences a phase change while crystallizing, heat is released into the system, otherwise known as the latent heat of crystallization. Therefore, an additional heat source is defined that expresses the temperature increment ∆T as a function of crystallization such as:

T g (0) -73.0 ( • C) [113] b 6 0.175 (K/MPa) a 0.175 (K/MPa) b 0.161 (K/MPa)
∆T = α(t)∆t ∆H c c p (4.9)
where ∆t is the time step, ∆H c is the crystallization enthalpy and c p is the specific heat. If the temperature of the melt is higher than the material's melting temperature (T > T m (P )), remelting of the crystals takes place leading to a negative ∆T = -α ∆H c c p .

The crystallization enthalpy is given in Table 4.2. As for the specific heat, it is typically defined using Moldflow R ⃝ material cards as a function of temperature including a crystallization peak. However, since in this implementation a crystallization-dependent a heat source is included, the specific heat is defined using a linear relation of the values presented in Table 4.3.

DESCRIPTION OF THE IMPLEMENTED MODELS

Table 4.3 -The specific heat data used in the Moldflow simulations along with self-defined models.

Temperature ( • C) Specific Heat (J kg -1 K -1 )
32 1286 220 1932

Viscosity Model

A modified Cross-WLF model is used to describe the viscosity in this implementation. This model takes into account the crystallization dependency in addition to the temperature, shear rate and pressure dependencies of the viscosity. The melt viscosity η is given by: η(T, P, α, γ) = η 0 (T, P, α)

1 + ( η 0 (T,P,α) γ τ ⋆ ) 1-λv (4.10)
where γ is the shear rate, τ ⋆ is the critical stress level at the transition to shear thinning, λ v is the power law index in the high shear rate regime and η 0 is the zero shear viscosity defined as:

η 0 (T, P, α) = D 1 exp ( -A 1 (T -T ref ) A 2 + (T -T ref ) ) ϑ(α) (4.11)
where T ref = T g (P ) with a pressure-dependency defined in equation 4.7 and D 1 , A 1 and A 2 are data-fitted coefficients. As for ϑ(α), it is a function describing the crystallization-dependency of the viscosity and is an extension of Kitano et al.'s relation [START_REF] Kitano | An empirical equation of the relative viscosity of polymer melts filled with various inorganic fillers[END_REF] for concentrated suspension of particles.

This function is defined as:

ϑ(α) = 1 (1 -α A ) B (4.12)
where A represents the maximum relative crystallinity reached by the material before it solidifies and B is a data-fitted exponent. This coupling function differs from the one defined in equation 1.27 since the function ϑ(α) = η 0 (T, P, α)/η 0 (T, P ). This choice is made to increase the effect of crystallization on viscosity at low shear rates as done by Pantani et al. [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF].

By excluding the crystallization-dependency ϑ(α), this Cross-WLF model still differs from the one implemented by default in Moldflow R ⃝ , shown in Appendix A, since the A 2 coefficient is not dependent on pressure in this implementation. Including this dependency as done in default Moldflow R ⃝ leads to a linear relation between η 0 and P for all temperatures. However, Rudolph et al. [START_REF] Rudolph | WLF model for the pressure dependence of zero shear viscosity of polycarbonate[END_REF] showed a nonlinear dependency between them for low temperatures approaching the glass transition temperature of a polycarbonate material. Therefore, in this work, the main connection between the viscosity and pressure is modeled using the glass transition temperature's dependency.

DESCRIPTION OF THE IMPLEMENTED MODELS

To identify the parameters of the WLF model, the dynamic temperature ramp test performed at a constant frequency of 3 rad/s and a 0.03 % strain with a heating/cooling rate of 3 • C/min shown in Figure 2.7 is used. The viscosity results before the onset of crystallization at 152.4 • C are utilized to determine the temperature dependency of the viscosity by fitting the following linear relation:

1 ln(η ⋆ /D 1 ) = - A 2 A 1 1 (T -T ref ) - 1 A 1 (4.13)
where η ⋆ is the measured complex viscosity from the rheological test. Additionally, to describe the shear rate dependency of the viscosity using the Cross model given in equation 4.10, Moldflow R ⃝ 's material database is used. The model parameters for the Cross-WLF model are summarized in Table 4.4.

PVT Model

Since this work deals with a semi-crystalline thermoplastic material, the Pressure-Volume-Temperature (PVT) model is defined in relation to the relative crystallinity α. For these purposes, a simple two-phase system is assumed. This makes it possible to describe the specific volume v using a mixing law of the molten and solidified phases' specific volumes, represented respectively as v m and v s . This law is written as:

v = αv s + (1 -α)v m . ( 4.14)
The specific volumes v m and v s are described using the empirical Tait equation given in summarized in Table 4.5.

Solidification Model

A solidification criterion is used in an injection molding simulation to determine whether the polymer is solidified or still flowing. This is usually defined by a no-flow temperature characteristic to a material. However, the solidification of the material is highly dependent on the crystallization degree present in it. Therefore, in this implementation, a crystallization-dependent solidification model is utilized such as:

• if (α > A): solidified polymer • if (α < A): flowing polymer
where A is the same parameter used to describe the dependency of the viscosity on crystallization using the function given in equation 4.12. The Solver Application Programming Interface (API) in Moldflow R ⃝ enables the user to create their own C++ functions which the solver can call to calculate a property during an analysis. In AMI2021.1, this feature allows user-defined models for viscosity, PVT, core shift, solidification and fiber orientation. Additionally, due to the newly-introduced Advection API functionality, it is possible to define material derivatives that are solved in time and space by the solver while providing the user's code access to the solution of the derivative(s).

To be able to use the Solver API feature, a shared library object needs to be compiled and copied into the binary directory of the AMI installation folder. The project should be compiled with all the header and C++ files needed to get a successful run when the Solver API option is enabled in Moldflow R ⃝ . Most importantly, this should include the header files defining the user and utility functions which represent the main link between Moldflow R ⃝ 's solver and the user's code:

• Utility functions are functions that are already defined inside the Moldflow R ⃝ solver which can be called by the user's code. They are used in order to gain access to information available to

• solidification criterion: SolverUserHb3dSolidificationAtNode.

These functions necessitate access to the model parameters, this is done by defining the path of the text files containing these parameters in the Moldflow R ⃝ simulation. Example codes for these functions can be found in the installed directory folder of Moldflow R ⃝ Insight and code explanations on the Autodesk website [START_REF]Autodesk Knowledge Network[END_REF].

To sum up, the flowchart presented in Figure 4.2 is provided as an overview of the implemented code which is compiled to a shared library file in the form of a dynamic-link library (DLL) to be used during the various Moldflow R ⃝ simulation runs presented in this chapter.

Simulation Environment

Similarly to the simulation model used in Chapter 3, the half length 3 mm thick geometry is again used to perform the simulation runs in this implementation. The 3D meshed model along with the feed system and cooling channels are presented in Figure 3.1. For the material data, the POM material card from Moldflow R ⃝ 's database is edited to have a linear specific heat relation with temperature as given in Table 4.3. The process settings of the base simulation are similarly defined using imported data from the ENGEL sim link software tool.

Additionally, since in this work multiple user-defined models are implemented, the Solver API option should be enabled as well as selecting a user-defined viscosity, PVT and solidification models.

Since the crystallization model is defined using the Advection functionality, it does not have its own drop-down option but if the SolverUserHb3dUserNodeFieldsEvolution user function is compiled in the code it will be taken into account. For the model parameters, text file paths need to be provided by the user in the parameter text strings in Process Settings > Advanced options... > Solver API window in Moldflow R ⃝ .

In this implementation, a single simulation requires around 110 minutes to be completed on a workstation with a 4.10 GHz processor and 32 GB RAM. In order to generate the different surrogate models, the automation utilities available in Moldflow R ⃝ are utilized to run efficiently the multiple simulations needed, refer to Section 3.2.1 for more information concerning this topic.

β (Pa -1 s -1 ) -10 12 -10 10 ln(K 0 ) (-) 23 83 A (-) 0.3 0.6 B (-) 2 5 HT C p (W m -2 • C -1 )
500 5000

Surrogate Models

Input Variables and Output

A total of three surrogate models are generated for three different processing conditions corresponding to the ones used during the experimental injection molding runs. These process settings are presented in Table 4.7 along with the surrogate model's name according to which experimental DoE it represents. For each surrogate model, five input parameters are varied: two crystallization model parameters (β, K 0 ), two viscosity model parameters (A, B) and the heat transfer coefficient (HT C p ) used during the packing analysis. These inputs are summarized in Table 4.8 along with their lower and upper limits.

Some of these ranges were determined with the help of heuristics while others are based on experimental results. For K 0 , the DSC data obtained from the measurements presented in Section 2.2 are used. Since only a couple of cooling rates were utilized for these measurements, the temperature range covered is limited making it difficult to obtain one constant value for K 0 . Therefore, it is part of the surrogate model variables and the range fits the six data points calculated using the Avrami-Ozawa relation given in equation 4. 4 and is centralized at ln(K 0 ) = 53; a value determined by Plummer

SURROGATE MODELS

and Kausch [127] for POM. As for the solidification criterion A, in the case of B = 2, the viscosity coupling is analogous to the one given by Metzner [START_REF] Metzner | Rheology of Suspensions in Polymeric Liquids[END_REF] where he specified A = 0.68 for smooth spheres and A ≈ 0.44 for rough compact crystals. Therefore, for this implementation, since there should exist spherical and thread-like morphologies during injection molding, we postulate a rough ±0.15 to the value given for compact crystals. As for the B viscosity exponent parameter, Kitano [START_REF] Kitano | An empirical equation of the relative viscosity of polymer melts filled with various inorganic fillers[END_REF] and

Metzner [START_REF] Metzner | Rheology of Suspensions in Polymeric Liquids[END_REF] defined it to be equal to 2 according to suspension theory. However, the B is varied between 2 and 5 in order to take into account the sudden increase in the POM viscosity as it starts to crystallize as already seen in Figure 2.7 during the dynamic temperature ramp tests.

The β parameter has a strong effect on the contribution of the flow-induced crystallization to the total relative crystallinity in the part, therefore, this parameter interval was not straightforward to define as its high value can lead to short shots due to the fast solidification of the thin gate.

Nevertheless, the β interval was assumed after performing multiple sensitivity studies and generating some test surrogate models in order to obtain a DoE that produces a representative amount of short shots and fully filled parts in order to efficiently calibrate the parameters later on using the optimization routine.

The output result used to train the surrogate models is the pressure signal at a surface node corresponding to the location of the second sensor P2 located directly after the gate shown in Figure 2.9.

Generation Methodology

A Latin Hypercube sampling method is chosen to generate the surrogate models. The number of runs chosen for the training of the models is 132 simulations and for the testing is 20 simulations.

This choice is postulated to be appropriate as five parameters are being varied in this implementation in comparison to the six parameters changed in the feasibility study performed in Chapter 3 where it was shown that the pressure predictions are highly accurate using 120 training simulations.

The generation of the three surrogate models is done in MATLAB R2019b using an analogous POD-NLR methodology to that presented in Section 3.4.1. However, in this implementation, the following differences exist:

1. Pre-processing: the pressure results P i (i = 1, • • • , S) are differently pre-processed as in Section 4.5. CALIBRATION ALGORITHM 3.4.1. In this implementation, the original time series are normalized between 0 and 1 according to the determined t start and t end by:

t norm = t original -t start t end (4.15)
where t original is the imported unprocessed simulation time, t start represents the time at which the flow front reaches the sensor node and produces a non-zero pressure signal and t end is the time at which the pressure signal goes to zero at the end of the packing stage. This normalization step is followed by the resampling of the pressure data using a common time vector from 0 to 1 s with 0.001 s time steps.

2. Model order reduction: the POD basis functions φ are calculated similarly as in Section 3.4.1 using equation 3.1. The only difference is that in this study a truncation error value of ϵ = 5 × 10 -4 is used.

3. Model fitting: a least-squares regression of a second-order polynomial is used to train a model to predict the POD parameters Γ, the time shift values t start and t end . The first two are determined using equations 3.2 and 3.3, respectively, with i = j. Whereas, t end is determined as follows:

t (s) end = g + h i X (s) i + z ij X (s) i X (s) j (4.16)
where X i and X j are the surrogate model input variables with i = j representing the number of variables (1, • • • , 5) and g, h, z are additional surrogate parameters to those present in equations 3.2 and 3.3. The implicit Einstein summation convention is used for indexes i and j in equation 4.16.

Calibration Algorithm

To identify the five modeling parameters presented in Table 4.8, a multi-objective optimization routine is performed using the lsqnonlin built-in MATLAB function from the optimization toolbox.

The experimental pressure signals are used as reference results. The performed optimization utilizes 4.5. CALIBRATION ALGORITHM all three surrogate models to obtain one set of optimized parameters using the following formulation:

find V = (β, K 0 , A, B, HT C p ) argmin | Ŷ l max -Y exp,k max | (l = 1, 2, 3 & k = 6, 14, 27) | tl end -t exp,k end | s.t. -10 12 ≤ β ≤ -10 10 , 23 ≤ ln(K 0 ) ≤ 83, 0.3 ≤ A ≤ 0.6, 2 ≤ B ≤ 5, 500 ≤ HT C p ≤ 5000.
In total, the optimization consists of six objective functions. Three of them are minimizing the difference between the maximum pressure value obtained experimentally Y exp,k max where k = 6, 14, 27 represents the experimental DoE number (V6, V14, V27) and that predicted by the corresponding surrogate model Ŷ l max with l = 1, 2, 3 the number of the surrogate model (SM1, SM2, SM3). The other three objective functions are minimizing the difference between the predicted tl end defined using equation 4.16 and the experimentally-determined t exp,k end . This value is obtained by fitting a tangent line to the decreasing pressure slope during the end of packing and checking the intercept with the x-axis of that line. This method is used to not take into account the probable over-compacting that was observed experimentally in Figure 2.13 at the end of packing and discussed in Section 2.4.4 since in the simulation it was not possible to recreate this. Table 4.9 summarizes the experimental values used to define the six objective functions. However, it still struggles for some cases to accurately predict the time shift t end as seen in the curve represented in dark blue.

Since it is difficult to compare visually all the training and testing prediction results, the error metrics for the five POD Γ parameters and time shift parameters are calculated and presented in for the testing set pointing that the chosen regression is not a good fit for the parameters in question.

Moving on to the SM2 model, Figure 4.6 provides four examples of the pressure predictions using the SM2 model in comparison to those obtained by the high fidelity simulation from both the training and testing sets. At the first glance, there is some minor differences between the predictions and the function representing the surrogate model as it performs badly using unseen input data during testing. This is also highlighted by the negative R 2 for Γ 3 as well as for t start .

Modeling Parameters Calibration

The three generated surrogate models are used to calibrate the implemented material models in Moldflow R ⃝ by identifying the five parameters β, K 0 , A, B, HT C p using the algorithm described in 4.10.

Starting by the optimized viscosity parameters, the identified solidification criterion parameter A = 0.35 seems reasonable since it is a parameter that takes into account the geometrical effects of the crystallites. According to Kennedy and Zheng [START_REF] Kennedy | Flow analysis of injection molds[END_REF] and Metzner [START_REF] Metzner | Rheology of Suspensions in Polymeric Liquids[END_REF], the criterion is approximately equal to 0.44 for rough compact crystals and since during injection molding non-spherical crystals develop due to shear, an A < 0.44 is expected. As for the exponent parameter B, the optimized value is nearly equal In the next sections, the identified parameters are used in the Moldflow R ⃝ implementation proposed in Section 4.2 and its results are compared to the default Moldflow R ⃝ simulation that does not take crystallization explicitly into account.

Pressure Results

As shown in Table 4.7, each surrogate model is based on one experimental DoE. Additionally, the optimization used the experimental pressure signals at the sensor located after the gate (P2) to identify the modeling parameters given in Table 4. By observing these comparison plots, it is clear that the simulation that takes into account polymer (c) V27

Figure 4.14 -The RMSE of the pressure predictions using the simulation with the calibrated models (Optimization) and the one using the default Moldflow models (Default) for each of the three processing conditions used to generate a surrogate model. crystallization is performing much better than the one with the default models. For all three conditions, the pressure predictions improved drastically at the three sensor locations using the proposed method especially at P2 and P3. The estimations are following the same trend as the experimental results using the optimized parameters which is not the case while using the default settings. One interesting observation is the ability of the proposed models to recreate the constant pressure seen after the switch-over in Figures 4.12a and 4.12b at P1 and P2, respectively, which leads to a delay in the signal reaching P3 as predicted in Figure 4.12c for the V14 condition.

To quantify the improvement in pressure predictions, the normalized RMSE are calculated for the simulations using the default models and the ones using the proposed models. The averaged normalized RMSE of the pressure predictions at P2 using the simulation with the calibrated models (Optimization) and the one using the default Moldflow models (Default) for all 27 experimental processing conditions for the half length 3 mm thick geometry.

As expected, the prediction error has significantly decreased by including the crystallization model along with the crystallization-dependent viscosity, PVT and solidification models with their identified parameters. Major improvement is apparent in the pressure estimation of P2 and P3 for conditions V6 and V14 where the error decreases three folds in the case of P2 and four times for P3. These sensors are located in the part cavity and therefore predicting the pressure at these locations more accurately can lead to more accurate shrinkage and warpage predictions. In the case of the V27 processing condition, the pressure predictions became closer to reality as seen in Figure 4.13, however, quantitatively, the error decreased by less than 5 %.

The above-discussed results seem quite promising, however, since the optimization is based on the three presented experimental results, it is still not yet certain that the identified parameters are able to predict other experimental conditions as well as the previously-shown ones. Therefore, signals at P2 obtained experimentally in comparison to the ones predicted by the default simulation and by the simulation using the proposed models with the optimized parameters for the three conditions in consideration. By observing the presented results, it is apparent that both simulations are not able to recreate the experimental signals as the latter is more complex than the ones used to calibrate the material models. In the experimental data, during the packing phase around 10 s, the pressure signal changes slopes creating a second bump in the signal. This was hypothesized in Section 2.4.4 to be due to a probable partial remelting of the gate leading to additional polymer melt entering the cavity.

Even though the model includes the effect of remelting, the current calibrated models are not able to recreate the shown bump in pressure. This shows the possibility that some model assumptions need to be revised or that the model is not fully calibrated yet. In Figure 4.16b, a slight change of slope is seen at around 7 s but not as prominent as in reality. As a final remark on this topic, the calibrated models are able to better predict the maximum pressure needed during filling as observed in Figure 4.16 and they predict the first bump in the pressure curve during packing more accurately than the default simulation.

Fill Predictions

As briefly mentioned in Section 2.4.4, some experimental DoEs for the HL3 part did not produce a fully filled part therefore a short shot was produced. Figure 4.17a summarizes the filling status of the 27 experimental DoEs. One of the main goals of implementing a crystallization model and crystallization-dependent material models is to be able to predict accurately when a part is fully filled

Skin Layer Thickness

Since the proposed simulation approach includes a crystallization model, one of the simulation results is the relative crystallinity evolution. Therefore, this result is used to determine the skin layer thickness that is solidified during the filling phase and compacting under shear strain. As mentioned by Pantani et al. [START_REF] Pantani | Modeling of morphology evolution in the injection molding process of thermoplastic polymers[END_REF][START_REF] Pantani | Analysis of flow induced crystallization through molecular stretch[END_REF], there exists a critical shear rate value above which the nucleation density and growth rate are affected due to flow. Since this critical value is hard to determine experimentally, a threshold value of 50 s -1 is postulated in this work below which the flow-induced crystallization is considered minimal. This value is obtained by analyzing the simulation results after switch-over 4.6. RESULTS AND DISCUSSION Table 4.11 -The skin layer thicknesses determined using the simulation utilizing the calibrated models across the thickness at the three sensor locations for four processing conditions (V1, V3, V19, V21).

Simulated Skin Layer Thickness (mm)

V1 V3 V19 V21
Cross-section P1 0.40-0.68 0.21 0.23 0.10-0.21 Cross-section P2 0.35 0.08 0.17 0.02 Cross-section P3 0.20 0.08 0.09where the shear rate decreases by two orders of magnitude to reach this threshold value. To determine the skin layer thickness solidified under the high shear conditions, the time at which the shear rate through the thickness at the sensor locations is below this threshold is first determined. Then the relative crystallinity result α is plotted at that time along the thickness of the part at the three sensor locations P1, P2 and P3 specified in Figure 2.9. Following this, the thickness where α > A is read from the plot. The simulation results are summarized in Table 4.11 and compared to the ones obtained using optical microscopy (refer to Section 2.4.5) in Figure 4.18.

By comparing the thicknesses determined experimentally and those from the simulation, it is important to keep in mind that both results are prone to errors and uncertainties. The experimental values could include some observational errors whereas the simulation ones are based on the assumption that the skin layer does not grow further when the shear rate is lower than 50 s -1 . Nevertheless, the simulation predictions are quite accurate with a mean absolute error of 7.5 % and give a good approximation of the skin layer thickness obtained using optical microscopy. Some discrepancies are present, however, the trend of the results is correct such as the thickness is bigger at the cross-section of P1 and it becomes smaller as we move further to P2 and P3. This is consistent with the experimental results.

Calibrated Models Used in Other Geometries

In this section, the calibrated parameters are tested on the three other geometries that are possible to produce with the same mold presented in Figure 2.9. The simulation results obtained using the default Moldflow R ⃝ models are compared to the ones obtained from the implemented models along with the identified parameters. injected with the full length 3 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental). 4.7. CONCLUSION is a relatively high average volume of around 14 V% of material already solidified and the switching to pressure-control filling may have led to a change of pressure distribution in the mold cavity that can explain the discontinuity seen in Figure 4.21a around 0.5 s. However, at the end of the filling stage, both simulations' pressure drops to zero following a linear line without applying the holding pressure of 50 MPa. This behavior and the previous observation can also be due to an artifact in Moldflow R ⃝ and since it appears in both the default and optimized results, it indicates that it is not due to the user's code. This behavior at P1 is encountered only in the optimized simulation for the W10 condition leading to a short shot unlike using the default simulation which predicts incorrectly a completely filled part as seen in Figure 4.21c. An interesting observation is the shape of the pressure signal in the default simulation which resembles the one measured experimentally but does a poor job in correctly recreating the end of the filling phase. As for the pressure predictions at the P3 sensor location for these two conditions presented in Figures 4.21b and 4.21d, the simulation using the proposed models along with the optimized parameters seems to estimate the experimental signals more accurately than the default simulation. Additionally, the simulation using the optimized parameters is able to correctly predict 24 out of the 27 short shots observed experimentally.

In general, the results for this geometry seem to be unreliable making it hard to assess the performance of the simulation using the calibrated model parameters. Therefore, they are not taken into account in the final assessment given in the following section.

Conclusion

In this chapter, the implementation of a crystallization model along with PVT, viscosity and solidification models is performed in Moldflow R ⃝ . To calibrate the implemented models, three surrogate models are generated using model order reduction and least-squares regression. After studying the performance of the generated surrogates, it became apparent that the fitting method using the same second-order polynomial for all POD parameters and time shifts is not optimal. However, the generated surrogate models were still able to predict the pressure signals with high accuracy. Therefore, a multiobjective optimization scheme is performed to identify the five uncertain model parameters and the calibrated models are then tested in the simulation in comparison to the default Moldflow R ⃝ simulation This comparison showed great improvements in the ability of the simulation using the calibrated models to predict the pressure signals more accurately where the prediction error decreased by more than half from the default simulation for the majority of the studied DoEs. Using the simulation with user-defined models, it became possible to predict short shots produced during processing whereas this was not possible with the default simulation. Additionally, to determine the skin layer thickness developed during the filling phase and compacting stage, a threshold shear rate is assumed under which the FIC is considered minimal. The determined thicknesses showed a good agreement with the ones measured experimentally and presented in Chapter 2.

Following this, the calibrated models are tested with the three other available mold geometries. The pressure signals are observed for a couple of processing settings and the simulation predictions using the user models and default ones are compared. As a general observation, for all three geometries, both simulation results performed poorly. By focusing only on the simulation with user-defined models, the HL15 geometry shed the light on a possible low heat transfer coefficient when the mold thickness is divided by half. Whereas the FL3 geometry suffered from fast FIC kinetics leading to early gate solidification which indicates a possible under-estimation of the β exponent parameter. But it can also underline the complexity of the physical phenomena happening in the gate making the original calibration scheme in the HL3 geometry highly sensitive to them and leading to non-generalized model parameters. As the FL15 geometry showed probable Moldflow R ⃝ artifacts, it is not taken into account in this final discussion.

difficult to predict full filled parts. These behaviors could be due to various reasons. One of them is the use of global heat transfer coefficients for the filling and packing phases instead of local ones. This can lead to inaccurate temperature profile predictions affecting the whole heat balance calculations and in turn the crystallization and solidification predictions. Another reason could be that the chosen parameter intervals, for which the surrogate models are generated, are not big enough since some of the optimized parameters lied on the extremities of these intervals. An additional explanation could be that the calibration process was limited by the gate geometry and the complex physical phenomena happening in it such as early solidification due to high shears and fast cooling.

For future works, the redesign of the mold can be considered in which the gate is removed and the experimental runs are performed using a pressure-controlled filling. This could mitigate the probable errors contributed by the inclusion of the gate during the calibration of the models. Additionally, in the case of short shots, the flow length could be included as an extra output result and used in the models' calibration as it can give important information concerning the material solidification.

In the presented work, the measured pressure signal during the injection molding process is used to calibrate the models as it is accurately and easily measured during processing and can be directly related to the simulated pressure. The temperature measurements, on the other hand, are not as straightforward as the previous ones such as in the best case scenario the temperature sensors measure the mold surface temperature. However, this does not correlate well with the simulated temperature.

Therefore, it would be interesting to use machine learning to bridge the gap between these two results in order to be able to use such measured data in a similar application as presented in this thesis.

As an outlook into the topic of surrogate modeling in the field of injection molding simulation, the present work can be extended to include a full time and space model. This could be used to visualize and assess solutions instantaneously while varying the input parameters. Such surrogate models provide the key to performing uncertainty quantification and efficient sensitivity analyses.

Another interesting application is the creation of a digital twin of the injection molding process using these methodologies.

INTRODUCTION

Introduction

Le moulage par injection est l'une des principales techniques de transformation des polymères avec laquelle plus d'un tiers des matériaux thermoplastiques sont produits [START_REF] Osswald | Material science of polymers for engineers[END_REF]. Le processus de moulage par injection est cyclique et idéal pour la production en série de géométries complexes avec des tolérances serrées [2]. En cas de défauts de la pièce pendant la production, un changement des conditions du processus peut ne pas être suffisant et des modifications du moule doivent être appliquées pour résoudre le problème [START_REF] Kennedy | Flow analysis of injection molds[END_REF]. Cependant, cette démarche est coûteuse et prend du temps. Il est donc important de pouvoir prévoir ces problèmes dès la phase de conception de la pièce à l'aide de techniques de Une approche similaire pourrait être utilisée pour calibrer d'autres types de modèles de matériaux complexes tels qu'un modèle de cristallisation. En raison de la complexité de l'identification expérimentale des paramètres de cristallisation dans des conditions de moulage par injection, de nombreux logiciels commerciaux n'en tiennent pas compte, même lorsqu'ils traitent de l'injection de matériaux semi-cristallins. Cela conduit à des inexactitudes dans les prédictions de remplissage ainsi que dans l'estimation du gauchissement et du rétrécissement dans la simulation [START_REF] Kennedy | Flow analysis of injection molds[END_REF].

Au cours des dernières décennies, de grands efforts ont été faits pour modéliser la cristallisation des polymères en tenant compte de l'histoire de l'écoulement de la matière fondue. Pour atteindre cet objectif, la théorie et l'état actuel de l'art dans le domaine du moulage par injection, de la cristallisation des polymères et de la modélisation de substitution ont d'abord été présentés dans ce chapitre. Le chapitre 1 se termine par la présentation de quelques travaux de la littérature qui ont utilisé la modélisation de substitution pour représenter des résultats de simulation de moulage par injection similaires à ceux qui seront présentés dans ce travail.

Ensuite, la caractérisation du matériau ainsi que les essais de moulage par injection réalisés sont décrits dans le chapitre 2. De plus, dans ce chapitre, certains des signaux de pression expérimentaux sont analysés pour les différentes géométries considérées dans cette thèse ainsi que certaines micrographies produites pour observer la morphologie cristalline dans l'épaisseur des pièces produites.

Dans le chapitre 3, une étude de faisabilité est réalisée pour tester l'utilisation de la modélisation de substitution dans le but de calibrer certains paramètres de modélisation dans la simulation de moulage par injection. Cette étude compare trois approches différentes de modélisation de substitution et teste l'effet du nombre de simulations d'entraînement sur la précision des résultats de pression prédits. Cette comparaison a montré de grandes améliorations dans la capacité de la simulation utilisant 5.6. CONCLUSIONS ET PERSPECTIVES les modèles calibrés à prédire les signaux de pression avec plus de précision où l'erreur de prédiction a diminué de plus de la moitié par rapport à la simulation par défaut pour la majorité des DoEs étudiées.

En utilisant la simulation avec des modèles définis par l'utilisateur, il est devenu possible de prédire les coups courts produits pendant le traitement alors que cela n'était pas possible avec la simulation par défaut. En outre, pour déterminer l'épaisseur de la couche cristalline développée pendant la phase de remplissage et l'étape de compactage, on suppose un taux de cisaillement seuil en dessous duquel le FIC est considéré comme minimal. Les épaisseurs déterminées ont montré un bon accord avec celles mesurées expérimentalement et présentées dans le chapitre 2.

Ensuite, les modèles calibrés sont testés avec les trois autres géométries de moule disponibles.

Les signaux de pression sont observés pour quelques paramètres de traitement et les prédictions de simulation utilisant les modèles de l'utilisateur et ceux par défaut sont comparés. D'une manière générale, pour les trois géométries, les deux résultats de simulation ont donné de mauvais résultats.

En se concentrant uniquement sur la simulation avec les modèles définis par l'utilisateur, la géométrie HL15 a mis en lumière un possible faible coefficient de transfert de chaleur lorsque l'épaisseur du moule est divisée par deux. La géométrie FL3, quant à elle, a souffert d'une cinétique FIC rapide conduisant à une solidification précoce de l'entrée, ce qui indique une possible sous-estimation du paramètre exposant β. Mais cela peut également souligner la complexité des phénomènes physiques qui se produisent dans la grille, rendant le schéma de calibration original dans la géométrie HL3 très sensible à ces phénomènes et conduisant à des paramètres de modèle non généralisés. Comme la géométrie FL15 a montré de probables artefacts Moldflow R ⃝ , elle n'est pas prise en compte dans cette discussion finale. Sandra Saad Towards the use of surrogate modeling in model parameter calibration in injection molding process simulation Résumé : Les coûts de calcul des simulations du procédé de moulage par injection ont augmenté au cours des dernières années en raison de la complexité accrue des modèles intégrés. Ceci est particulièrement problématique pour des calculs d'optimisation de pièces ou pour les analyses de sensibilité de paramètres. Une façon de surmonter ce problème est d'implémenter des métamodèles pour réaliser ces simulations haute-fidélités. Ces métamodèles peuvent jouer un rôle important dans le cas du moulage par injection de polymères semi-cristallins pour modéliser le processus de cristallisation induit par l'écoulement. À ce jour, la plupart des logiciels commerciaux ne prennent pas explicitement en compte la cristallisation des polymères, ce qui conduit à diverses erreurs dans les prédictions de remplissage, d'évolution des champs de contraintes et dans les prédictions dimensionnelles. La faisabilité de l'utilisation des métamodèles pour identifier les paramètres d'un modèle physique est d'abord présentée. Il s'en suit la mise en oeuvre d'un modèle de cristallisation thermo-mécanique afin de décrire la cristallisation induite par l'écoulement d'un matériau thermoplastique semi-cristallin non renforcé pendant le procédé de moulage par injection. Le modèle de cristallisation est couplé aux modèles de viscosité et PVT dans le logiciel commercial Autodesk R ⃝ Moldflow R ⃝ Insight 2021 en utilisant la fonction Solver API. Les paramètres du modèle sont identifiés à l'aide d'un schéma de recalage qui utilise trois différents métamodèles représentant les résultats de pression simulés pour effectuer une optimisation multi-objectifs. Les prédictions de remplissage ainsi que les pressions calculées sont présentées en utilisant les paramètres du modèle calibré en comparaison avec celles mesurées pendant le moulage par injection réel d'une pièce en polyoxyméthylène avec différentes conditions de processus. Les résultats montrent des améliorations majeures dans les prédictions des signaux de pression ainsi que dans l'état de remplissage des pièces produites et les épaisseurs estimées de la couche de peau formée dans des conditions de cisaillement élevé. Les modèles recalés sont testés en utilisant différentes géométries de moules pour évaluer leurs performances. Mots clés : modélisation de substitution, moulage par injection, cristallisation, simulation, thermoplastique, paramètre de modèle, métamodèle.

Conclusions et Perspectives

Abstract: The computational costs of injection molding simulations have been increasing in the past years due to the higher complexity of the embedded models. This is especially problematic in case of using such simulation models for optimization routines or sensitivity analyses. One way to overcome this challenge is by having a surrogate model, also known as a metamodel, of these high-fidelity simulations, which provides a cheaper way to perform these types of analyses. These surrogate models can play an important role in the case of the injection molding of semi-crystalline polymers to model the flow-induced crystallization process. To date, most commercial software do not explicitly take polymer crystallization into account leading to various errors in the fill predictions as well as the calculation of warpage and shrinkage. This is mainly due to the common complexity of the models used to describe crystallization and the challenging respective model parameter identification process under injection molding conditions. To close this gap, in this thesis, the feasibility of using surrogate modeling to identify modeling parameters is first studied. This is then followed by the implementation of a thermo-mechanical crystallization model in order to describe the flow-induced and quiescent crystallization of an unreinforced semi-crystalline thermoplastic material during injection molding. The crystallization model is defined alongside crystallization-dependent viscosity, PVT and solidification models in the commercial software Autodesk R ⃝ Moldflow R ⃝ Insight 2021 using the Solver API feature. The model parameters are identified using a calibration scheme that employs three surrogate models representing the simulated pressure results to perform a multi-objective optimization. The fill predictions as well as the calculated pressure fields are presented using the calibrated model parameters in comparison to those measured during the actual injection molding of a polyoxymethylene part with different process conditions. The results show major improvements in the predictions of the pressure signals as well as the filling status of the produced parts and the estimated skin layer thicknesses formed under high-shear conditions. Additionally, the calibrated models are tested using various mold geometries to assess the calibrated models' performance. Keywords: surrogate modeling, injection molding, crystallization, simulation, thermoplastic, model parameter, metamodel.

1. 1 103 4. 6 105 4. 8 105 4. 9

 1103610581059 Overview of the literature in the field of modeling the flow-induced crystallization during polymer processing. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 List of the commonly used models to couple the viscosity and the crystallization using a normalized rheological function. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Overview of the literature in the field of surrogate modeling based optimization in injection molding simulation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1 Overview of the DSC measurement runs and the studied samples' weight. . . . . . . . 2.2 Overview of the material characteristics of the studied POM determined using the DSC measurements. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3 Determined shot volumes and holding times from the performed filling and packing studies. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4 The upper and lower bounds of the variable process settings used in the injection molding trials. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.5 DoE of the injection molding experiments for the full length 3 mm thick part. . . . . . 2.6 DoE of the injection molding experiments for the half length 3 mm thick part. . . . . 2.7 DoE of the injection molding experiments for the full length 1.5 mm thick part. . . . . 2.8 DoE of the injection molding experiments for the half length 1.5 mm thick part. . . . 2.9 The measured skin layer thickness obtained from analyzing the micrographs from the polarized optical microscope. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xv LIST OF TABLES 3.1 Varied parameters in the simulation runs used to generate the surrogate models along with their upper and lower bounds. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Definition of the testing simulations for evaluating the performance of the surrogate models. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 An overview of the proposed DoE and surrogate modeling techniques. . . . . . . . . . 3.4 The model and boundary condition parameters used in the default simulation and the ones obtained in the optimization use case. . . . . . . . . . . . . . . . . . . . . . . . . 4.1 The model parameters of the implemented crystallization model. . . . . . . . . . . . . 4.2 Material characteristics of the studied POM homopolymer obtained from the DSC measurements. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3 The specific heat data used in the Moldflow simulations along with self-defined models. 4.4 The model parameters of the implemented viscosity model. . . . . . . . . . . . . . . . 4.5 The model parameters of the implemented PVT model. . . . . . . . . . . . . . . . . . 4.6 Comparison between the implemented models and the ones used by default in a Moldflow simulation for the injection molding of a semi-crystalline thermoplastic material. . 4.7 The process settings used in the basis simulations for the three generated surrogate models. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.8 The surrogate model variables with their upper and lower bounds. . . . . . . . . . . . 4.9 The experimentally-determined parameters used to define the six objective functions used in the optimization routine to identify the modeling parameters. . . . . . . . . . 4.10 The identified modeling parameters from the multi-objective optimization algorithm. . 4.11 The skin layer thicknesses determined using the simulation utilizing the calibrated models across the thickness at the three sensor locations for four processing conditions (V1, V3, V19, V21). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xvi List of Figures 1 Model parameter optimization scheme using a surrogate model. . . . . . . . . . . . . . 1.1 Simple sketch of a reciprocating-screw injection molding machine. . . . . . . . . . . . . 1.2 Scheme of the hierarchical structure of an injection molded part with a semi-crystalline material. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 An illustrative example plot of the specific volume as a function of temperature at a fixed pressure using a crystallinity-dependent Tait model for five different cooling rates Ṫ = 1, 5, 10, 20, 50 • C/min (solid lines) in comparison to using the Tait model (dashed line). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Overview of the most commonly used methods during the generation of a surrogate model. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.5 Latin Hypercube sampling in a 2-dimensional design space with a set of sampling sites denoted by integers from which the shaded cells are the selected sampling sites. . . . . 1.6 Comparison of 2-dimensional design spaces with 1000 sample points obtained using (a) QMC, (b) LHS. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1 DSC thermograms using five different cooling rates. . . . . . . . . . . . . . . . . . . . 2.2 Crystallization peaks of the different DSC measurements. . . . . . . . . . . . . . . . .2.3The dependency of the crystallization onset temperature T c,onset as a function of the different cooling rates Ṫ along with two linear fits, Linear fit 1: data from Ṫ = 1, 3, 5, 10, 40• C/min, Linear fit 2: data from Ṫ = 1, 3, 5, 10 • C/min. . . . . . . . . . . . . . . . . . . xvii LIST OF FIGURES 2.4 DSC thermogram obtained for POM using a cooling rate of 10 • C/min showing the 4-points baseline determination technique used to calculate the relative crystallinity evolution. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.5 The calculated relative crystallinity α versus temperature for the five cooling rates used in the DSC measurements. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.6 The calculated relative crystallinity α versus time for the five cooling rates used in the DSC measurements. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.7 The complex viscosity measured during two dynamic temperature ramp tests as a function of temperature under constant frequency specified in rad/s and constant strain given in %. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.8 The mold produced to carry out the injection molding trials. . . . . . . . . . . . . . .2.9 Sketch of the injection-molded part including some characteristic dimensions in mm and the location of the three pressure and temperature combination sensors (p-T-sensors) as well as the insert's location. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.10 3D-view of the full length 3 mm thick part. . . . . . . . . . . . . . . . . . . . . . . . . 2.11 3D-view of the half length 3 mm thick part. . . . . . . . . . . . . . . . . . . . . . . . . 2.12 The pressure signals measured during the injection molding of the full length 3 mm thick mold cavity for two different processing conditions: (a) S3, (b) S6. . . . . . . . . 2.13 The pressure signals measured during the injection molding of the half length 3 mm thick mold cavity for four different processing conditions: (a) V1, (b) V9, (c) V15, (d) V22. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.14 The pressure signals measured during the injection molding of the 1.5 mm thick mold cavity: (a) using the full length geometry for condition W3, (b) using the half length geometry for condition T20. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.15 Short shoted parts produced during the injection molding of the half length 3 mm thick mold cavity under two processing conditions: (a) V1, (b) V11. . . . . . . . . . . . . . xviii LIST OF FIGURES 2.16 The location of the microtome cuts used to observe the crystal morphology using polirized optical microscopy in the half length 3 mm thick parts: (a) in the case of a short shoted part, (b) in the case of a fully filled part. . . . . . . . . . . . . . . . . . . . . . 2.17 Micrographs of the microtome cut made at the second pressure sensor located after the gate (Cut 2) for the V1 processing condition under three different magnifications. . . . 2.18 Micrographs of Cut 1 using the same magnification of 25x for the four studied processing conditions: (a) V1, (b) V3, (c) V19, (d) V21. . . . . . . . . . . . . . . . . . . . . . . . 2.19 Micrographs of Cut 2 using the same magnification of 25x for the four studied processing conditions: (a) V1, (b) V3, (c) V19, (d) V21. . . . . . . . . . . . . . . . . . . . . . . . 2.20 Micrographs of Cut 3 using the same magnification of 25x for the four studied processing conditions: (a) V1, (b) V3, (c) V19, (d) V21. . . . . . . . . . . . . . . . . . . . . . . . 3.1 The meshed simulation model including the part (dark green), runner and sprue (light green), cooling channels (blue) and feed system (red). . . . . . . . . . . . . . . . . . . 3.2 Workflow for generating an adaptive surrogate model using the pyUQ active learning algorithm. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Predicted pressure signals by the three proposed surrogate models after training with 30 simulations in comparison to the true high-fidelity simulation results for four test cases. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 RMSE of the three studied surrogate modeling methods after training with 30 simulation runs. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.5 Predicted pressure signals by the three proposed surrogate models after training with 60 simulations in comparison to the true high-fidelity simulation results for four test cases. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.6 RMSE of the three studied surrogate modeling methods after training with 60 simulation runs. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.7 Predicted pressure signals by the three proposed surrogate models after training with 120 simulations in comparison to the true high-fidelity simulation results for four test cases. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xix LIST OF FIGURES 3.8 RMSE of the three studied surrogate modeling methods after training with 120 simulation runs. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.9 Average normalized RMSE over the five test runs for each surrogate model generated as a function of the number of training simulations used. . . . . . . . . . . . . . . . . . 3.10 The first order Sobol sensitivity indices obtained using the generated adaptive-PCE-LS surrogate model trained with 120 simulation runs. . . . . . . . . . . . . . . . . . . . . 3.11 The second order Sobol sensitivity indices obtained using the generated adaptive-PCE-LS surrogate model trained with 120 simulation runs. . . . . . . . . . . . . . . . . . . 3.12 Comparison between the experimental pressure signals and those obtained by the simulation using the default HT C f illing , HT C packing and D 3 values (Simulation: Default) and the optimized ones (Simulation: Optimization) for four example DoEs. . . . . . . 3.13 RMSE of the simulated pressure results obtained using the default HT C f illing , HT C packing and D 3 values (Simulation: Default) and the optimized ones (Simulation: Optimization) for the 27 experimental DoEs. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1 Fitting surface relating the different pressure dependencies of the melting temperature (b 6 ), equilibrium melting temperature (a) and glass transition temperature (b). . . . . 4.2 Flowchart representing the link between the Moldflow solver (turquoise nodes) and the implemented C++ code used to define the user models (gray nodes), DLL=Dynamiclink library. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3 The POD basis functions of the three generated surrogate models. . . . . . . . . . . . 4.4 Comparison between the pressure predictions using the SM1 surrogate model and those obtained by four example high-fidelity simulations used for the model's (a) training and (b) testing. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.5 The error metrics of the calculated POD parameters Γ and time shifts t start and t end using the SM1 surrogate model when using input data from the training and testing sets.Comparison between the pressure predictions using the SM2 surrogate model and those obtained by four example high-fidelity simulations used for the model's (a) training and (b) testing. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xx LIST OF FIGURES 4.7 The error metrics of the calculated POD parameters Γ and time shifts t start and t end using the SM2 surrogate model when using input data from the training and testing sets.Comparison between the pressure predictions using the SM3 surrogate model and those obtained by four example high-fidelity simulations used for the model's (a) training and (b) testing. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . The error metrics of the calculated POD parameters Γ and time shifts t start and t end using the SM3 surrogate model when using input data from the training and testing sets.106 4.10 The norm of the residuals as a function of the number of iterations performed by lsqnonlin during the optimization routine. . . . . . . . . . . . . . . . . . . . . . . . . . 107 4.11 The pressure signals for the V6 processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for this condition. . . . . . . . . . . . . . . . . . . 109 4.12 The pressure signals for the V14 processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for this condition. . . . . . . . . . . . . . . . . . . 109 4.13 The pressure signals for the V27 processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for this condition. . . . . . . . . . . . . . . . . . . 110 4.14 The RMSE of the pressure predictions using the simulation with the calibrated models (Optimization) and the one using the default Moldflow models (Default) for each of the three processing conditions used to generate a surrogate model. . . . . . . . . . . . . . 110 4.15 The averaged normalized RMSE of the pressure predictions at P2 using the simulation with the calibrated models (Optimization) and the one using the default Moldflow models (Default) for all 27 experimental processing conditions for the half length 3 mm thick geometry. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 xxi LIST OF FIGURES 4.16 The P2 pressure signals for three processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for these conditions. . . . . . . . . . . . 4.17 The experimental fill results (a) in comparison to the ones predicted using the (b) default simulation, (b) simulation with the calibrated models. . . . . . . . . . . . . . . 4.18 The comparison between the skin layer thicknesses determined experimentally and those obtained using the optimized simulation along the part thickness of the three sensor locations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.19 Two example P1 and P3 pressure predictions for two processing conditions T2 and T25 injected with the half length 1.5 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental). . . . 4.20 Two example P1 and P3 pressure predictions for two processing conditions S1 and S5 injected with the full length 3 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental). . . . 4.21 Two example P1 and P3 pressure predictions for two processing conditions W2 and W10 injected with the full length 1.5 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental). . . . 5.1 Schéma d'optimisation des paramètres du modèle utilisant un modèle de substitution.

  of a crystallization model implemented in an injection molding simulation. In order to reach this INTRODUCTION goal, a feasibility study is first performed to assess the use of different surrogate modeling techniques and design of experiments in calibrating modeling parameters in the injection molding simulation of a polyoxymethylene thermoplastic part. This is followed by the implementation of a modified Poitou et al.'s crystallization model along with crystallization-dependent viscosity, PVT and solidification models in the commercial software Autodesk R ⃝ Moldflow R ⃝ Insight 2021 using the Solver API feature.

Content 1 . 1

 11 Injection Molding Process . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Injection Molding Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.1 Main Approximations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.2 Governing Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Polymer Crystallization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.1 Semi-crystalline Thermoplastics . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.2 Crystallization from the Melt . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.3 Experimental Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.4 Modeling Approaches . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.5 Effect on Other Material Properties . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Surrogate Modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.1 General Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.2 Design of Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.3 Model Order Reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.4 Surrogate Modeling Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.5 Fitting Techniques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.6 Some Applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.5 Surrogate Modeling in Injection Molding Simulation . . . . . . . . . . . . . 1.6 Summary and Thesis Structure . . . . . . . . . . . . . . . . . . . . . . . . . .

1. 2 .Figure 1 . 1 -

 211 Figure 1.1 -Simple sketch of a reciprocating-screw injection molding machine, based on [18].

Figure 1 . 2 -

 12 Figure 1.2 -Scheme of the hierarchical structure of an injection molded part with a semi-crystalline material, based on [1].

( b 1 α + b 2 α b 3 )

 13 EmpiricalTitomanlio et al.[START_REF] Titomanlio | On the Simulation of Thermoplastic Injection Moulding Process: II Relevance of Interaction Between Flow and Crystallization[END_REF], Guo and Narh[START_REF] Guo | Computer simulation of stress-induced crystallization in injection molded thermoplastics[END_REF]

Figure 1 .

 1 Figure 1.3 using a self-produced plot by showing the difference in the specific volume predictions using a cooling rate or crystallization dependent pressure-volume-temperature (PVT) model in comparison to the widely used Tait model. As the cooling rate changes between 1 and 50 • C/min, the modified Tait model predicts different transition zone depending on the cooling rate whereas the Tait model estimates the same specific volume profile for all cooling rates.

2 . 3 . 4 . 6 .Figure 1 . 4 -

 234614 Figure 1.4 -Overview of the most commonly used methods during the generation of a surrogate model [6, 83].

Figure 1 .

 1 Figure 1.4 provides a more detailed look into the most commonly used techniques to perform the second and third steps in the surrogate model generation process presented above. After defining the design space, an additional optional step is to perform a model order reduction on the output result(s)

1. 4 .Figure 1 . 5 -

 415 Figure 1.5 -Latin Hypercube sampling in a 2-dimensional design space Ω = [0, 1] 2 with a set of sampling sites denoted by integers from which the shaded cells are the selected sampling sites.

Figure 1 . 6 -

 16 Figure 1.6 -Comparison of 2-dimensional design spaces Ω = [0, 1] 2 with 1000 sample points obtained using (a) QMC, (b) LHS.

1. 5 .

 5 SURROGATE MODELING IN INJECTION MOLDING SIMULATION their injection molding simulation. This was done by generating an ANN-based surrogate model representing the fiber orientation evolution result while varying the two model parameters using a full factorial design space.

Content 2 .

 2 1. MATERIAL DESCRIPTION

3 .

 3 Ramp 10 • C/min to 220 • C. 4. Isothermal for 2 min. 5. Ramp x • C/min to 25 • C with x = 1, 3, 5, 10, 40 • C/min.

Figure 2 .

 2 Figure 2.1 shows the DSC thermograms that were obtained for the five different cooling rates. The main focus of the performed DSC experiments is to study the crystallization of the polymer and determine the evolution of the crystallinity as a function of time and temperature. Therefore, Figure 2.2 presents the crystallization peaks for the different cooling rates evaluated. It can be seen that as the cooling rate increases the peak and area underneath it increase while shifting to lower temperatures. From the peaks shown in Figure 2.2, it is possible to determine the crystallization temperature (T c ), the crystallization onset temperature (T c,onset ) and the crystallization enthalpy (∆H c ). Additionally,

Figure 2 . 1 -

 21 Figure 2.1 -DSC thermograms using five different cooling rates.

Figure 2 . 2 -

 22 Figure 2.2 -Crystallization peaks of the different DSC measurements.

Figure 2 . 3 -3

 23 Figure 2.3 -The dependency of the crystallization onset temperature T c,onset as a function of the different cooling rates Ṫ along with two linear fits, Linear fit 1: data from Ṫ = 1, 3, 5, 10, 40 • C/min, Linear fit 2: data from Ṫ = 1, 3, 5, 10 • C/min.

Figure 2 . 4 -

 24 Figure 2.4 -DSC thermogram obtained for POM using a cooling rate of 10 • C/min showing the 4-points baseline determination technique used to calculate the relative crystallinity evolution.

Figure 2 .

 2 Figure 2.7 presents the results of the dynamic temperature ramp tests described in the previous section. By observing the two curves, it is apparent that as the samples cool down from 210 • C until the onset of crystallization their viscosity gradually increases. Whereas, after the start of crystallization, a singularity appears and the viscosity increases more than two orders of magnitudes between 152

Figure 2 . 8 -

 28 Figure 2.8 -The mold produced to carry out the injection molding trials.

Figure 2 . 9 -

 29 Figure 2.9 -Sketch of the injection-molded part including some characteristic dimensions in mm and the location of the three pressure and temperature combination sensors (p-T-sensors) as well as the insert's location.

Figure 2 .

 2 Figure 2.11 -3D-view of the half length 3 mm thick part.

and 2 .

 2 11, showing the full length geometry and the half length one, respectively.2.4.2 Filling and Packing StudiesBefore starting the official experimental runs, filling and packing studies are performed on the four different geometries of the produced mold. A filling study is done to determine the shot volume for each geometry by trying to fill the cavity up to 99% without applying any holding pressure. Whereas a packing study, it is used to find the appropriate holding time by specifying a suitable packing pressure and varying the holding time until the weight of the part becomes constant indicating the freezing of the gate. For this work, the studies are made on an electrical injection molding machine (ENGEL E-Motion 440/220 T) with a material melt temperature of 220 • C. For the filling studies, a mold temperature of

Figures 2 .

 2 Figures 2.12a and 2.12b show two example pressure signals measured during the injection of the full length 3 mm thick cavity under the conditions S3 and S6 (refer to Table 2.5), respectively. It is

Figure 2 .Figure 2 .Figure 2 .Figure 2 . 12 -

 222212 Figure 2.13 shows four representative results of the pressure signals observed during the injection of the half length 3 mm thick part.Figure 2.13a presents the pressure signals of the V1 run (refer to Table2.6) for which the part shorted e.g. not filled fully. The resultant parts are presented in Figure2.15a as well as another example of a short shot produced using the V11 process settings in Figure2.15b. These short shots are due to the early freezing of the gate which takes place after the switch-over. The low holding pressure of 20 MPa which is applied after the switch-over to finish the filling of the part is not sufficient to produce a complete part. Figure2.13b shows the pressure signals for the V9 run which produced a full part. An interesting phenomenon is observed in the presented V9

Figure 2 . 13 -

 213 Figure 2.13 -The pressure signals measured during the injection molding of the half length 3 mm thick mold cavity for four different processing conditions: (a) V1, (b) V9, (c) V15, (d) V22.

Figure 2 . 14 -

 214 Figure 2.14 -The pressure signals measured during the injection molding of the 1.5 mm thick mold cavity: (a) using the full length geometry for condition W3, (b) using the half length geometry for condition T20.

Figure 2 . 15 -

 215 Figure 2.15 -Short shoted parts produced during the injection molding of the c mold cavity under two processing conditions: (a) V1, (b) V11.

Figure 2 . 16 -

 216 Figure 2.16 -The location of the microtome cuts used to observe the crystal morphology using polarized optical microscopy in the half length 3 mm thick parts: (a) in the case of a short shoted part, (b) in the case of a fully filled part.

Figure 2 . 17 -

 217 Figure 2.17 -Micrographs of the microtome cut made at the second pressure sensor located after the gate (Cut 2) for the V1 processing condition under three different magnifications.

Figure 2 . 18 -

 218 Figure 2.18 -Micrographs of Cut 1 using the same magnification of 25x for the four studied processing conditions: (a) V1, (b) V3, (c) V19, (d) V21.

Figure 2 . 19 -

 219 Figure 2.19 -Micrographs of Cut 2 using the same magnification of 25x for the four studied processing conditions: (a) V1, (b) V3, (c) V19, (d) V21.

Figure 2 . 20 -

 220 Figure 2.20 -Micrographs of Cut 3 using the same magnification of 25x for the four studied processing conditions: (a) V1, (b) V3, (c) V19, (d) V21.

Chapter 3 Feasibility

 3 Study: Use of Surrogate Modeling for Modeling Parameter Calibration Content 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Software Tools . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.1 Injection Molding Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.2 Python Uncertainty Quantification Library . . . . . . . . . . . . . . . . . . . . 3.2.3 MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Simulation Environment . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 Proposed DoE and Surrogate Modeling Techniques . . . . . . . . . . . . . . 3.4.1 POD-NLR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4.2 Regular-PCE-LS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4.3 Adaptive-PCE-LS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.5 Evaluation Strategy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.6 Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.7 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.7.1 Comparison of Surrogate Models . . . . . . . . . . . . . . . . . . . . . . . . . . 3.7.2 Sensitivity Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.7.3 Modeling Parameter Calibration . . . . . . . . . . . . . . . . . . . . . . . . . . 3.8 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

  1. a non-linear regression model of proper orthogonal decomposition coefficients trained using a fixed Latin Hypercube sampling (POD-NLR), 2. a polynomial chaos expansion model fitted with a least-squares optimization technique using the same fixed Latin Hypercube sampling points (Regular-PCE-LS), and 3. a polynomial chaos expansion model fitted with a least-squares optimization technique but using an adaptive DoE that is enriched while generating the model (Adaptive-PCE-LS).

  1 using the respective POM material card given in the Moldflow R ⃝ database. The simulations are recreating the injection molding experiments for the half length 3 mm thick part presented in Section 2.4 along with example pressure signals in Figure 2.13. The simulation consists of a cooling, fill and pack analyses. Therefore, the simulation model included cooling channels meshed as beam elements and the part as well as the feed system (nozzle & flange) meshed using tetrahedral elements with 24 layers through the thickness. The feed system is defined as hot runner in the simulation model. Lateral and top views from the meshed model are shown in Figure 3.1a and Figure 3.1b, respectively. The process settings of the base simulation are defined according to imported data from the ENGEL sim link software tool. These include the filling and packing profiles along with the switch-over ram position and the machine settings. A single simulation requires around 75 minutes to be completed on a workstation with a 4.10 GHz processor and 32 GB RAM.

3. 3 .

 3 SIMULATION ENVIRONMENT (a) Lateral view (b) Top view

Figure 3 . 1 -

 31 Figure 3.1 -The meshed simulation model including the part (dark green), runner and sprue (light green), cooling channels (blue) and feed system (red).

2 .

 2 Polynomial chaos expansion: since the surrogate model input variables have a uniform distribution, after applying the Gram-Schmidt orthogonalization to equation 1.43, the Legendre polynomial family is obtained and thus used as basis functions to represent the pressure signal results from the simulation. A hyperbolic truncation of 1 is utilized to generate the model which 3.4. PROPOSED DOE AND SURROGATE MODELING TECHNIQUES

3 .

 3 Model fitting: the algorithm uses the pressure time-dependent training data to solve an optimization problem which minimizes the squared residuals between the simulation's pressure signals and the predicted output of the PCE model. The solution of the least-squares optimization problem is the vector containing the polynomial coefficients.

Figure 3 . 2 -

 32 Figure 3.2 -Workflow for generating an adaptive surrogate model using the pyUQ active learning algorithm.

  6 are utilized. The optimization problem aims at minimizing the difference between the experimental results Y exp and the surrogate model predictions Ŷ by varying the respective uncertain input parameters of the surrogate model. It can be represented as follows: find V ar = (HT C f illing , HT C packing , D 3 ), -Ŷi (X)∥, with V = (V inj , T c,in , P hold ), X = (V inj , T c,in , P hold , HT C f illing , HT C packing , D 3 ), s.t. 5000 ≤ HT C f illing ≤ 10000 (W/m 2• C), 1000 ≤ HT C packing ≤ 5000 (W/m 2• C), 0 ≤ D 3 ≤ 0.40 (K/MPa).

Figure 3 . 3 -Figure 3 . 4 -

 3334 Figure 3.3 -Predicted pressure signals by the three proposed surrogate models after training with 30 simulations in comparison to the true high-fidelity simulation results for four test cases.

Figure 3 .

 3 Figure 3.3 presents the models' predictions obtained after using 30 training simulation runs in comparison to the true high-fidelity simulation results. At first glance, it seems that all three surrogate models are able to capture the main features of the pressure signal after only 30 training simulations. However, one can distinguish some particular discrepancies such as the inability of the POD-based model to predict the t start for Test 2 and error of the adaptive-PCE for estimating the pressure at the end of packing in Tests 3 and 4. Additionally, there are deviations by predicting the pressure peak at the end of the filling phase for all surrogate models. A quantitative analysis of the performance of the different surrogate models is given in Figure 3.4, where the normalized RMSE is plotted for the four test simulations. The average normalized RMSEs of 0.084, 0.093 and 0.119 for POD-NLR, regular-PCE-LS and adaptive-PCE-LS, respectively, show that the POD-based surrogate model is slightly more accurate in recreating the simulation results after training with 30 high-fidelity simulations than the PCE-based models.Looking to increase the prediction capability, the number of simulations used to train the surrogate models is doubled to 60. Overall, this increase in training simulations led to noticeable improvements in the pressure signal predictions of the test runs as shown in Figure3.5. For example, by comparing Figures 3.3b and 3.5b, the POD-based model is now able to better predict the t start for Test 2 as well as both PCE-based models can estimate the pressure result more accurately during the packing phase. In general, there is still room for improvements when estimating the pressure for Test 4 as well as predicting the pressure peak at the end of filling and the time at which the pressure goes to zero

Figure 3 . 5 -Figure 3 . 6 -

 3536 Figure 3.5 -Predicted pressure signals by the three proposed surrogate models after training with 60 simulations in comparison to the true high-fidelity simulation results for four test cases.

Figures 3 .

 3 10 and 3.11 show the resultant first order and some of the second order sensitivity indices, respectively.The first order indices underline the contribution of a certain surrogate model parameter to the variance of the pressure output as explained in Section 1.4.6.1. The sensitivity of the different process settings shown in Figure3.10a are in agreement with the common understanding of the injection molding process. The injection velocity is the more sensitive factor during the filling stage because

Figure 3 . 7 -Figure 3 . 8 -

 3738 Figure 3.7 -Predicted pressure signals by the three proposed surrogate models after training with 120 simulations in comparison to the true high-fidelity simulation results for four test cases.

Figure 3 . 9 -Figure 3 . 10 -

 39310 Figure 3.9 -Average normalized RMSE over the five test runs for each surrogate model generated as a function of the number of training simulations used.

Figure 3 . 3 Figure 3 . 11 -

 33311 Figure 3.11 presents some of the second order sensitivity indices representing the most importantinteractions between the surrogate model parameters as a function of time in relation to their combined contribution to the output pressure. During the filling phase, the main interaction is given by the injection velocity and D 3 . As for the beginning of the packing phase, the most important interactions are given by the holding pressure along with both the injection velocity and D 3 . Whereas, along the packing phase, it turns out that there is not significant two-parameter relations affecting the

Figure 3 . 9 .D 3 .

 393 Figure 3.9. The POD-based surrogate model is used for the parameter optimization as the coupling of this particular model implementation with the optimization algorithm presented in Section 3.6 was straightforward in comparison with the other surrogate model implementations. The goal of this calibration is to identify the optimal uncertain simulation model parameters: HT C f illing , HT C packing

Figure 3 . 12 -Figure 3 . 13 -

 312313 Figure 3.12 -Comparison between the experimental pressure signals and those obtained by the simulation using the default HT C f illing , HT C packing and D 3 values (Simulation: Default) and the optimized ones (Simulation: Optimization) for four example DoEs.

3. 8 .

 8 CONCLUSION training data. By increasing the number of training data, the adaptive-PCE-LS technique improves steadily the accuracy of the surrogate model and reaches a comparable performance with the fixed-DoE techniques for 120 training simulations. The advantages of the adaptive-PCE-LS should appear more evident when having a high number of variables.Based on a Sobol sensitivity analysis using the generated surrogate model, we found that the heat transfer coefficient during filling has a low impact on the pressure evolution at a sensor node in front of the heat transfer coefficient during packing or the pressure-dependent viscosity coefficient.By optimizing the HT C f illing , HT C packing and the pressure-dependent viscosity model parameter D 3 using measured pressure signals, it was feasible to calibrate the high-fidelity simulation in order to decrease the deviations in the pressure estimation when compared to those obtained using default simulation values. However, there is still a gap for getting an accurate estimation of the pressure fields in injection molding, which could be closed with this methodology by including additional simulation model parameters in addition to a refinement of the models implemented in the simulation as shown in the next chapter by implementing a flow-induced crystallization model along with crystallizationdependent PVT, viscosity and solidification criterion for semi-crystalline thermoplastic polymers.4.1. INTRODUCTION

Figure 4 . 1 -

 41 Figure 4.1 -Fitting surface relating the different pressure dependencies of the melting temperature (b 6 ), equilibrium melting temperature (a) and glass transition temperature (b).

  equations 1.29 to 1.31. The PVT model parameters are taken from Moldflow R ⃝ 's material database for POM and are 4.2. DESCRIPTION OF THE IMPLEMENTED MODELS Table 4.5 -The model parameters of the implemented PVT model. Model Parameter Value (units) b 1m 8.452 × 10 -4 (m 3 /kg) b 2m 5.65 × 10 -7 (m 3 /kgK) × 10 -4 (m 3 /kg) b 2s 2.944 × 10 -7 (m 3 /kgK) b 3s 2.876 63 × 10 8 (Pa) b 4s 4.784 × 10 -3 (1/K)

4. 3 .Figure 4 . 2 -

 342 Figure 4.2 -Flowchart representing the link between the Moldflow solver (turquoise nodes) and the implemented C++ code used to define the user models (gray nodes), DLL=Dynamic-link library.

Figure 4 . 5 -

 45 Figure 4.5 -The error metrics of the calculated POD parameters Γ and time shifts t start and t end using the SM1 surrogate model when using input data from the training and testing sets.

Figure 4 .

 4 4b shows the pressure signals of four simulations used to test the model in comparison to the pressure results obtained by the fitted model. It is apparent that the model is able to recreate the form of the pressure signal.

Figure 4 . 5 .

 45 Figure 4.5. Figures 4.5a and 4.5b show the RMSE and R 2 score, respectively, of these parameters obtained during the training and testing steps of the SM1 surrogate model. The parameter errors for the training set are, as expected, lower than those obtained for the testing set. Using unseen input data on the trained SM1 model doubles the RMSE of the POD parameter Γ 1 and Γ 4 as well as the time shift t start . Overall, the surrogate model's parameters are predicted with a RMSE lower than 26 % which is considered for these purposes as acceptable. By analyzing the coefficient of determination results, the calculated Γ 3 and Γ 4 parameters by the fitted surrogate model produce a negative R 2 score

Figure 4 . 7 -

 47 Figure 4.7 -The error metrics of the calculated POD parameters Γ and time shifts t start and t end using the SM2 surrogate model when using input data from the training and testing sets.

Figure 4 . 8 -

 48 Figure 4.8 -Comparison between the pressure predictions using the SM3 surrogate model and those obtained by four example high-fidelity simulations used for the model's (a) training and (b) testing.

Figure 4 . 9 -

 49 Figure4.9 -The error metrics of the calculated POD parameters Γ and time shifts t start and t end using the SM3 surrogate model when using input data from the training and testing sets.

Section 3 . 6 .

 36 The optimization routine necessitated 15 iterations to reach a local minimum with a residuals norm of 11.72 as shown in Figure4.10. In order to check if there exist other local minima for this multi-objective optimization and to find the global minimum, the built-in MultiStart MATLAB function is used to find the possible other local solutions to the optimization problem by starting from various initial points. This resulted in finding the same combination of parameters suggesting that the found local minimum is in fact the global one. The final optimized parameters are given in Table

Figure 4 . 10 - 107 4. 6 .

 4101076 Figure 4.10 -The norm of the residuals as a function of the number of iterations performed by lsqnonlin during the optimization routine.

  10. Therefore, to assess the proposed implementation's results, the pressure signals obtained using the simulation with the calibrated models (Simulation: Optimization) are plotted in comparison to the ones obtained by the default Moldflow R ⃝ simulation (Simulation: Default) and as a reference the experimental results (Experimental). The plots presented in Figures 4.11

, 4 .

 4 12 and 4.13 correspond to the previously-mentioned results for the three processing conditions V6, V14 and V27, respectively and they show these predictions at the three sensor locations P1, P2 and P3 (refer to Figure2.9).

Figure 4 . 11 -

 411 Figure 4.11 -The pressure signals for the V6 processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for this condition.

Figure 4 . 12 -

 412 Figure 4.12 -The pressure signals for the V14 processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for this condition.

Figure 4 . 13 -

 413 Figure 4.13 -The pressure signals for the V27 processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for this condition.

  Figures 4.14a

  , 4.14b and 4.14c present the obtained results for the three processing conditions V6, V14 and V27, respectively.

Figure 4 .

 4 Figure4.15 -The averaged normalized RMSE of the pressure predictions at P2 using the simulation with the calibrated models (Optimization) and the one using the default Moldflow models (Default) for all 27 experimental processing conditions for the half length 3 mm thick geometry.

Figure 4 .

 4 Figure 4.15 presents the averaged normalized RMSE at P2 for all experimental processing conditions (Experimental DoE) summarized in Table 2.6 while using Moldflow R ⃝ with default models or with the calibrated models of this work. By observing Figure 4.15, for the majority of the different DoEs, the prediction error is decreased by more than half from the default results. The only exceptions are seen for the experimental DoEs V7, V8 and V9. To analyze the probable reasons behind this discrepancy, Figure 4.16 presents the pressure

Figure 4 .

 4 Figure 4.16 -The P2 pressure signals for three processing condition obtained using the simulation with the calibrated models (Simulation: Optimization) in comparison to the ones obtained by the default Moldflow simulation (Simulation: Default) along with the corresponding experimental results (Experimental) for these conditions.

  (a) Experimental (b) Simulation: Default (c) Simulation: Optimization

Figure 4 . 17 -

 417 Figure 4.17 -The experimental fill results (a) in comparison to the ones predicted using the (b) default simulation, (b) simulation with the calibrated models.

Figure 4 . 18 -

 418 Figure 4.18 -The comparison between the skin layer thicknesses determined experimentally and those obtained using the optimized simulation along the part thickness of the three sensor locations.

Figure 4 . 19 -

 419 Figure 4.19 -Two example P1 and P3 pressure predictions for two processing conditions T2 and T25 injected with the half length 1.5 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental).

Figure 4 . 20 -

 420 Figure 4.20 -Two example P1 and P3 pressure predictions for two processing conditions S1 and S5injected with the full length 3 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental).

Figure 4 . 21 -

 421 Figure 4.21 -Two example P1 and P3 pressure predictions for two processing conditions W2 and W10 injected with the full length 1.5 mm thick geometry using the simulation with the calibrated models (Simulation: Optimization) and the one with the default models (Simulation: Default) in comparison to the experimental results (Experimental).

4. 7 .

 7 CONCLUSIONresults for the half length 3 mm thick part.

Figure 5 . 1 -

 51 Figure 5.1 -Schéma d'optimisation des paramètres du modèle utilisant un modèle de substitution [4].

Enfin, dans le chapitre 4 5 . 3

 453 , l'implémentation d'un modèle de cristallisation basé sur la thermomécanique est décrite ainsi que les modèles de viscosité, de PVT et de solidification dépendant de la cristallisation dans Moldflow R ⃝ pour simuler le moulage par injection d'un matériau thermoplastique non renforcé. Cinq paramètres de modélisation sont identifiés à l'aide de trois modèles de substitution générés par une routine d'optimisation multi-objectifs. Les résultats de la simulation utilisant les modèles calibrés sont ensuite analysés pour différentes conditions de traitement et géométries de 5.3. CARACT ÉRISATION DES MAT ÉRIAUX ET ESSAIS DE MOULAGE PAR INJECTION moule. Caractérisation des Matériaux et Essais de Moulage par Injection Le matériau polymère utilisé dans ce travail est un homopolymère de polyoxyméthylène (POM) non renforcé de qualité industrielle sous sa forme granulaire. Ce polymère, également connu sous le nom de polyacétal, est un matériau thermoplastique semi-cristallin d'ingénierie largement utilisé dans l'industrie automobile en raison de sa grande rigidité, de sa faible friction et de ses excellentes propriétés anti-usure [21, 113]. Deux types d'expériences sont réalisées afin de caractériser le matériau utilisé. Le premier type était les mesures DSC qui ont fourni des informations sur le comportement de cristallisation du matériau sous différentes vitesses de refroidissement. Ceci a permis de déterminer la température de fusion du matériau à 178,88 • C et l'enthalpie de cristallisation de 155,35 J/g. L'évolution de la cristallinité relative avec la température et le temps sous les différents taux de refroidissement a également été déterminée et sera utilisée pour quantifier la cinétique de cristallisation du POM dans le chapitre 4. Le deuxième type d'expérience est un test dynamique de rampe de température utilisé pour identifier la dépendance de la viscosité à la température et observer son comportement après le début de la cristallisation. Ces données seront également utilisées pour déterminer les paramètres du modèle de viscosité plus tard dans le travail. Un moule a été conçu pour permettre la production de quatre géométries différentes variant en longueur et/ou en épaisseur, ainsi que la mesure des signaux de pression pendant le processus de moulage par injection à quatre endroits de la cavité. Pour mettre en place les essais de moulage par injection, des études de remplissage et d'emballage ont été réalisées, suivies des essais officiels où, pour chaque géométrie, trois paramètres du processus sont modifiés. Les paramètres variables sont la température du moule, la vitesse d'injection et la pression de maintien. Certains signaux de pression caractéristiques ont été analysés dans ce chapitre pour donner un aperçu des résultats obtenus qui sont utilisés dans ce travail pour la calibration des paramètres du modèle. Certains des échantillons moulés par injection ont été étudiés à l'aide d'un microscope optique polarisé afin de mieux comprendre les morphologies des cristaux et les différentes couches qui se 5.4. ÉTUDE DE FAISABILIT É : UTILISATION DE LA MOD ÉLISATION DE SUBSTITUTION POUR LA CALIBRATION DES PARAM ÈTRES DE MOD ÉLISATION forment au cours du processus. Des coupes au microtome ont été réalisées à l'emplacement des capteurs pour différents paramètres du processus afin d'étudier l'effet de la température du moule et de la vitesse d'injection sur la morphologie et l'épaisseur de la couche de peau. En observant les différentes micrographies, il est apparu que plus la température du moule et la vitesse d'injection sont élevées, plus la couche cristalline est fine.5.4 Étude de Faisabilité : Utilisation de la Modélisation de Substitution pour la Calibration des Paramètres de Modélisation La faisabilité de l'utilisation de la modélisation de substitution pour calibrer les paramètres dans une simulation de moulage par injection est étudiée. Trois paramètres de modélisation sont choisis pour la calibration dans cette étude. Il s'agit des coefficients de transfert de chaleur pendant le remplissage et le compactage ainsi que du coefficient de dépendance à la pression dans le modèle de viscosité Cross-WLF [117,118]. Les paramètres susmentionnés sont connus pour affecter les résultats de pression dans la simulation et leur détermination expérimentale est difficile et prend du temps [119-121]. Pour cette calibration, la géométrie de la demi-longueur de 3 mm d'épaisseur présentée sur la figure 2.11 est utilisée et les trois paramètres de modélisation sont étalonnés à l'aide des données expérimentales du capteur de pression P2 (voir la figure 2.9). Dans cette étude de faisabilité, on étudie la précision de trois techniques de modélisation de substitution basées sur la régression et de deux méthodes de plan d'expérience dans leur capacité à approximer le signal de pression à un seul noeud de la simulation. Les approches de modélisation proposées sont les suivantes : 1. un modèle de régression non linéaire des coefficients de la décomposition orthogonale propre entraîné à l'aide d'un échantillonnage fixe de l'Hypercube latin (POD-NLR), 2. un modèle de polynomial chaos expension ajusté par une technique d'optimisation des moindres carrés utilisant les mêmes points d'échantillonnage fixes de l'Hypercube latin (Regular-PCE-LS), et 3. un modèle de polynomial chaos expension ajusté par une technique d'optimisation des moindres carrés, mais utilisant une DoE adaptative qui est enrichie pendant la génération du modèle (Adaptive-PCE-LS).

5. 5 .

 5 UTILISATION DE LA MOD ÉLISATION DE SUBSTITUTION POUR CALIBRER UN MOD ÈLE DE CRISTALLISATION INDUITE PAR L' ÉCOULEMENT DANS MOLDFLOW PVT, de la viscosité et le critère de solidification dépendant de la cristallisation pour les polymères thermoplastiques semi-cristallins. 5.5 Utilisation de la Modélisation de Substitution pour Calibrer un Modèle de Cristallisation Induite par l' Écoulement dans Moldflow L'étude de faisabilité a montré qu'il est possible d'optimiser ces paramètres en utilisant les résultats de pression expérimentale à l'aide d'un modèle de substitution représentant les résultats de pression simulée à un noeud de surface. Par conséquent, dans ce chapitre, une méthodologie similaire est suivie afin de calibrer un modèle de cristallisation induit par l'écoulement ainsi que des modèles PVT, de la viscosité et de solidification dépendant de la cristallisation. Ces modèles sont implémentés dans Moldflow R ⃝ à l'aide d'un code utilisateur. Au total, trois modèles de substitution sont utilisés pour identifier les cinq paramètres de modélisation incertains et ces substituts sont générés à l'aide d'une méthode POD-NLR analogue à celle présentée au chapitre 3. Par exemple, les signaux de pression sont d'abord réduits à l'aide de la technique de décomposition orthogonale appropriée, suivie d'une régression des moindres carrés d'un polynôme du second ordre. Cette méthodologie a été choisie car elle a montré une bonne performance en utilisant un faible nombre de données d'apprentissage et un faible nombre de paramètres d'entrée. Après avoir obtenu les paramètres de modélisation optimisés, les performances de la simulation sont évaluées en étudiant la capacité à prédire les signaux de pression, l'état de remplissage et l'épaisseur de la couche de peau à l'aide des modèles définis par l'utilisateur par rapport aux résultats expérimentaux. En outre, les modèles calibrés sont testés en utilisant trois géométries de moule différentes afin de vérifier la généralité des paramètres identifiés. Après avoir étudié les performances des modèles de substitution générés, il est apparu que la méthode d'ajustement utilisant le même polynôme du second ordre pour tous les paramètres POD et les décalages temporels n'est pas optimale. Cependant, les modèles de substitution générés étaient toujours capables de prédire les signaux de pression avec une grande précision. Par conséquent, un schéma d'optimisation multi-objectif est exécuté pour identifier les cinq paramètres incertains du modèle et les modèles calibrés sont ensuite testés dans la simulation en comparaison avec les résultats de simulation par défaut de Moldflow R ⃝ pour la demi-longueur de la pièce de 3 mm d'épaisseur.

  Dans cette thèse, après avoir prouvé la faisabilité de l'utilisation de modèles de substitution pour identifier les paramètres de modélisation par un schéma d'optimisation simple, trois modèles de substitution sont utilisés pour calibrer des modèles définis par l'utilisateur et implémentés dans une simulation de moulage par injection. Le logiciel commercial Autodesk R ⃝ Moldflow R ⃝ Insight 2021 est utilisé pour simuler le moulage par injection d'un matériau thermoplastique semi-cristallin. Dans Moldflow, quatre modèles utilisateur sont définis : un modèle de cristallisation décrivant à la fois la cristallisation induite par l'écoulement et la cristallisation naturelle, ainsi que des modèles de viscosité, de PVT et 5.6. CONCLUSIONS ET PERSPECTIVES de solidification dépendant de la cristallisation. Ces modèles sont mis en oeuvre afin de prendre en compte le processus de cristallisation qui se produit pendant le traitement du matériau et qui affecte la solution d'écoulement ainsi que les prédictions de déformation et de rétraction dans la simulation. Une optimisation multi-objectifs est ensuite réalisée pour identifier les paramètres du modèle en utilisant les modèles de substitution générés avec l'aide de signaux de pression mesurés expérimentalement. La simulation utilisant les modèles proposés avec leurs paramètres identifiés a montré des améliorations majeures dans les prédictions des signaux de pression dans diverses conditions pour la géométrie du moule de base. L'un des progrès cruciaux réalisés grâce à cette mise en oeuvre est la capacité de la simulation à prédire avec précision l'état de remplissage des pièces produites (remplissage partiel ou remplissage complet). En outre, grâce à l'inclusion d'un modèle de cristallisation, il est devenu possible d'observer l'évolution de la cristallinité relative au cours des différentes phases du processus de moulage par injection, ce qui a permis d'estimer l'épaisseur des couches cristallines formées dans des conditions de cisaillement élevé. 5.6. CONCLUSIONS ET PERSPECTIVES pourrait atténuer les erreurs probables dues à l'inclusion du seuil lors de la calibration des modèles. En outre, dans le cas de remplissages partiels, la longueur d'écoulement pourrait être incluse comme un résultat de sortie supplémentaire et utilisée dans la calibration des modèles, car elle peut donner des informations importantes sur la solidification du matériau. Dans le travail présenté, le signal de pression mesuré pendant le processus de moulage par injection est utilisé pour calibrer les modèles car il est mesuré avec précision et facilement pendant le traitement et peut être directement relié à la pression simulée. Les mesures de température, en revanche, ne sont pas aussi simples que les précédentes. Dans le meilleur des cas, les capteurs de température mesurent la température de surface du moule. Cependant, cette mesure n'est pas bien corrélée avec la température simulée. Par conséquent, il serait intéressant d'utiliser l'apprentissage machines basés sur les données pour combler l'écart entre ces deux résultats afin de pouvoir utiliser ces données mesurées dans une application similaire à celle présentée dans cette thèse. En guise de perspective sur le sujet de la modélisation de substitution dans le domaine de la simulation du moulage par injection, le présent travail peut être étendu pour inclure un modèle temporel et spatial complet. Celui-ci pourrait être utilisé pour visualiser et évaluer les solutions instantanément tout en faisant varier les paramètres d'entrée. De tels modèles de substitution fournissent la clé pour effectuer la quantification des incertitudes et des analyses de sensibilité efficaces. Une autre application intéressante est la création d'un jumeau numérique du processus de moulage par injection à l'aide de ces méthodologies.

Table 1 .

 1 2 -List of the commonly used models to couple the viscosity and the crystallization using a normalized rheological function Γ = η/η 0 .

	Author(s)	Equation	Approach
	Kitano et al. [67], Metzner [68]	Γ = 1/ (1 -α/b 1 ) -2	Suspension
	Ziabicki [69]	Γ = 1/ (1 -α/b 1 ) b 2	Empirical, based on suspension
	Zuidema et al. [29], Doufas et al.	Γ = exp
	[46]

Table 1 .

 1 3 provides an overview of some of the literature available in the field of surrogate modeling in injection molding simulation. More comprehensive reviews can be found in[START_REF] Cheng | Multiobjective optimization of injection molding parameters based on soft computing and variable complexity method[END_REF][START_REF] Gao | Intelligent methods for the process parameter determination of plastic injection molding[END_REF].1.5. SURROGATE MODELING IN INJECTION MOLDING SIMULATION

Table 1 .

 1 3 -Overview of the literature in the field of surrogate modeling based optimization in injection molding simulation.

	Design Variables		Mold and melt temperatures, injection and packing	times, packing pressure	Mold and melt temperatures, packing time, packing	pressure		Packing and cooling times, packing pressure profile	Runner diameters, packing and cooling times, packing	pressure	Mold temperature, packing time			Mold and melt temperatures, injection, packing and	cooling times, packing pressure	Mold and melt temperatures, injection, packing and	cooling times, packing pressure	Mold temperature, ram speed, packing time and pres-	sure	12 process settings		Cooling temperature, injection, packing and cooling	times, packing pressure	Mold and melt temperatures, injection, packing and	cooling times, packing pressure	Mold and melt temperatures, injection, packing and	cooling times, packing pressure, injection velocity	Two fiber orientation model parameters
	Optimization Objective		Minimize warpage		Minimize deflection along	length direction and the maxi-	mum injection pressure	Minimize shrinkage	Minimize shrinkage, cycle time,	runner volume	Minimize shrinkage and part	weight		Minimize warpage		Minimize shrinkage, part	weight, flash	Minimize warpage and opti-	mize gate design	Minimize warpage		Minimize warpage		Minimize warpage, shrinkage	and weldline	Minimize warpage and shrink-	age	Minimize fiber orientation dis-	tribution prediction error
	Surrogate Mod-	eling Method	Kriging		Kriging			RBF	ANN		Linear Regression,	Gaussian Process		ANN		ANN		Kriging		Kriging		RBF		ANN-RBF		Kriging	ANN
	DoE Method		Optimal midpoint	Latin Hypercube	Latin Hypercube			Latin Hypercube	Adaptive DoE		Central composite,	maximin, Latin	Hypercube	Latin Hypercube		Taguchi		Adaptive DoE		Optimal Latin	Hypercube	Latin Hypercube		Orthogonal Array	Latin Hypercube	Central composite	Full factorial
	Author(s)		Gao and	Wang [5]	Chen et al.	[6]		Li et al. [7]	Cheng et	al. [8]	Villarreal-	Marroquin et	al. [110]	Shi et al. [10]		Xu and Yang	[111]	Wang et al.	[106]	Kang et al.	[107]	Kitayama et	al. [108]	Feng and	Zhou [2]	Mukras et	al. [112]	Ivan et al.	[11]

  DIFFERENTIAL SCANNING CALORIMETRY 2. at 10 • C/min for the measurement runs with cooling rates of 3, 5, 10 • C/min, 3. at 40 • C/min for the measurement run with a cooling rate of 40 • C/min.

• C/min, 2.2.

Table 2 .

 2 1 -Overview of the DSC measurement runs and the studied samples' weight.

	Run Cooling rate Ṫ ( • C/min)	Sample weight m (mg)
	1 2	1	20.67 20.74
	3 4	3	10.34 10.93
	5 6	5	10.99 9.38
	7 8	10	5.06 5.93
	9 10	40	3.79 3.69

Table 2 .

 2 3 -Determined shot volumes and holding times from the performed filling and packing studies.

	Geometry (abbrev.)	Shot Volume (cm 3 ) Holding time (s)
	Full length, 3 mm thick (FL3)	42.0	16
	Half length, 3 mm thick (HL3)	28.3	16
	Full length, 1.5 mm thick (FL15) -	10
	Half length, 1.5 mm thick (HL15) 26.8	10

• C is specified and for the packing studies the mold temperature is set at 110 • C and the holding pressure at 80 MPa. A summary of the determined values is given in Table 2.3. One important thing 2.4. INJECTION MOLDING EXPERIMENTS

Table 2 .

 2 4 -The upper and lower bounds of the variable process settings used in the injection molding trials.

	Process Settings	Lower Bound			Upper Bound
	FL3 HL3 FL15 HL15	FL3 HL3 FL15 HL15
	Injection Velocity, V inj (cm 3 /s) 5	10	10	10	15	50	50	50
	Mold Temperature, T mold ( • C) 80	80	80	80	110 110 110	110
	Holding Pressure, P hold (MPa) 80	20	50	50	110 80	110	110

Table 2 . 5

 25

	S1	80	5	80
	S2	80	15	110
	S3	95	5	80
	S4	95	15	110
	S5	110	5	80
	S6	110	15	110

-DoE of the injection molding experiments for the full length 3 mm thick part.

DoE T mold ( • C) V inj (cm 3 /s) P hold (MPa)

Table 2 .

 2 6 -DoE of the injection molding experiments for the half length 3 mm thick part.

	V1	80	10	20
	V2	80	30	20
	V3	80	50	20
	V4	80	10	50
	V5	80	30	50
	V6	80	50	50
	V7	80	10	80
	V8	80	30	80
	V9	80	50	80
	V10 95	10	20
	V11 95	30	20
	V12 95	50	20
	V13 95	10	50
	V14 95	30	50
	V15 95	50	50
	V16 95	10	80
	V17 95	30	80
	V18 95	50	80
	V19 110	10	20
	V20 110	30	20
	V21 110	50	20
	V22 110	10	50
	V23 110	30	50
	V24 110	50	50
	V25 110	10	80
	V26 110	30	80
	V27 110	50	80

DoE T mold ( • C) V inj (cm 3 /s) P hold (MPa)

Table 2 .

 2 [START_REF] Li | A modified global optimization method based on surrogate model and its application in packing profile optimization of injection molding process[END_REF] -DoE of the injection molding experiments for the full length 1.5 mm thick part.DoE T mold ( • C) V inj (cm 3 /s) P hold (MPa)

	W1	80	10	50
	W2	80	30	50
	W3	80	50	50
	W4	80	10	80
	W5	80	30	80
	W6	80	50	80
	W7	80	10	110
	W8	80	30	110
	W9	80	50	110
	W10 95	10	50
	W11 95	30	50
	W12 95	50	50
	W13 95	10	80
	W14 95	30	80
	W15 95	50	80
	W16 95	10	110
	W17 95	30	110
	W18 95	50	110
	W19 110	10	50
	W20 110	30	50
	W21 110	50	50
	W22 110	10	80
	W23 110	30	80
	W24 110	50	80
	W25 110	10	110
	W26 110	30	110
	W27 110	50	110

Table 2 .

 2 8 -DoE of the injection molding experiments for the half length 1.5 mm thick part.DoE T mold ( • C) V inj (cm 3 /s) P hold (MPa)

	T1	80	10
	T2	80	20
	T3	80	30
	T4	80	10
	T5	80	20
	T6	80	30
	T7	80	10
	T8	80	20
	T9	80	30
	T10 95	10
	T11 95	30
	T12 95	50
	T13 95	10
	T14 95	30
	T15 95	50
	T16 95	10
	T17 95	30
	T18 95	50
	T19 110	10
	T20 110	30
	T21 110	50
	T22 110	10
	T23 110	30
	T24 110	50
	T25 110	10
	T26 110	30
	T27 110	50

Table 2

 2

	(a) 25x	(b) 50x	(c) 100x

.9, Figures 2.18, 2.19 and 2.20 present the micrographs of all mi-2.4. INJECTION MOLDING EXPERIMENTS

Table 2 .

 2 9 -The measured skin layer thickness obtained from analyzing the micrographs from the polarized optical microscope.

		Skin Layer Thickness (mm)
	V1	V3	V19	V21
	Cut 1 0.2-0.4 0.1-0.2	0.14-0.2 0.06-0.16
	Cut 2 0.2-0.3 0.11-0.13 0.12	0.09
	Cut 3 0.1	0.07-0.09 -	0.04

Table 3 .

 3 2 -Definition of the testing simulations for evaluating the performance of the surrogate models.

	DoE	V inj	T c,in ( • C) P hold	HT C f illing	HT C packing	D 3
	Name	(cm 3 /s)		(MPa)	(W/m 2• C)	(W/m 2• C)	(K/MPa)
	Test 1	36.0	98.2	55.1	5175	4522	0.33
	Test 2	15.5	87.3	74.7	9751	2143	0.07
	Test 3	45.6	109.4	47.5	8092	3945	0.21
	Test 4	27.2	106.2	25.9	7338	1640	0.12

Table 3 .

 3 3 -An overview of the proposed DoE and surrogate modeling techniques.

	Surrogate model	DoE Technique	Surrogate model	Fitting Method
	Acronym			Technique
	POD-NLR	Latin Hypercube	Proper orthogonal	Least-squares	re-
				decomposition		gression
	Regular-PCE-LS	Latin Hypercube	Polynomial	chaos	Least-squares	re-
				expansion		gression
	Adaptive-PCE-LS	Sequential	DoE
		(LHS & QMC)

Table 4 .

 4 2 -Material characteristics of the studied POM homopolymer obtained from the DSC measurements.

	Melting Temperature, T m (0)	Crystallization Enthalpy, ∆H c
	178.9 • C	155.35 J/g
	4.2.1.1 Heat Source

Table 4 .

 4 4 -The model parameters of the implemented viscosity model.

	Model Parameter Value (units)
	λ v	0.1608 (-)
	τ ⋆	398 000 (Pa)
	D 1	5.55 × 10 13 (Pa.s)
	A 1	26.596 (-)
	A 2	28.875 (K)

Table 4 .

 4 6 presents a summary of the implemented models in comparison to the ones used by the

	default Moldflow R ⃝ solver.

Table 4 .

 4 6 -Comparison between the implemented models and the ones used by default in a Moldflow simulation for the injection molding of a semi-crystalline thermoplastic material.

	Modeling	Default Moldflow R ⃝	Proposed Implementation
	Crystallization Not included	Poitou et al. [17] thermo-mechanical
			model: FIC & quiescent crystallization
	Viscosity	Cross-WLF model without crystalliza-	Modified Cross-WLF model with crys-
		tion dependency and a linear pressure	tallization coupling and a non-linear
		dependency	pressure dependency
	PVT	2-domain Tait model with a transition	Tait model coupled with crystallization
		temperature linearly dependent on pres-	and pressure
		sure
	Heat Source	Not included	Heat created due to crystallization and
			removed due to remelting of crystals
	4.3 Implementation in Moldflow using Solver API
	4.3.1 Solver API Feature

Table 4 .

 4 [START_REF] Li | A modified global optimization method based on surrogate model and its application in packing profile optimization of injection molding process[END_REF] -The process settings used in the basis simulations for the three generated surrogate models.

	Experimental	Surrogate	T c,in ( • C)	V inj (cm 3 /s)	P hold (MPa)
	DoE Name	Model Name
	V6	SM1	84	50	50
	V14	SM2	99	30	50
	V27	SM3	114	50	80

Table 4 .

 4 8 -The surrogate model variables with their upper and lower bounds.

	Surrogate Model Variables (units) Lower Bound Upper Bound

Table 4 .

 4 [START_REF] Gao | Intelligent methods for the process parameter determination of plastic injection molding[END_REF] -The experimentally-determined parameters used to define the six objective functions used in the optimization routine to identify the modeling parameters.

	k	Y exp,k max (MPa) t exp,k end (s)
	6	38.09	12.8
	14 35.25	15.2
	27 63.65	20.47

  .2. CONTENU ET STRUCTURE DE LA TH ÈSE L'objectif principal de cette thèse est l'utilisation de la modélisation de substitution dans la calibration des paramètres d'un modèle de cristallisation implémenté dans une simulation de moulage par injection. Afin d'atteindre cet objectif, une étude de faisabilité est d'abord réalisée pour évaluer l'utilisation de différentes techniques de modélisation de substitution et de plans d'expériences dans la calibration des paramètres de modélisation dans la simulation de moulage par injection d'une pièce thermoplastique en polyoxyméthylène. Cette étape est suivie de l'implémentation d'un modèle de cristallisation modifié de Poitou et al. ainsi que de modèles de viscosité, de PVT et de solidification dépendant de la cristallisation dans le logiciel commercial Autodesk R ⃝ Moldflow R ⃝ Insight 2021 en utilisant la fonction Solver API. Les modèles sont mis en oeuvre en 3D et calibrés à l'aide de trois modèles de substitution représentant la pression simulée dans différentes conditions de traitement. La capacité des modèles calibrés à prédire avec précision la pression, l'état de remplissage et l'épaisseur de la couche cristalline est testée dans différentes conditions en utilisant une géométrie de moule de base et trois variations de celle-ci. De plus, afin d'améliorer les prédictions de remplissage global ainsi que les estimations de gauchissement et de retrait, le critère de solidification doit être couplé à la cristallisation au lieu d'utiliser une valeur de température constante sans écoulement pour déterminer la solidification. Néanmoins, comme nous l'avons vu dans la section 1.3, le processus de cristallisation dans des conditions d'écoulement est complexe et la physique qui le sous-tend n'est pas encore totalement comprise, ce qui rend ce problème difficile à résoudre. De plus, la plupart des modèles disponibles sont des modèles phénoménologiques avec un grand nombre de paramètres qui sont difficiles à déterminer expérimentalement en raison de l'impossibilité de recréer les conditions extrêmes rencontrées pendant le processus de moulage par injection en utilisant les dispositifs expérimentaux actuels. Une façon de surmonter 5.2. CONTENU ET STRUCTURE DE LA TH ÈSE cet inconvénient est d'utiliser un modèle de substitution pour représenter un résultat de simulation qui peut être mesuré expérimentalement sur la machine de moulage par injection et de l'utiliser pour identifier les paramètres du modèle comme Ivan et al. [11] l'ont fait pour les paramètres du modèle d'orientation des fibres. Par conséquent, dans cette thèse, l'objectif principal est d'implémenter un modèle de cristallisation dans une routine de simulation d'injection 3D tout en prenant en considération l'effet de la cristallisation sur les modèles de viscosité, PVT et solidification. Afin d'identifier les paramètres expérimentaux difficiles à déterminer, des modèles de substitution sont générés pour représenter les résultats de pression obtenus par la simulation et pour effectuer une calibration à l'aide de mesures expérimentales de signaux de pression obtenues lors d'essais de moulage par injection.

	Une approche com-
	mune utilisée pour décrire la cristallisation induite par l'écoulement, observée pendant le moulage par
	injection, est basée sur la théorie de nucléation et de croissance d'Avrami-Kolmogorov [12,13] combinée
	aux équations de vitesse de Schneider [14]. Cette méthode nécessite la détermination expérimentale
	d'un grand nombre de paramètres du modèle et la résolution d'un système complexe et interdépendant
	d'équations différentielles [15]. En revanche, Poitou et al. [16,17] ont proposé une approche thermomé-
	canique décrivant à la fois la cristallisation induite par l'écoulement et la cristallisation naturelle avec

un nombre beaucoup plus faible de paramètres de modèle, ce qui rend son utilisation plus attrayante dans les logiciels de simulation. 55.2 Contenu et Structure de la Thèse La simulation du moulage par injection constitue un moyen rentable de concevoir et d'améliorer la qualité des pièces avant de produire un moule, ainsi que d'identifier les paramètres optimaux du processus pour la production. Toutefois, pour y parvenir, les simulations doivent être très précises et comporter de moins en moins d'approximations. Dans le cas des thermoplastiques semi-cristallins, cela signifie qu'il faut prendre en compte le processus complexe de cristallisation qui a lieu pendant le traitement du matériau, ainsi que son effet sur d'autres propriétés telles que la viscosité et le volume spécifique.

1.6. SUMMARY AND THESIS STRUCTURE

Les modèles sont entraînés à l'aide de 30, 60 et 120 cycles de simulation et leurs prédictions sont évaluées en conséquence. Les trois modèles sont construits en faisant varier un total de six paramètres d'entrée, y compris les paramètres du processus et les variables de modélisation. Les paramètres du processus sont la vitesse d'injection, la température d'entrée du liquide de refroidissement, la pression de maintien et les paramètres de modélisation, comme déjà mentionné, sont les coefficients de transfert de chaleur pendant le remplissage et le conditionnement ainsi que le coefficient dépendant de la pression dans le modèle de viscosité Cross-WLF. Après avoir obtenu les modèles de substitution, une analyse de sensibilité est effectuée pour étudier la contribution de chaque paramètre varié au signal de pression ainsi que leur interaction les uns avec les autres. Enfin, en utilisant les données de pression expérimentales mesurées et l'un des modèles de substitution générés, les paramètres de modélisation sont calibrés et l'erreur entre la simulation optimisée et la simulation par défaut est comparée. Cette étude a montré que le modèle de substitution POD-NLR et le modèle de substitution PCE-LS régulier présentent fondamentalement les mêmes performances, en particulier lorsque le nombre de données d'entraînement est faible. En augmentant le nombre de données d'entraînement, la technique PCE-LS adaptative améliore régulièrement la précision du modèle de substitution et atteint une performance comparable à celle des techniques à DoE fixe pour 120 simulations d'entraînement. Les avantages de la technique adaptive-PCE-LS devraient apparaître plus clairement lorsque le nombre de variables est élevé.Sur la base d'une analyse de sensibilité 'Sobol' utilisant le modèle de substitution généré, nous avons constaté que le coefficient de transfert de chaleur pendant le remplissage a un faible impact sur l'évolution de la pression au niveau d'un noeud de capteur devant le coefficient de transfert de chaleur pendant le conditionnement ou le coefficient de viscosité dépendant de la pression. En optimisant les paramètres HT C f illing , HT C packing et le paramètre D 3 du modèle de viscosité dépendant de la pression à l'aide des signaux de pression mesurés, il a été possible de calibrer la simulation haute-fidélité afin de diminuer les déviations dans l'estimation de la pression par rapport à celles obtenues en utilisant les valeurs de simulation par défaut. Cependant, il existe encore une lacune dans l'obtention d'une estimation précise des champs de pression dans le moulage par injection, qui pourrait être comblée avec cette méthodologie en incluant des paramètres de modèle de simulation supplémentaires en plus d'un raffinement des modèles mis en oeuvre dans la simulation, comme le montre le chapitre suivant en mettant en oeuvre un modèle de cristallisation induite par l'écoulement ainsi qu'un diagramme

La performance des modèles calibrés est également évaluée en utilisant trois variations de la géométrie du moule de base. Les prédictions de pression pour ces géométries utilisant les modèles calibrés n'étaient pas satisfaisantes. Cela a mis en évidence différentes lacunes des modèles mis en oeuvre ainsi que leurs paramètres identifiés. En modifiant l'épaisseur du moule, il n'était plus possible de prédire avec précision le remplissage partiel. En revanche, en ayant un canal plus long avec la même épaisseur de base, il est devenu difficile de prédire les pièces entièrement remplies. Ces comportements peuvent être dus à diverses raisons. L'une d'entre elles est l'utilisation de coefficients de transfert de chaleur globaux pour les phases de remplissage et de post-remplissage au lieu de coefficients locaux.Cela peut conduire à des prédictions de profil de température inexactes qui affectent l'ensemble des calculs de bilan thermique et, par conséquent, les prédictions de cristallisation et de solidification. Une autre raison pourrait être que les intervalles de paramètres choisis, pour lesquels les modèles de substitution sont générés, ne sont pas assez grands puisque certains des paramètres optimisés se situent aux extrémités de ces intervalles. Une autre explication pourrait être que le processus de calibration a été limité par la géométrie de la porte et les phénomènes physiques complexes qui s'y produisent, comme la solidification précoce due à des cisaillements élevés et à un refroidissement rapide.Pour les travaux futurs, on peut envisager de revoir la conception du moule en supprimant le seuil d'injection et en effectuant les essais expérimentaux avec un remplissage contrôlé par la pression. Cela

Acknowledgments

First of all, I would like to take this opportunity to express my deepest gratitude to my university supervisors Prof. Amine Ammar and Prof. Gilles Régnier. Their support and guidance throughout my PhD time was instrumental to my work and its success. They were always ready to offer me their time and knowledge to help me advance in my research. I consider myself very lucky to have had them

xiii CONTENTS the various micrographs, it became apparent that the higher the mold temperature and the injection velocity the thinner the skin layer is.

IMPLEMENTATION IN MOLDFLOW USING SOLVER API

the solver during an analysis or to perform specific actions.

• User functions are defined by the user's own code and are called by the solver during an analysis when needed. Function calls depend on the options selected in the Solver API window in the process settings of the current running study. Called functions need to be defined in the compiled shared library file or the solver will exit with failure.

Implementation Code

The models presented in Section 4.2 are implemented using a C++ code in a Visual Studio project. This project file includes a class for each user-defined model in which the class member functions are used by the main source file to calculate the viscosity, crystallization and specific volume. The main source file plays a major role in the implementation of the various models using the Solver API feature as it defines the user functions called by the solver during an analysis. The structure of this file resembles any other source file, it begins with #include directives which are read and interpreted by the preprocessor when the Visual Studio project is compiled. In this work, the following C++ Standard Library headers are used:

• input/output library: <iostream>, <sstream>, <fstream>,

• error handling: <cassert>,

• numerics library: <cmath>, <math.h>.

To be able to obtain the relative crystallinity α at a specific node in the simulation at time t, the crystallization model equation given by equation 4.1 needs to be solved in time and space. This is possible using the new Solver API advection functionality through the new user functions recognized as UserNodeFields. The user's code can access the integration results of the defined equation using the GetUserNodeField utility function.

In total, five main user functions are defined in the C++ code to compute:

• relative crystallinity: SolverUserHb3dUserNodeFieldsEvolution,

• viscosity: SolverUserHb3dViscosityAtNode,

• specific volume: SolverUserHb3dSpecificVolumeAtNode,

• heat generation/loss: SolverUserHb3dTemperatureIncrementAtNode,

Results and Discussion

Surrogate Models Performance

The performance of each generated surrogate model is assessed according to how well it is able to recreate the full pressure signal results of the high-fidelity simulation at one sensor location. Therefore, in this section, we present some example predicted pressure results from both the training and testing DoE sets as well as the error metrics of the POD parameters and time shifts obtained using the trained surrogate model, in other words, using the fitted polynomials defined by equations 3.2, 3.3 and 4.16.

The error metrics presented in this chapter are the normalized RMSE defined in equation 3.5 and the coefficient of determination also known as the R 2 score. It is a statistical measure that indicates how well the data fit the regression model and how well the unseen samples are likely to be predicted by the model. The best score is 1.0 and has a range between -∞ < R 2 ≤ 1. A value far from 1 means that the surrogate model response Ŷ is unable to predict accurately the true output results Y . R 2 is calculated as follows:

In order to support the presentation of the previously-mentioned results, the POD basis functions of each surrogate model are presented in Figure 4.3. The SM1 and SM2 models are reduced using five basis functions whereas the SM3 model is reduced using only three. The reason behind this difference is the use of the same truncation error ϵ for all models. Additionally, since the pressure signal for the V27 process settings is less complex than the one for V6 and V14, the SM3 model necessitates less basis functions to represent it.

Starting off by the SM1 surrogate model, Figure 4.4a presents the pressure results of four training simulations and the approximated results for these same simulations obtained by the model after its training. Even though the model is trained using these simulations, there still exists some discrepancy between the model predictions and that of the simulation especially in the case of the green curve.

These differences are expected since a non-interpolation technique is used to train the model as the signal is first reduced using POD and then its parameters are fitted to generate the final model. It is also clear that the form of the pressure result in the case of the green curve is different than the other reference simulation results. For the training simulations, the amplitude of the pressure during packing in the turquoise curve is under-predicted whereas for the time shift t end in the green and pink curves is over-predicted. As for the testing simulations, it seems that the SM2 model is performing very well on the four presented examples in Figure 4.6b. This is confirmed by analyzing the error metrics given by Figure 4.7 where the overall RMSE for all parameters is lower than 0.2. Additionally, there is very minimal difference between the model's ability in predicting seen and unseen data when comparing the error between the training and testing results.

As for the SM3 model, Figures 4.8a

mm Thick Geometry

The first geometry considered is the one with the same cavity length as the main geometry (HL3) but with half of its thickness with 1.5 mm. This geometry is referred to as HL15. 4.19c, both simulation results are nearly the same and seem to recreate accurately the experimental pressure signal at P1. However, it is important to remember that the models behind these simulations are quite different and it is not expected that they produce the same result. This highlights that the implemented models with the crystallization-dependent solidification criterion performs similarly to the no-flow temperature criterion before the gate where the same runner system is used for both the HL3 and HL15 geometries. This shows that the proposed models are able to reproduce the no-flow temperature criterion for some special cases such as in this case where the incomplete solidification of the runner due to a lower cooling time led to the same pressure results for both simulations.

At the P3 sensor location seen in Figures 4.19b and 4.19d, the pressure prediction similarities of the two simulations is not anymore observed and the material solidifies earlier using the no-flow condition than using the crystallization-dependent criterion. Nevertheless, as shown in Figure 4.19b, the experimental P3 signal is zero for the T2 condition which signifies a shorted part. However, both simulations are not able to reproduce this result. Even when no short shot is produced experimentally such as it is the case for T25 shown in Figure 4.19d, both simulations are unable to predict the signal accurately. For this geometry, the pressure signals in the part along with the fill predictions are not as well-predicted as the ones in the main HL3 geometry. This is due to the under-prediction of the crystallization evolution which could result form multiple sources. One possible explanation is the low identified value of the HT C p at 500 W m -2 • C -1 . As the mold thickness decreases, the heat exchange between the metal and the polymer melt is more prominent. This can directly relate to the bad prediction of short shots as the temperature profile is probably wrongly predicted leading to a more delayed solidification than observed experimentally. To summarize, the parameters calibrated using the HL3 geometry are unable to estimate the pressure signals correctly in the part cavity as well as the final filling status of the half length 1.5 mm thick geometry.

RESULTS AND DISCUSSION

Full Length 3 mm Thick Geometry

The second geometry considered is the one with the same cavity thickness as the main geometry (HL3) but with the whole length in order to make a U-shaped part. This geometry is referred to as FL3. Figure 4.20 presents two experimental pressure results in comparison to the two different simulation results at P1 and P3 for the two processing conditions S1 and S5. At the first glance, it is apparent that both simulations are unable to predict accurately the experimental pressure results especially during the packing phase of the process. By focusing on the simulation predictions using the implemented models with the optimized parameters presented in Table 4.10, it is clear that the filling phase is well-predicted but due to the early solidification of the gate the packing pressure is not delivered to the part cavity. This phenomenon is observed due to the prolonged shearing experienced along the 1 mm thick gate leading to higher flow-induced crystallization and therefore solidification followed by the restriction of the gate's cross-section. It is important to point that the injection velocity used for S1 and S5 is 5 cm 3 /s which is slightly outside the bounds of the process settings of the HL3 geometry. With such a low velocity, the filling of the FL3 geometry takes four times longer than when using the lowest velocity of 10 cm 3 /s for the HL3 geometry. The calibration using the HL3 geometry seems to be sensitive to the various complex physical phenomena happening at the gate leading to early solidification due to the fast FIC kinetics.

As for the filling predictions for this geometry, the simulation using the calibrated model predicted two short shots for the S2 and S4 conditions, however, experimentally no short shots were produced for all six conditions given in Table 2.5.

Full Length 1.5 mm Thick Geometry

The third and final geometry considered is the one with the full length U-shaped cavity having a 1.5 mm thickness. This geometry is referred to as FL15. For all experimental conditions presented in Table 2.7, the part produced was a short shot and the machine reached the maximum injection pressure of 240 MPa during the filling phase. Figure 4.21 presents two experimental pressure results in comparison to the two different simulation results at P1 and P3 for the two processing conditions W2 and W10. Both of the simulation pressure predictions at P1 for the W2 condition seen in Figure 4.21a predict correctly the short shot. However, after the machine maximum pressure is reached in both simulations around 0.5 s, the filling is switched to being pressure-controlled. In both these cases, there

Conclusions and Outlook

In this thesis, after having proved the feasibility of using surrogate modeling to identify modeling parameters by a simple optimization scheme, three surrogate models are used to calibrate userdefined models implemented in an injection molding simulation. The commercial software Autodesk R ⃝ Moldflow R ⃝ Insight 2021 is used to simulate the injection molding of a semi-crystalline thermoplastic material. In Moldflow R ⃝ , four user models are defined; a crystallization model describing both the flow-induced and quiescent crystallization along with crystallization-dependent viscosity, PVT and solidification models. These models are implemented in order to take into account the crystallization process happening during the processing of the material as it affects the flow solution as well as the predictions of warpage and shrinkage in the simulation. A multi-objective optimization is then performed to identify the model parameters using the generated surrogate models with the help of experimentally-measured pressure signals.

The simulation using the proposed models along with their identified parameters showed major improvements in the predictions of pressure signals in various conditions for the base mold geometry.

One of the crucial advances made by using this implementation was the ability of the simulation to accurately predict the filling status of the produced parts (short shot or fully filled) . Additionally, due to the inclusion of a crystallization model, it became possible to observe the relative crystallinity evolution along the different phases of the injection molding process which helped in estimating the skin layer thicknesses formed under high-shear conditions.

The performance of the calibrated models is additionally assessed using three variations of the base mold geometry. The pressure predictions for these geometries using the calibrated models were not satisfactory. This highlighted different shortcomings of the implemented models along with their identified parameters. By changing the mold's thickness, it was not anymore possible to accurately predict the short shots. Whereas, by having a longer channel with the same base thickness, it became Appendix A

Cross-WLF Viscosity Model

The Cross-WLF viscosity model [START_REF] Cross | Rheology of non-Newtonian fluids: A new flow equation for pseudoplastic systems[END_REF] describes the temperature, shear rate, and pressure dependency of the viscosity for thermoplastic materials. This model is used in Autodesk Moldflow Insight 2021.1 to calculate the viscosity of the polymer during its injection molding.

where:

• η is the viscosity of the melt,

• η 0 is the zero shear viscosity,

• γ is the shear rate,

• τ * is the critical stress at the transition to shear thinning,

• n is the power law index in the high shear rate regime.

The zero shear viscosity parameter, η 0 , in the above equation is given by the WLF model [START_REF] Williams | The Temperature Dependence of Relaxation Mechanisms in Amorphous Polymers and Other Glass-forming Liquids[END_REF]:

where:

• T is the temperature,

• T * = D 2 + D 3 P is the glass transition temperature,

• A 2 = A 3 + D 3 P ,

• P is the pressure,

• A 1 , A 2 , D 1 , D 2 and D 3 are data-fitted coefficients.