# Efficient identification of a flow-induced crystallization model for injection molding simulation

Sandra Saad, Camilo Cruz, Gilles Régnier, Amine Ammar

<figure><img src="image_1.png" /></figure>

## To cite this version:

Sandra Saad, Camilo Cruz, Gilles Régnier, Amine Ammar. Efficient identification of a flow-induced crystallization model for injection molding simulation. International Journal of Advanced Manufacturing Technology, 2024, 133 (11-12), pp.5243-5270. 10.1007/s00170-024-13961-6. hal-05228688

HAL Id: hal-05228688

https://hal.science/hal-05228688v1

Submitted on 1 Sep 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Efficient identification of a flow-induced crystallization model for injection molding simulation

Sandra Saad

sandra.saad@de.bosch.com

Robert Bosch GmbH https://orcid.org/0000-0001-6978-5430

Camilo Cruz

Gilles Régnier

Amine Ammar

## Research Article

**Keywords:** Crystallization, Injection Molding, Surrogate Modeling, Process Simulation, Multi-objective Optimization, Thermoplastics, Parameter Calibration, Polymers

**Posted Date:** March 22nd, 2024

**DOI:** https://doi.org/10.21203/rs.3.rs-4044458/v1

**License:** ©⑈ This work is licensed under a Creative Commons Attribution 4.0 International License.
Read Full License
# Efficient identification of a flow-induced crystallization model for injection molding simulation

Sandra Saad¹,²*, Camilo Cruz¹, Gilles Régnier³ and Amine Ammar²

¹*Corporate Sector Research and Advance Engineering, Robert Bosch GmbH, Robert-Bosch-Campus 1, Renningen, 71272, Germany.

²LAMPA Laboratory, Arts et Métiers Institute of Technology, Boulevard du Ronceray 2, Angers, 49035, France.

³PIMM Laboratory, Arts et Métiers Institute of Technology, Boulevard de l'Hôpital 151, Paris, 75013, France.

*Corresponding author(s). E-mail(s): sandra.saad@de.bosch.com;

Contributing authors: camilo.cruz@de.bosch.com; gilles.regnier@ensam.eu;

amine.ammar@ensam.eu;

## Abstract

Most commercial software used to simulate the injection molding process of semi-crystalline thermo-plastic polymers do not explicitly take into account the polymer crystallization, which could lead to errors in the estimations of filling as well as 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 under injection molding conditions. To close this gap, in this work we use a simple thermo-mechanical crystallization model to describe the flow-induced and quiescent crystallization of an unreinforced semi-crystalline thermoplastic material during injection molding. The crystallization model is implemented in the commercial software Autodesk® Moldflow® Insight 2021 using the Solver API feature alongside crystallization-dependent viscosity, PVT and solidification models. The model parameters were identified using a calibration workflow that employs surrogate models representing the simulated pressure results to perform a multi-objective optimization. The filling 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 (POM) part using different process conditions. The results show major improvements in the estimations of the time-dependent pressure field as well as the level of filling of the produced parts.

**Keywords:** Crystallization, Injection Molding, Surrogate Modeling, Process Simulation, Multi-objective Optimization, Thermoplastics, Parameter Calibration, Polymers

# 1 Introduction

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 is a numerical implementation to solve a
d crystallization model for injection molding simulation

set of conservation equations using material mod-
els [1]. Many commercial and academic software
have been developed to solve this problem, which
involves flow, diverse heat transfer mechanisms
in addition to a phase change as well as time-
dependent boundary conditions of the polymer
material domain during the various stages of the
process. Taking into account all these physical
phenomena in the simulation is not straightfor-
ward and necessitates some simplifying assump-
tions. An overview of those hypotheses with a
focus on the injection molding of an unfilled
semi-crystalline thermoplastic [1] is given in the
following:

* Specific heat capacity: Data taken from differen-
tial scanning calorimetry (DSC) tests performed
at constant cooling rates, typically quite lower
than the ones encountered close to the mold
wall in the actual process. The data includes the
crystallization peak measured during the DSC
test.

* Heat sources: Latent heat of crystallization is considered indirectly by means of the specific heat capacity data mentioned above, i.e. from a single constant cooling rate. Therefore, those solvers do not implement explicitely a heat source term related to crystallization in the energy equation.

* Solidification: This phenomenon is approximated numerically by using a single no-flow temperature (or transition temperature) below which the local melt viscosity is set to a very high value (~ $10^6$ Pa.s). This implementation ignores the changes in crystallization temperature induced by different local cooling rates occurring in the actual process.

* Pressure-Volume-Temperature (PVT): Data taken from experiments using a constant cooling rate, typically far below than the ones in the actual process, and therefore neglecting the changes in crystallization kinetics triggered by different cooling rates.

Due to the complexity of experimentally identi-
fying crystallization parameters under injection
molding conditions and the effect of crystallization
on material properties during processing, many
commercial software do not take it into account
by default, even when dealing with the injection
of semi-crystalline materials. This could lead to
inaccuracies in the fill predictions as well as in

the estimation of warpage and shrinkage in the
simulation [1].

To overcome the assumptions mentioned above
and, at the same time, to maintain the effort
in material characterization as low as possible,
in this work, we chose a simple crystallization
model able to consider flow-induced and quiescent
phenomena [2] to simulate explicitly crystalliza-
tion during injection molding. From this point
on, this work has two main contributions to the
state-of-the-art in injection molding simulation.
The first one is the numerical implementation of
a flow-induced crystallization model in a com-
mercial software for injection molding simulation.
Concretely, we present the implementation of a
flow-induced crystallization model coupled with
melt viscosity, PVT and solidification models in
Autodesk® Moldflow® Insight 2021 by using its
solver Application Programming Interface (known
as Solver API). The second contribution is a
novel material-parameter identification method-
ology, which uses a surrogate model-based opti-
mization taking as reference experimental sensor
data from an actual injection molding process.
One innovation here is the material-parameter
identification under actual processing conditions,
which differs from traditional approaches where
the crystallization models are calibrated based on
offline experiments with thermomechanical condi-
tions far away from industrial injection molding
ones.

In-mold pressure signals are reliable state-of-
the-art sensing equipment in industrial context.
That is one of the reasons why we chose in-mold
pressure sensor signals as reference data for the
model calibration. Another important reason is
because this physical quantity is strongly depen-
dent on the crystallization process of the injected
semi-crystalline thermoplastic material. In other
words, the simulated pressure field is highly sen-
sitive to the crystallization model parameters
that we look for calibrating. A subjacent working
hypothesis here is that the rest of parameters in
the simulation (of material or numerical nature)
are correctly calibrated and that other computed
physical quantities (e.g., temperature field) are
accurate.

This work looks for elucidating two research
hypotheses. On the one hand, that the implemen-
tation of a flow-induced crystallization model in
injection molding simulation increases the robust-
ness of the mold filling predictions, particularly
in terms of pressure-field evolution and incom-
plete fillings (or short-shots). On the other hand,
that the proposed method of material parame-
ter identification is especially advantageous for
the modeling of crystallization in injection mold-
ing. In fact, we propose that the generation of a
surrogate model of a high-fidelity injection mold-
ing simulation model (where the surrogate input
variables are the unknown material parameters in
the high-fidelity simulation model) and the subse-
quent identification by calibrating the generated
surrogate model is more efficient than the tradi-
tional approaches based exclusively on lab exper-
iments. In a previous work dealing with default
injection molding simulation [3], this methodology
performed well for a low number of training data
and input parameters.

This paper is structured in four main sections.
In the next section, we revisit the literature about
polymer crystallization and its modeling as well
as the numerical methods associated to the meta-
model generation and its application in injection
molding simulation. In the section 3, we show the
implementation of the flow-induced crystallization
model in Autodesk® Moldflow® Insight as well
as the respective material parameter identification
process by means of a pilot case using a poly-
oxymethylene material. In the subsequent section,
we discuss the performance of the identified mod-
els at the level of the high-fidelity simulation in
terms of pressure estimation and prediction of
short shots. Finally, we close the paper with a
conclusion.

2 Theoretical Background

## 2.1 Polymer Crystallization

2.1.1 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 pro-
cess and therefore only the crystallization from
the melt is considered. This process involves two
stages:

1. Nucleation or the formation of active nuclei in the liquid phase acting as starting points for the appearance of crystals [1]. Two types of nucleation can be distinguished [4, 5]:

*   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 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 microstructures having different spatial growing mechanisms [1, 6]:

* Spherulites: Grow radially in space and form spherical structures known as spherulites, typically seen under no-flow conditions (quiescent crystallization).

*   Thread-like: Oriented semi-crystalline microstructures in which crystals mainly grow perpendicularly to flow direction, leading to deformed spherulites, shish-kebab structures and even fibrils for high strain rates (flow-induced crystallization).

The flow conditions during injection molding
affects the crystallization from the melt in several
ways [7]:

* Increase the nucleation density,

* Raise the crystallization kinetics,

* Induce changes in the semi-crystalline morphological structures, and

* Increase the crystallization and melting temperature.

Comprehensive research work has been done in the
last few decades to better understand the physics
behind polymer crystallization whether under qui-
escent or flow conditions. An extensive review of
the various theories postulated is given by Zhang
et al. [8]. The current understanding of this pro-
cess under flow conditions is that there exists a
Weissenberg number threshold below which no
changes in crystal morphology are observed and
above which the nucleation density and growth
d crystallization model for injection molding simulation

rate are altered by the intensity of flow in addition to a second threshold where thread-like (fibrillar) morphologies start appearing [9, 10].

### 2.1.2 Modeling Approaches

The first modeling approaches used to describe the kinetics of quiescent crystallization are based on the Kolmogorov-Avrami-Evans (KAE) model [11–13] which was developed considering isothermal conditions. Such model describes the nucleation and growth of spherulites using integral equations while taking into account crystal impingement. To include a temperature dependence into the modeling, it is common to use the Hoffman-Lauritzen model which describes the temperature-dependent growth rate constant [14]. Additionally, Ozawa [15] extended the phenomenological KAE model to adhere for the non-isothermal case but it was found that it applies for a limited range of cooling rates. Nakamura et al. [16] introduced an isokinetic approach to describe the kinetic constant in the KAE model. Thereafter, Ziabicki [17] developed a generalized theory to predict non-isothermal crystallization kinetics using an empirical relation for the rate constant. As for the description of the temperature-dependent nucleation evolution, a relation based on Koscher & Fulchiron’s experimental work [18] is commonly used in literature. An alternative approach for the modeling of non-isothermal crystallization is proposed by Schneider et al. [19]. This approach uses a set of four first-order differential equations derived from the KAE model to solve for the relative crystallinity. Additionally, the solution of Schneider’s equations provides morphological information concerning the formed crystalline structures.

The pioneering work of Eder et al. [20] set the stage for the modeling of flow-induced crystallization (FIC). They based their model on the shear rate as a driving force for crystallization and proposed a model dependent on the shear rate to account for the effects of flow. Other researchers like Doufas et al. [21] took into account the effect of flow by defining a stress intensification factor. Zuidema et al. [22] proposed a strain-based model describing the spherulite evolution using the Schneider equations and a modified Eder model to account for the evolution of shish-kebabs.

Guo and Narh [23] and Kim et al. [24] developed a model to describe the molecular orientation based on Nakamura’s model. As for Coppola et al. [25] and Titomanlio and Lamberti [26], who correlated the increase of free energy to the raise of melting temperature and decrease of induction time. Zinet et al. [27] considered the first invariant of the extra-stress tensor as driving force and extended the Schneider equations accordingly. A thermo-mechanical based approach is used by Poitou et al. [2] to describe the flow-induced crystallization in the framework of generalized standard materials. Roozemond et al. [28] used the modified Schneider equations to account for the growth of both spherulites and shish-kebab crystallites. Most recently, Spina et al. [29, 30] used a multi-scale approach to model flow-induced crystallization.

To sum up, in the literature we can find a broad range of 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:

* The model family based on Nakamura’s approach: The increase of crystallization kinetics is modelled by multiplying the kinetic function by an intensification 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.

* The family based on KAE model and/or Schneider’s rate equations: The increase of crystallization kinetics is coupled with an enhancement of the nucleation rate and/or growth rate of the crystalline structures. In some particular cases, the increase of nucleation rate is connected to an higher volumetric free energy difference between molten and crystalline phases.

Table 1 presents an overview of the literature in the field of FIC modeling categorized according to the two above-mentioned clusters and mentioning the proposed driving force describing the effect of the flow on the crystallization kinetics. These two modeling approaches are not strictly separated as some of the functions describing the driving force in the Nakamura-based approaches can be integrated into the Schneider rate equations to obtain morphological information.
<table><caption>**Table 1** Overview of the literature in the field of modeling the flow-induced crystallization during polymer processing, based on [9]</caption><thead><tr><th>Crystallization model based on</th><th>Effect of flow on crystallization kinetics using</th><th>Author(s)</th></tr></thead><tbody><tr><td rowspan="5">Nakamura's model [16]</td><td>Multiplying factor function of the extra stress tensor</td><td>Doufas et al. [21], Zinet [27]</td></tr><tr><td>Multiplying factor function of the strain and the shear rate</td><td>Kulkarni et al. [31], Tanner et al. [32], Brahmia [33], Mu et al. [34, 35]</td></tr><tr><td>Melting temperature increase function of stress</td><td>Haas and Maxwell [36], Titomanlio et al. [37], Guo and Narh [23]</td></tr><tr><td>Melting temperature increase function of molecular strain</td><td>Titomanlio and Lamberti [26], Kim et al. [24], Pantani et al. [10, 38]</td></tr><tr><td>Irreversible thermodynamics based on the standard material formalism</td><td>Poitou et al. [2, 39]</td></tr><tr><td rowspan="3">KAE model [11–13], Schneidner's rate equations [19] (Nucleation and growth models)</td><td>Enhancement of nucleation rate</td><td>Eder et al. [20], Zuidema et al. [22], Koscher and Fuchiron [18], Roozemond et al. [28], Kim et al. [40]</td></tr><tr><td>Enhancement of growth rate</td><td>van Meerveld et al. [41], Zinet et al. [42], Roozemond et al. [43], Troisi and Arntz [44]</td></tr><tr><td>Enhancement function of the free energy contribution affecting the nucleation rate</td><td>Acierno et al. [45], Zheng and Kennedy [46], Zheng et al. [47], Laschet et al. [48], Schrank et al. [49]</td></tr></tbody></table>

### 2.1.3 Effect on Other Material Properties

During processing and forming of semi-crystalline thermoplastics, the development of crystalline structures influences other material properties such as viscosity, specific volume, thermal capacity and thermal conductivity. Therefore, substantial research effort has been made to develop models that describe the effect of crystallization on these material properties.

As it is difficult to measure the viscosity and the crystallinity simultaneously, separate testing methodologies are typically performed to assess the impact of the crystallization phenomenon on

the effective viscosity. Rheological measurements are done to describe the viscosity evolution as a function of temperature and shear rate. Whereas, DSC 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 [50]. The experimental data obtained from the previously mentioned methods are used to obtain models describing the viscosity enhancement due to crystallization. Those coupling functions can be divided in two main classes:

* Models taken from the framework of particles suspension rheology, and
d crystallization model for injection molding simulation

* Empirical equations recreating the abrupt increase in viscosity induced by crystallization.

Pantani et al. [9] and Lamberti et al. [50] present concise reviews on the different developed viscosity functions found in literature that take into account the crystallinity content in relation to viscosity.

As for the specific volume, several works have been presented by Luyé et al. [51], Fulchiron et al. [52], Zheng et al. [47] and Zhao et al. [53]. They included the relative crystallinity into the calculation of the specific volume by assuming a simple mixing rule of the molten and solidified phases' specific volumes and using Tait-based forms of equation of state. An alternative approach is to extend the classical modified 2-domain Tait pvT model by including a cooling-rate dependency as proposed by Cook et al. [54], Wang et al. [55] and Hopmann et al. [56].

## 2.2 Surrogate Modeling

Despite the advances in numerics and high performance computing, there is still a large number of engineering problems that are difficult to solve accurately and in an acceptable time frame [57]. This is mainly due to the strong coupling and multi-dimensionality of the subjacent physical models [58]. In response to this fact, there has been a growing interest in surrogate modeling techniques to address these challenges. Such methods approximate the response of complex models using a surrogate model also known as a metamodel, which is basically a supervised machine learning model trained on a finite set of input-output data from the complex (high-fidelity) model. 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, and sensitivity analysis.

### 2.2.1 Generation of a Surrogate Model

The construction of a surrogate model consists of multiple steps. The basic process can be summarized as follows [59–61]:

1. *Choice of input (design) variables:* Selection of input parameters, which presumably have a non-negligible impact on the model output. This choice is usually supported by preliminary experiments, whether physical or numerical ones.

2. *Sampling in the design space:* Definition of a sampling strategy also known as design of experiments and evaluation of the respective design points by means of a high-fidelity simulation or actual experiments.

3. *Supervised learning:* Selection of an emulator model in accordance with the problem at hand and model identification by fitting the input-output data obtained in the previous step.

4. *Model validation:* Assessing the performance of the surrogate model by calculating diverse statistical criteria for training and testing sets. In case of unsatisfactory results, the identification of new design points for further model enrichment is triggered.

5. *(Optional) Model updating:* Building an updated surrogate model using the additional design points along with the previous ones.

6. *Model exploitation:* Use of the generated surrogate model for further analyses such as parameter sensitivity analysis, optimization routines or uncertainty quantification purposes.

Figure 1 provides a summary of the mathematical/numerical techniques most commonly used 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), which has to be approximated by the surrogate model. Following this, the type of emulator is chosen along with the most suitable fitting (or regression) method. In the following sections, a closer look is taken into each of the four steps shown in Figure 1.

### 2.2.2 Surrogate Modeling in Injection Molding Simulation

In the last years, 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 [64] 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 profile. Similar works were
<figure><img src="image_3.png" /><figcaption>Fig. 1 Overview of the most commonly used methods during the generation of a surrogate model [62, 63]</figcaption></figure>

performed by Chen et al. [63], Wang et al. [65] and Kang et al. [66]. Other authors used radial basis function [67–69], artificial neural networks [70, 71], Gaussian process [72] as fast emulator to optimize process parameters for controlling shrinkage and warpage in the final part. Additional applications for surrogate models have been used for the optimization of cycle time [68] and part weight [73].

All of the previously mentioned publications use a surrogate model to find optimal process parameters. However, another interesting utilization of surrogate modeling in injection molding was recently published by Ivan et al. [74], 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 used in their high-fidelity injection molding simulation. This was done by generating an ANN-based surrogate model recreating the fiber orientation tensors across the thickness of several regions of interest in a plate as function of two fiber orientation model parameters. In this work, we use an analogous approach to calibrate a crystallization model implemented in a high-fidelity injection molding simulation.

Table 2 provides a brief overview of the literature dealing with surrogate modeling in injection molding simulation. More comprehensive reviews can be found in [71, 75].

# 3 Model Implementation and Parameter Identification

A workflow of the proposed material parameter identification is sketched in Figure 2. It includes the implementation of the crystallization model into a high-fidelity injection molding simulation software, the generation of a respective surrogate model (having as input variables the unknown material model parameters), the generation of reference data (in-mold pressure sensor signals) from actual injection molding as well as the optimization loop over the surrogate model for identifying the unknown material model parameters. In the following, we will present in detail those different blocks by using a polyoxymethylene (POM) homopolymer as pilot case.

## 3.1 Experimental Data Generation

### 3.1.1 Material

The polymer material used in this work is an injection molding grade of an unreinforced POM homopolymer. To characterize the quiescent crystallization, we carried out differential scanning calorimetry (DSC) measurements using a TA Instruments, Inc. DSC Q1000 with a heating rate of 10 °C/min and several cooling rates (1, 3, 5, 10 °C/min). The DSC thermograms of the different cooling rates are presented in Figure 3.
d crystallization model for injection molding simulation

<table><caption>**Table 2** Overview of the literature in the field of surrogate modeling based optimization in injection molding simulation</caption><thead><tr><td>Author(s) (Year)</td><td>DoE Method</td><td>Surrogate Modeling Method</td><td>Optimization Objective</td><td>Design Variables</td></tr></thead><tbody><tr><td>Gao and Wang (2009) [64]</td><td>Optimal mid-point Latin Hypercube</td><td>Kriging</td><td>Minimize warpage</td><td>Mold and melt temperatures, injection and packing times, packing pressure</td></tr><tr><td>Chen et al. (2010) [63]</td><td>Latin Hypercube</td><td>Kriging</td><td>Minimize deflection along length direction and the maximum injection pressure</td><td>Mold and melt temperatures, packing time, packing pressure</td></tr><tr><td>Li et al. (2010) [67]</td><td>Latin Hypercube</td><td>RBF</td><td>Minimize shrinkage</td><td>Packing and cooling times, packing pressure profile</td></tr><tr><td>Cheng et al. (2013) [71]</td><td>Adaptive DoE</td><td>ANN</td><td>Minimize shrinkage, cycle time, runner volume</td><td>Runner diameters, packing and cooling times, packing pressure</td></tr><tr><td>Villarreal-Marroquin et al. (2013) [73]</td><td>Central composite, maximin, Latin Hypercube</td><td>Linear Regression, Gaussian Process</td><td>Minimize shrinkage and part weight</td><td>Mold temperature, packing time</td></tr><tr><td>Shi et al. (2013) [70]</td><td>Latin Hypercube</td><td>ANN</td><td>Minimize warpage</td><td>Mold and melt temperatures, injection, packing and cooling times, packing pressure</td></tr><tr><td>Xu and Yang (2015) [76]</td><td>Taguchi</td><td>ANN</td><td>Minimize shrinkage, part weight, flash</td><td>Mold and melt temperatures, injection, packing and cooling times, packing pressure</td></tr><tr><td>Wang et al. (2015) [65]</td><td>Adaptive DoE</td><td>Kriging</td><td>Minimize warpage and optimize gate design</td><td>Mold temperature, ram speed, packing time and pressure</td></tr><tr><td>Kang et al. (2016) [66]</td><td>Optimal Latin Hypercube</td><td>Kriging</td><td>Minimize warpage</td><td>12 process settings</td></tr><tr><td>Kitayama et al. (2017) [68]</td><td>Latin Hypercube</td><td>RBF</td><td>Minimize warpage</td><td>Cooling temperature, injection, packing and cooling times, packing pressure</td></tr><tr><td>Feng and Zhou (2019) [69]</td><td>Orthogonal Array Latin Hypercube</td><td>ANN-RBF</td><td>Minimize warpage, shrinkage and weldline</td><td>Mold and melt temperatures, injection, packing and cooling times, packing pressure</td></tr><tr><td>Mukras et al. (2019) [77]</td><td>Central composite</td><td>Kriging</td><td>Minimize warpage and shrinkage</td><td>Mold and melt temperatures, injection, packing and cooling times, packing pressure, injection velocity</td></tr><tr><td>Ivan et al. (2022) [74]</td><td>Full factorial</td><td>ANN</td><td>Minimize fiber orientation distribution prediction error</td><td>Two fiber orientation model parameters</td></tr></tbody></table>
<figure><img src="image_4.png" /><figcaption>**Fig. 2** Sketch of the surrogate modeling-based optimization for parameter identification in injection molding simulation</figcaption></figure>

<table><caption>**Table 3** Material characteristics of the studied POM homopolymer obtained from DSC measurements</caption><thead><tr><th>Crystallization<br>Peak<br>Temperature,<br>T<sub>c</sub>(0)</th><th>Melting Peak<br>Temperature,<br>T<sub>m</sub>(0)</th><th>Sample Crys-<br>tallization<br>Enthalpy, ΔH<sub>c</sub></th></tr></thead><tbody><tr><td>149.4 °C</td><td>178.9 °C</td><td>155.3 J/g</td></tr></tbody></table>

Using the DSC data, we determined the crystallization (peak) and melting (peak) temperatures along with the crystallization enthalpy for the studied POM material. Table 3 summarizes those characteristical values.

Additionally, rheological measurements are performed using an ARES rheometer from TA Instruments, Inc. with a plate-plate geometry. The test is carried out with a constant frequency of 3 rad/s and a maximal shear strain of 0.03 % with a cooling rate of 3 °C/min. Figure 4 presents the results of the rheological dynamic tests with temperature ramp.

### 3.1.2 Injection Molding Experiments

Injection molding experiments are performed on an electrical injection molding machine (ENGEL

E-Motion 440/220 T). We used a constant cross-section cavity mold geometry, where four pressure sensors are mounted in the cavity to measure pressure signals during the injection molding process. A technical drawing of the part along with the design dimensions are shown in Figure 5. The sensors’ location is specified by the gray circles in Figure 5. By means of a filling study with a melt temperature of 220 °C, a mold temperature of 80 °C and an injection velocity of 10 cm³/s, we determined a shot volume of 28.3 cm³ that is used as fixed initial ram position for all subsequent injection molding experiments. On the other hand, a packing study is done with a mold temperature of 110 °C and a holding pressure of 80 MPa for defining a holding time of 16 s.

A design of experiments consisting of 27 sampling points is performed where the mold temperature, the injection velocity and the holding pressure are varied. The upper and lower bounds of these variables are presented in Table 4. All experimental runs are carried out with a melt temperature of 220 °C. To import the machine and process settings along with the sensor signals during the different injection molding runs, we used the ENGEL sim link interface.
<table><caption>**Table 4** Variable process settings in the injection molding experimental runs with their upper and lower bounds</caption><thead><tr><td>Process Settings</td><td>Units</td><td>Lower Bound</td><td>Upper Bound</td></tr></thead><tbody><tr><td>Injection Velocity, V<sub>inj</sub></td><td>cm<sup>3</sup>/s</td><td>10</td><td>50</td></tr><tr><td>Mold Temperature, T<sub>mold</sub></td><td>°C</td><td>80</td><td>110</td></tr><tr><td>Holding Pressure, P<sub>hold</sub></td><td>MPa</td><td>20</td><td>80</td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**Fig. 3** DSC thermograms of the studied POM under different constant cooling rates (1, 3, 5, 10 °C/min)</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 4** Complex viscosity measured under a cooling rate of 3 C/min with a constant frequency of 3 rad/s and a maximal shear strain of 0.03%</figcaption></figure>

Some of the experimental runs presented short shots due to the incomplete filling of the mold cavity caused by the early freezing of the gate.

## 3.2 Implementation of User-defined Models

### 3.2.1 Simulation Software

We selected Autodesk® Moldflow® Insight 2021.1 (AMI2021.1) as software to perform the high-fidelity injection molding simulations because it offers the possibility to implement user-defined functions by means of the Solver API feature.

The Solver Application Programming Interface (API) in Moldflow® enables the user to create their own C++ functions, which are called by the solver to calculate an user-defined variable (e.g., physical property, material descriptor, etc.) 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.

### 3.2.2 Crystallization Model

The implemented crystallization model is a thermo-mechanical based model developed by Poitou et al. [2] in the framework of irreversible thermodynamics. The model is derived using the standard material formalism commonly used in solid mechanics to describe various coupled phenomena [2]. Such formalism necessitates two potentials, a thermodynamic potential and a pseudo-potential, 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 [2]. More detailed information concerning this method can be found in [2, 39].
<figure><img src="image_7.png" /><figcaption>Fig. 5 Technical drawing of the injection-molded part including some characteristic dimensions in mm and the location of the four pressure sensors (P0, P1, P2 and P3)</figcaption></figure>

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 [16] and the potential referring to the mechanical constitutive behavior [2]. 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. In Poitou et al.'s model [2], the crystallization's kinetics are modeled using a temperature-dependent function $\chi(T)$ defined by Hieber [78]. However, in this work, this relation is substituted by a rate equation described by Lauritzen and Hoffman [14], known as the Hoffman-Lauritzen equation, and related to the previous by $\chi(T) = 1/K(T)^{1/n}$.

The crystallization model, which only applies for temperatures lower than the material's melting temperature $T_m(P)$, reads as follows:

$$ \dot{\alpha}(t) = (1 - \alpha) \left[ -\frac{1}{\beta} \frac{\partial \eta}{\partial \alpha} \text{tr}\mathbf{D}^2 + nK(T)^{\frac{1}{n}}(-\ln(1 - \alpha))^{1 - \frac{1}{n}} \right] \quad (1) $$

$$ K(T) = K_0 \exp \left( -\frac{3K_g}{T(T_m^0(P) - T)} \right) \times \exp \left( -\frac{3U^*}{R(T - T_\infty(P))} \right) \quad (2) $$

where $\beta$, $n$, $K_0$ and $K_g$ are data-fitted parameters, $\mathbf{D}$ is the strain rate tensor, $T_m^0(P)$ is the pressure-dependent 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_\infty = T_g(P) - 30$ K with $T_g(P)$ as the pressure-dependent glass transition temperature. The first term in equation 1 represents the contribution of the flow-induced crystallization, whereas the second term in equation 1 is the contribution of the quiescent one. This last term is a function of the quiescent crystal growth rate given in equation 2, which corresponds phenomenologically to the Avrami kinetic model parameter too.

To decrease the amount of unknown material model parameters, we determined a linear relationship between the kinetic constants of the Hoffman-Lauritzen equation, which is particular for the POM used in this work:

$$ K_g = 6732 \ln(K_0) + 11296. \quad (3) $$

This relationship defines the allowed couplets of $(K_g, K_0)$, which fit equation 2 to the values of $K_{Avrami}$ identified in the reduced range of temperatures that is typically accessible by experimental means. We determined, in turn, the Avrami
ed crystallization model for injection molding simulation

<figure><img src="image_8.png" /><figcaption>**Fig. 6** Fitting surface relating the different pressure dependencies of the melting temperature ($b_6$), equilibrium melting temperature (a) and glass transition temperature (b)</figcaption></figure>

kinetic constants (crystal growth rate at constant temperature) from the Ozawa kinetic constants (crystal growth rate at constant cooling rates) using the following expression:

$$K_{Avrami} = K(T) = \left[ -\frac{d}{dT} \left( K_{Ozawa} (T)^{1/n} \right) \right]^n \quad (4)$$

where *n* is the Avrami exponent and also appears in equation 1.

We needed this Avrami-Ozawa correlation because our DSC experiments were performed under constant cooling rates, so we identified first directly the Ozawa kinetic constants.

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

$$T_m(P) = T_m(0) + b_6 P, \quad (5)$$

$$T_m^0(P) = T_m^0(0) + aP, \quad (6)$$

$$T_g(P) = T_g(0) + bP \quad (7)$$

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

<table><caption>**Table 5** The model parameters of the implemented crystallization model</caption><thead><tr><td>Model Parameter</td><td>Value (units)</td></tr></thead><tbody><tr><td>n</td><td>3 (-)</td></tr><tr><td>U*</td><td>6270 (J/mol)</td></tr><tr><td>T<sub>m</sub><sup>0</sup>(0)</td><td>199.7 (°C) [79]</td></tr><tr><td>T<sub>m</sub>(0)</td><td>178.9 (°C)</td></tr><tr><td>T<sub>g</sub>(0)</td><td>-73.0 (°C) [80]</td></tr><tr><td>b<sub>6</sub></td><td>0.175 (K/MPa)</td></tr><tr><td>a</td><td>0.175 (K/MPa)</td></tr><tr><td>b</td><td>0.161 (K/MPa)</td></tr></tbody></table>

given in equation 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 6 shows the surface representing the relation between the parameters, which seems to follow a linear relationship:

$$b_6 = f \cdot a + g \cdot b \quad (8)$$

where $f = 0.629$ and $g = 0.4015$ are the obtained fitted parameters. In this work, we suppose that the melting temperature and the equilibrium melting temperature have the same pressure dependency such as $b_6 = a = 0.175$ K/MPa (value provided by the material supplier). As consequence, the pressure dependency of the glass transition temperature was found to be $b = 0.161$ K/MPa by using equation 8. The crystallization model parameters are summarized in Table 5.

### 3.2.3 Heat of Crystallization

In the standard Moldflow® solver, the implemented energy equation in terms of temperature *T* reads as follows:

$$\rho c_p \left( \frac{\partial T}{\partial t} + \mathbf{v} \cdot \nabla T \right) = \nabla \cdot (k \nabla T) + \tau : \nabla \mathbf{v} + \zeta T \left( \frac{\partial P}{\partial t} + \mathbf{v} \cdot \nabla P \right) \quad (9)$$

where $\rho$ is the polymer density, *t* is the time, $\mathbf{v}$ is the velocity vector, *P* is the pressure, $\tau$ is the viscous stress tensor, *k* is the polymer thermal conductivity, $c_p$ is the specific heat capacity of the melt and $\zeta = -\frac{1}{\rho} \frac{\partial \rho}{\partial T}$ is the coefficient of volume expansion. Notice that there is not explicit
<table><caption>**Table 6** The specific heat data used in the Moldflow® simulations along with self-defined models</caption><thead><tr><td>Temperature (°C)</td><td>Specific Heat (J kg<sup>-1</sup> K<sup>-1</sup>)</td></tr></thead><tbody><tr><td>32</td><td>1286</td></tr><tr><td>220</td><td>1932</td></tr></tbody></table>

term of heat source. In this work, we model however explicitely the crystallization heat release by adding a heat source term to the right of equation 9 as follows:

$$ \dot{Q} = \rho \Delta H_c \dot{\alpha}(t) \quad (10) $$

where $\Delta H_c$ is the latent heat of crystallization of a material sample with ultimate absolute crystallinity. The numerical implementation of the heat source term via the solver API feature is done by defining a local temperature increment $\Delta T$ as function of the crystallization increment in a given time step:

$$ \Delta T = \dot{\alpha}(t) \Delta t \frac{\Delta H_c}{c_p} \quad (11) $$

where $\Delta t$ is the solver time step. This model supposes that the latent heat is released uniformly along the whole crystallization process. On the other hand, if the local temperature of the melt becomes higher than the melting temperature ($T > T_m(P)$), an instantaneous remelting of the crystals should take place. This endothermic process should lead to a local temperature reduction, modeled as follows: $\Delta T = -\alpha \frac{\Delta H_c}{c_p}$.

The latent heat of crystallization is given in Table 3. As here we implement explicitely the latent heat of crystallization, we have to eliminate the crystallization peak from the specific heat data as typically found in Moldflow® material cards. In this work, we define the temperature-dependent specific heat using the values given in Table 6.

### 3.2.4 Viscosity Model

A modified Cross-WLF model is used to describe the viscosity in this implementation. In this model the shear viscosity depends on crystallization, temperature, shear rate and pressure. The melt shear viscosity $\eta$ is given by:

$$ \eta(T, P, \alpha, \dot{\gamma}) = \frac{\eta_0(T, P, \alpha)}{1 + \left(\frac{\eta_0(T, P, \alpha)\dot{\gamma}}{\tau^*}\right)^{1-\lambda_v}} \quad (12) $$

where $\dot{\gamma}$ is the shear rate, $\tau^*$ is the critical stress level at the transition to shear thinning, $\lambda_v$ is the power law index in the high shear rate regime and $\eta_0$ is the zero shear viscosity defined as:

$$ \eta_0(T, P, \alpha) = D_1 \exp \left( \frac{-A_1(T - T_{ref})}{A_2 + (T - T_{ref})} \right) \vartheta(\alpha) \quad (13) $$

where $T_{ref} = T_g(P)$ with the pressure-dependency defined in equation 7, $D_1$ is theoretically the shear viscosity at $T_{ref}$, and $A_1$ and $A_2$ are data-fitted coefficients approaching the universal constant values of the WLF theory. As for $\vartheta(\alpha)$, it is a function describing the crystallization-dependency of the viscosity and is an extension of Kitano et al.'s relation [81] for concentrated suspension of particles. This function is defined as:

$$ \vartheta(\alpha) = \frac{1}{\left(1 - \frac{\alpha}{A}\right)^B} \quad (14) $$

where A would represent a critical relative crystallinity at which the viscosity of the system tends to infinite (solidification) and B is a data-fitted exponent. The coupling function $\vartheta(\alpha) = \eta_0(T, P, \alpha)/\eta_0(T, P)$ is chosen to increase the effect of crystallization on viscosity at low shear rates as done by Pantani et al. [9].

Be aware that in our implementation of the Cross-WLF model the $A_2$ coefficient does not depend on pressure. This differs from Moldflow® implementation, which includes this pressure dependency and leads to a linear relation between $\eta_0$ and P for all temperatures. However, Rudolph et al. [82] showed a non-linear dependency between $\eta_0$ and P for low temperatures approaching the glass transition temperature of a polycarbonate material. Therefore, in this work, the main connection between shear viscosity and pressure is modeled using the glass transition temperature's dependency.

To identify the parameters of the WLF model, the rheological measurement presented in Figure 4 is used. The viscosity results before the onset of crystallization at 152 °C are utilized to identify the temperature dependency of the shear viscosity by fitting the following linear relation:

$$ \frac{1}{\ln(\eta^*/D_1)} = -\frac{A_2}{A_1} \frac{1}{(T - T_{ref})} - \frac{1}{A_1} \quad (15) $$
<table><caption>**Table 7** The model parameters of the implemented viscosity model</caption><thead><tr><td>Model Parameter</td><td>Value (units)</td></tr></thead><tbody><tr><td>λ<sub>v</sub></td><td>0.16 (-)</td></tr><tr><td>τ<sup>*</sup></td><td>398 000 (Pa)</td></tr><tr><td>D<sub>1</sub></td><td>5.55 × 10<sup>13</sup> (Pa.s)</td></tr><tr><td>A<sub>1</sub></td><td>26.59 (-)</td></tr><tr><td>A<sub>2</sub></td><td>28.87 (K)</td></tr><tr><td>T<sub>g</sub>(0)</td><td>-73.0 (°C) [80]</td></tr><tr><td>b</td><td>0.161 (K/MPa)</td></tr></tbody></table>

where $η^*$ is the magnitude of the measured complex viscosity from the dynamic rheological test. On the other hand, the Cross parameters describing the dependency on shear rate dependency are taken from the Moldflow®'s material database. The model parameters for the Cross-WLF model are summarized in Table 7.

### 3.2.5 PVT Model

Since this work deals with a semi-crystalline thermoplastic material, the Pressure-Volume-Temperature (PVT) relationship is defined by considering explicitely the relative crystallinity $α$. We assume that the material is a simple two-phase system. In this context, we can 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 = \alpha v_s + (1 - \alpha)v_m. \quad (16)$$

The specific volumes $v_m$ and $v_s$ are defined by using the respective empirical Tait equation of state:

$$v_x(T, P) = v_0(T) \left[ 1 - C \ln \left( 1 + \frac{P}{B(T)} \right) \right] \quad (17)$$

with $x = m, s$. In equation 17, $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), \quad (18)$$

<table><caption>**Table 8** The model parameters of the implemented PVT model</caption><thead><tr><td>Model Parameter</td><td>Value (units)</td></tr></thead><tbody><tr><td>b<sub>1m</sub></td><td>8.45 × 10<sup>-4</sup> (m<sup>3</sup>/kg)</td></tr><tr><td>b<sub>2m</sub></td><td>5.65 × 10<sup>-7</sup> (m<sup>3</sup>/kgK)</td></tr><tr><td>b<sub>3m</sub></td><td>1.57 × 10<sup>8</sup> (Pa)</td></tr><tr><td>b<sub>4m</sub></td><td>7.67 × 10<sup>-3</sup> (1/K)</td></tr><tr><td>b<sub>5</sub></td><td>448.15 (K)</td></tr><tr><td>b<sub>1s</sub></td><td>7.45 × 10<sup>-4</sup> (m<sup>3</sup>/kg)</td></tr><tr><td>b<sub>2s</sub></td><td>2.94 × 10<sup>-7</sup> (m<sup>3</sup>/kgK)</td></tr><tr><td>b<sub>3s</sub></td><td>2.88 × 10<sup>8</sup> (Pa)</td></tr><tr><td>b<sub>4s</sub></td><td>4.78 × 10<sup>-3</sup> (1/K)</td></tr></tbody></table>

$$B(T) = b_{3x} \exp(-b_{4x}(T - b_5)) \quad (19)$$

where $b_{1x}$, $b_{2x}$, $b_{3x}$, $b_{4x}$ and $b_5$ are data-fitted coefficients.

In this work, the respective Tait parameters are taken from Moldflow®'s material database and are summarized in Table 8.

### 3.2.6 Solidification Model

A solidification criterion is used in an injection molding simulation to mimic a melt-to-solid transition in the manufacturing process. Numerically it is treated as an abrupt increase of local viscosity to be able to solve the flow problem without shifting the system boundaries. Traditionally, the criterion is defined by a constant no-flow temperature, which should be characteristic of the material. However, the solidification of a semi-crystalline polymer is primarily dependent on its crystallization degree. Therefore, we propose a crystallization-dependent solidification model as follows:

* If $α ≥ A$ then polymer is treated as solid, and

* If $α < A$ then polymer is in melt state.

where $A$ is the same parameter used to describe the dependency of the viscosity on crystallization given in equation 14.

Table 9 presents a summary of the implemented models in comparison to the ones used by the default Moldflow® solver.
<table><caption>**Table 9** 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</caption><thead><tr><td>Model</td><td>Default Moldflow®</td><td>Proposed Implementation</td></tr></thead><tbody><tr><td>Crystallization</td><td>Not included</td><td>Poitou et al. [2] thermo-mechanical model: Quiescent and flow-induced crystallization</td></tr><tr><td>Viscosity</td><td>Cross-WLF model (linear pressure dependency) without crystallization dependency</td><td>Modified Cross-WLF model (non-linear pressure dependency) coupled with crystallization</td></tr><tr><td>PVT</td><td>2-domain Tait model with a transition temperature linearly dependent on pressure</td><td>2-phases mixing law coupled with relative crystallinity</td></tr><tr><td>Heat of Crystallization</td><td>Indirectly considered in the specific heat data (including crystallization peak)</td><td>Heat release due to crystallization and heat removal due to remelting of crystals</td></tr></tbody></table>

## 3.3 Surrogate Model Generation

### 3.3.1 High-fidelity Simulation Model

The simulations are set-up to recreate the injection molding experiments presented in Section 3.1.2. The injection-molded part corresponds to the geometry shown in Figure 5. The high-fidelity simulation emulates the in-mold cooling process of the injection-molded part, starting from the filling of the mold cavity up to the mold opening for part ejection. The thermal boundary conditions are given by a previous thermal simulation of the complete mold within a stable injection-molding cycle. In Moldflow® software the described simulation corresponds to a Cool(FEM) + Fill + Pack analysis. The geometrical model includes cooling channels meshed as beam elements and the part as well as the feed system (nozzle and flange) meshed using tetrahedral elements with 24 layers through the thickness, where the feed system is defined as a hot runner. Lateral and top views of the meshed model are shown in Figure 7a and Figure 7b, respectively. The computation time of a single simulation requires around 110 minutes in a PC with 4.10 GHz processor and 32 GB RAM.

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. Since in this work multiple user-defined models are implemented, the Solver API option should be enabled and the paths of the text files containing the various model parameters should be specified accordingly.

### 3.3.2 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 10 along with the respective labels of the experimental DOE and surrogate model. Since a cooling analysis is performed in the simulations, the inlet cooling temperature $T_{c,in}$ is varied such as $T_{c,in} = T_{mold} + 4°C$ of the $T_{mold}$ set experimentally. Each surrogate model is built by considering the following input variables: Cystallization model parameters $\beta$ and $K_0$, viscosity model parameters $A$ and $B$, and the heat transfer coefficient during the packing phase $HTC_p$. These input variables are summarized in Table 11 along with their lower and upper limits.

Some of these ranges were determined by making some numerical tests while others were
ed crystallization model for injection molding simulation

<table><caption>**Table 10** Process settings of the baseline simulations for the generation of surrogate models</caption><thead><tr><th>Experimental<br>DoE Label</th><th>Surrogate<br>Model Label</th><th>Inlet cooling<br>temperature,<br>(°C)</th><th>T<sub>c,in</sub></th><th>Injection Velocity,<br>V<sub>inj</sub> (cm<sup>3</sup>/s)</th><th>Holding Pressure,<br>P<sub>hold</sub> (MPa)</th></tr></thead><tbody><tr><td>V6</td><td>SM1</td><td>84</td><td></td><td>50</td><td>50</td></tr><tr><td>V14</td><td>SM2</td><td>99</td><td></td><td>30</td><td>50</td></tr><tr><td>V27</td><td>SM3</td><td>114</td><td></td><td>50</td><td>80</td></tr></tbody></table>

<table><caption>**Table 11** Surrogate model variables with respective upper and lower bounds</caption><thead><tr><th>Surrogate Model Variables</th><th>Units</th><th>Lower Bound</th><th>Upper Bound</th></tr></thead><tbody><tr><td>β</td><td>Pa<sup>-1</sup> s<sup>-1</sup></td><td>-10<sup>12</sup></td><td>-10<sup>10</sup></td></tr><tr><td>ln(K<sub>0</sub>)</td><td>-</td><td>23</td><td>83</td></tr><tr><td>A</td><td>-</td><td>0.3</td><td>0.6</td></tr><tr><td>B</td><td>-</td><td>2</td><td>5</td></tr><tr><td>HTC<sub>p</sub></td><td>W m<sup>-2</sup> °C<sup>-1</sup></td><td>500</td><td>5000</td></tr></tbody></table>

<figure><img src="image_9.png" /><figcaption>**Fig. 7** The meshed simulation model including the part (dark green), runner and sprue (light green), cooling channels (blue) and feed system (red)</figcaption></figure>

based on experimental results. For $K_0$, for example, we took as basis the the DSC measurements presented in Section 3.1.1. However, given the experimental limitations for covering a large range of cooling rates, we have to deal with an uncertain identification of $K_0$. In fact, all values within the proposed interval of variation for $K_0$ fit the available experimental data according to the Avrami-Ozawa relation given in equation 4. On top of that, the interval of variation is centered at $ln(K_0) = 53$, which is a value determined by Plummer and Kausch [83] 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 [84], where $A = 0.68$ is specified for smooth spheres and $A \approx 0.44$ for rough compact crystals. As we expect the formation of spherulites and thread-like crystal morphologies during the injection molding process of POM, we propose a rough variation of $\pm 0.15$ around the value given for compact crystals. As for $B$ in the viscosity model, Kitano [81] and Metzner [84] defined it to be equal to 2 according to suspension theory. However, in our work $B$ is varied between 2 and 5 to encompass larger viscosity increases with evolving crystallization, as observed experimentally for POM in Figure 4.
$\beta$ is a parameter that controls the contribution of the flow-induced crystallization to the total relative crystallinity in the model. It turns out that the overall filling simulation is extremely sensitive to this parameter because it can trigger artificially fast solidification of the gate and therefore incomplete cavity filling. The interval of variation for $\beta$ was defined after multiple high-fidelity simulation sensitivity studies and surrogate model tests, guaranteeing a DoE that produces a representative amount of short shots and fully filled parts to efficiently calibrate the material-dependent model parameters.

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

### 3.3.3 Generation Methodology

The training set for the generation of the different surrogate models is sampled using the Latin Hypercube method. 132 simulations were chosen for training and 20 simulations for testing purposes. This choice is postulated to be appropriate as five parameters are being varied in this work in comparison to the six parameters changed in our previous work [3], where an acceptable accuracy for the estimation of the pressure signal was reached using 120 training simulations.

Once the results of the high-fidelity simulation DoE are available, the generation of the surrogate models is done in MATLAB R2019b using proper orthogonal decomposition (POD) of the pressure signal and non-linear regression (NLR) of the POD basis coefficients. This methodology is partially analogous to the POD-NLR method presented in [3]. The following steps describe the process of generating the surrogate models using the POD-NLR technique used in this work:

1. Normalization: The pressure signal results $P_i$ ($i = 1, \dots, S$) are normalized between 0 and 1 as follows:

$$t_{norm} = \frac{t_{original} - t_{start}}{t_{end}} \qquad (20)$$

where $t_{original}$ is the imported unprocessed simulation time, $t_{start}$ represents the time at which the flow front reaches the sensor node producing a non-zero pressure value, and $t_{end}$ is

the time at which the pressure signal becomes zero, i.e. the local pressure equalizes the ambient pressure. In this work, $t_{start}$ is substantially smaller than $t_{end}$, that is why we have neglected it in the denominator of the normalization equation. Subsequently, the simulated pressure is resampled using 1000 points uniformly distributed in the normalized time interval between 0 and 1.

2. Model order reduction: The normalized simulated pressure signal is dimensionally reduced by means of proper orthogonal decomposition (POD). To obtain the POD basis functions, the eigenvalue problem $\mathbf{P}\mathbf{P}^T V_k = \lambda_k V_k$ needs to be first solved, where $\mathbf{P}$ is a snapshot matrix of different pressure signals and $V_1, V_2, V_3, \dots$ are the eigenvectors associated to the eigenvalues $\lambda_1, \lambda_2, \lambda_3, \dots$. The POD basis functions $\phi_k$ can be then determined by normalizing the eigenvectors by $\phi_k = V_k/\lambda_k$. These basis functions create an orthogonal system that captures the dominant energy modes within the data and their coefficients are determined by projecting the original data onto these functions. The reconstruction of the pressure signal can be obtained by:

$$P_{reconstructed}^{(s)} = \sum_{k=1}^{n} \Gamma_{k}^{(s)} \phi_{k} \qquad (21)$$

where $s$ is the high-fidelity simulation number, $n$ is the number of modes or basis functions obtained according to a specified error value and $\Gamma_k^{(s)}$ is the POD basis coefficient for a specific mode $k$ and a simulation $s$. The truncation criterion is done according to an error value of $\epsilon = 5 \times 10^{-4}$, which cut off the number of modes once the eigenvalue condition $\frac{\lambda}{\lambda_{max}} \ge \epsilon$ is not more satisfied.

3. Regression: A least-squares regression of a second-order polynomial is used to correlate the POD basis coefficients $\Gamma$, and the time shift values $t_{start}$ and $t_{end}$ with the metamodel variables. The regression equations read as follows, where the implicit Einstein summation convention is used for indexes $i$ and $j$:

$$\Gamma_{n}^{(s)} = a^{(n)} + b_{i}^{(n)} X_{i}^{(s)} + c_{ij}^{(n)} X_{i}^{(s)} X_{j}^{(s)} \qquad (22)$$

$$t_{start}^{(s)} = d + e_{i}X_{i}^{(s)} + f_{ij}X_{i}^{(s)}X_{j}^{(s)} \qquad (23)$$
$$t_{end}^{(s)} = g + h_i X_i^{(s)} + z_{ij} X_i^{(s)} X_j^{(s)} \quad (24)$$

where $X_i$ and $X_j$ are the surrogate model input variables and $a, b_i, c_{ij}, d, e_i, f_{ij}, g, h_i, z_{ij}$ are the surrogate model parameters.

### 3.3.4 Performance of the Surrogate Models

The performance of each surrogate model is assessed by measuring its capacity to recreate the pressure signal given by the high-fidelity simulation at one sensor location. In this section, we present some predicted pressure results from both the training and testing DoE sets as well as the error metrics for the estimation of the POD basis coefficients and time shifts. In other words, we measure the quality of the fitted polynomials defined in equations 22, 23 and 24 and their predictive capability.

One of the error metric in this work is the normalized root mean squared error (RMSE) that uses the min-max normalization method to facilitate the comparison between various surrogate modeling techniques:

$$ \text{Normalized RMSE}(Y, \hat{Y}) = \frac{\text{RMSE}(Y, \hat{Y})}{\frac{Y_{\max} - Y_{\min}}{\sqrt{\frac{1}{n} \sum_{i=1}^{n} (Y^i - \hat{Y}^i)^2}}} = \frac{Y_{\max} - Y_{\min}}{Y_{\max} - Y_{\min}} \quad (25) $$

where $Y$ and $\hat{Y}$ represents the target and estimated output, respectively. In this work, the estimated output is the surrogate model response.

The second metric is 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 unseen samples are likely to be predicted by the model. $R^2$ ranges between $-\infty < R^2 \le 1$, where 1.0 is the best score. $R^2$ is calculated as follows:

$$R^2(Y, \hat{Y}) = 1 - \frac{\sum_{i=1}^{n} (Y^i - \hat{Y}^i)^2}{\sum_{i=1}^{n} (Y^i - \bar{Y})^2} \quad (26)$$

where $\bar{Y} = \frac{1}{n} \sum_{i=1}^{n} Y^i$ is the mean value of $Y$.

To ease the analysis of the previously-mentioned metrics, the POD basis functions associated to each surrogate model are presented in Figure 8. For the molding conditions V6 and V14,

the simulation output is reduced using five basis functions. For the molding condition V27, data reduction required only three basis functions by using the same truncation error $\epsilon$. This fact means that the simulated pressure signal for the process settings V27 is less fluctuating than the ones for the other process conditions.

Starting with the surrogate model SM1, Figure 9a presents the pressure results of four training simulations and the approximated results given by the trained surrogate model. Even considering training data, there still exists some discrepancy between the target simulation and the surrogate model predictions, especially in the case of the green curve. These differences are somehow explainable since the metamodel is trained on a reduced form of the predicted signal and not directly on the full resolution simulated pressure response. On the other hand, the pressure signal depicted in the green curve is more wavy than the other three curves making it more challenging to approximate using the same number of basis functions. Figure 9b shows the pressure signals of four simulations used to test the capability of the surrogate model to predict unseen data. The surrogate model is able to recreate reasonably well the form of the different pressure signals but struggles to accurately predict $t_{end}$ as seen in the case colored in dark blue.

The error metrics of the surrogate model SM1 for estimating the five POD basis coefficients $\Gamma$ and the two time shift parameters are presented in Figure 10. Figures 10a and 10b show the normalized RMSE and $R^2$ score, respectively, differentiating between training and testing scenarios. As expected, error metrics for training are lower than for testing data sets. The normalized RMSE for estimating the POD basis coefficients $\Gamma_1$ and $\Gamma_4$ as well as the time shift $t_{start}$ doubles when passing from training to testing data sets. Overall, the surrogate model SM1 delivers a normalized RMSE lower than 26 %, which is considered as acceptable for the purposes in this work. By analyzing the coefficient of determination, the estimation of the POD basis coefficients $\Gamma_3$ and $\Gamma_4$ produce a negative $R^2$ score for the testing data set, which points out that the SM1 should not be used as single regressor for those parameters. Even if the $R^2$ values of the different parameters (basis coefficients and time shifts) are low, notice that
<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 8** The POD basis functions of the three generated surrogate models</figcaption></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /><figcaption>**Fig. 9** 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</figcaption></figure>

the reconstruction of the simulated pressure signals (Figure 9) is accurate enough for the purposes of this work. In fact, the reconstruction exhibits a mean averaged error of 2.5 MPa and 3.7 MPa for the training and testing data sets, respectively.

Moving on to the surrogate model SM2, Figure 11 shows high-fidelity simulation results in front of the surrogate model estimations for both training and testing data sets. Some minor differences between high-fidelity simulation and surrogate model are visible. In the training set, the maximum pressure of the simulation in turquoise is under-predicted whereas the end times $t_{end}$ for the
ed crystallization model for injection molding simulation

<figure><img src="image_15.png" /><figcaption>**Fig. 10** The error metrics of the calculated POD basis coefficients $\Gamma$ and time shifts $t_{start}$ and $t_{end}$ using the SM1 surrogate model when using input data from the training and testing sets</figcaption></figure>

simulations in green and pink are over-predicted. In the testing set, SM2 model performs reasonably well as seen in Figure 11b. This is confirmed by a normalized RMSE lower than 0.2 for all POD basis coefficients as shown in Figure 12. Additionally, there is a relative small difference of performance metrics between the training and testing data sets. In this case, the reconstruction of the simulated pressure signals using SM2 shows a mean average error of 2.5 MPa and 2.3 MPa for the training and testing data sets, respectively.

The SM3 model exhibits also a remarkable performance as observed in Figures 13a and 13b for training and testing data sets, respectively.

<figure><img src="image_16.png" /><figcaption>**Fig. 11** 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</figcaption></figure>

The most noticeable deviation occurs in the training set, where the maximum pressures of the high-fidelity simulations represented in green and turquoise (Figure 13a) are under-predicted by the metamodel. Interestingly, the model-order reduction in this case required only three POD basis functions (see Figure 8c) to satisfy the truncation error of $\epsilon = 5 \times 10^{-4}$. For this reason, Figure 14 presents the error metrics for only three $\Gamma$ POD basis coefficients. In the training set, the normalized RMSE for all parameters is lower than 12 %. However, the prognosis of the third POD basis coefficient $\Gamma_3$ as well as $t_{start}$ is inaccurate in
<figure><img src="image_17.png" /><figcaption>Fig. 12 The error metrics of the calculated POD basis coefficients Γ and time shifts t$_{start}$ and t$_{end}$ using the SM2 surrogate model when using input data from the training and testing sets</figcaption></figure>

the testing set (negative R²), pointing out a deficiency of the chosen polynomial regressor for those parameters. Nonetheless, the reconstruction of the simulated pressure signals using SM3 is acceptable as far as it displays a mean average error of 3.6 MPa and 4.6 MPa for the training and testing data sets, respectively.

## 3.4 Surrogate Model-based Parameter Calibration

### 3.4.1 Calibration Algorithm

To identify the five modeling parameters presented in Table 11, a multi-objective optimization

<figure><img src="image_18.png" /><figcaption>Fig. 13 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</figcaption></figure>

routine is performed using the `lsqnonlin` built-in MATLAB function from the optimization toolbox. The experimental pressure signals are used as reference for defining the objective functions. The calibration utilizes all three surrogate models to obtain one single set of optimized parameters using the following formulation:

$$ \text{find} \quad V = (\beta, K_0, A, B, HTC_p) $$

$$ \begin{array}{ll} \text{argmin} & | \hat{Y}_{max}^l - Y_{max}^{exp,k} | \quad (l = 1, 2, 3) \\ & | \hat{t}_{end}^l - t_{end}^{exp,k} | \quad \& k = 6, 14, 27) \end{array} $$

$$ \text{s.t.} \quad -10^{12} \le \beta \le -10^{10}, $$
ed crystallization model for injection molding simulation

<figure><img src="image_19.png" /><figcaption>**Fig. 14** The error metrics of the calculated POD basis coefficients $\Gamma$ and time shifts $t_{start}$ and $t_{end}$ using the SM3 surrogate model when using input data from the training and testing sets</figcaption></figure>

$$ 
\begin{array}{l}
23 \le \ln(K_0) \le 83, \\
0.3 \le A \le 0.6, \\
2 \le B \le 5, \\
500 \le HTC_p \le 5000.
\end{array}
$$

The optimization problem contains in total six objective functions. Three of them minimizes the difference between experiment and surrogate model for the maximal pressure value in three different processing conditions. The experimental maximal pressure is $Y_{max}^{exp,k}$, where $k = 6, 14, 27$ represents the experimental DoE labels V6, V14, and V27. On the other hand, the maximal pressure

<figure><img src="image_20.png" /><figcaption>**Fig. 15** The time intercept of a tangent straight line fitted to the pressure signal decay at end of packing for the P2 sensor signal of the DOE level V6</figcaption></figure>

<table><caption>**Table 12** The experimentally-determined parameters used to define the six objective functions used in the optimization routine to identify the modeling parameters</caption><thead><tr><td>k</td><td>Y<sub>max</sub><sup>exp,k</sup> (MPa)</td><td>t<sub>end</sub><sup>exp,k</sup> (s)</td></tr></thead><tbody><tr><td>6</td><td>38.1</td><td>12.8</td></tr><tr><td>14</td><td>35.2</td><td>15.2</td></tr><tr><td>27</td><td>63.6</td><td>20.5</td></tr></tbody></table>

determined by the surrogate model $l$ is $\hat{Y}_{max}^l$. The remaining objective functions minimize the difference between predicted end time $\hat{t}_{end}^l$, calculated according to equation 24, and the experimentally-determined end time $t_{end}^{exp,k}$. The experimental end time value is obtained as the time intercept of a tangent straight line fitted to the pressure signal decay at end of packing. A visual example of the end time determination of V6 is given in Figure 15. Table 12 summarizes the experimental reference values used to define the six objective functions.

### 3.4.2 Calibration Results

For calibrating the material models implemented in Moldflow®, we can use now the three generated surrogate models by identifying $\beta, K_0, A, B, HTC_p$ according to the algorithm described in Section 3.4.1. The optimization routine necessitated 15 iterations to reach a local minimum as shown in Figure 16. The optimization
<figure><img src="image_21.png" /><figcaption>**Fig. 16** The norm of the residuals as a function of the number of iterations performed by *lsqnonlin* during the optimization routine</figcaption></figure>

<table><caption>**Table 13** Identified model parameters after multi-objective optimization</caption><thead><tr><td>Model Parameter</td><td>Units</td><td>Identified Value</td></tr></thead><tbody><tr><td>β</td><td>Pa<sup>-1</sup> s<sup>-1</sup></td><td>-10<sup>11.37</sup></td></tr><tr><td>ln(K<sub>0</sub>)</td><td>-</td><td>83.00</td></tr><tr><td>A</td><td>-</td><td>0.35</td></tr><tr><td>B</td><td>-</td><td>4.96</td></tr><tr><td>HTC<sub>p</sub></td><td>W m<sup>-2</sup> °C<sup>-1</sup></td><td>500</td></tr></tbody></table>

routine stopped because the function tolerance fell below $10^{-4}$. To check if there exist other local minima in this multi-objective optimization and finally obtain the global minimum, the built-in *MultiStart* MATLAB function is used to find other possible local minima to the optimization problem by starting from different initial points. We did not find a different combination of optimized parameters suggesting that we found the global minimum from the first try. The optimized parameters are given in Table 13.

The identified solidification criterion parameter $A = 0.35$ seems reasonable since it depends theoretical on the geometry of the suspending solid crystallites. According to Kennedy and Zheng [1] and Metzner [84], $A$ is approximately equal to 0.44 for rough compact crystals and since during injection molding fibrillar crystal structures could develop due to intense shear strain, an $A < 0.44$ is plausible. The optimized value for the exponent parameter $B$ is nearly equal to the upper bound of the searching interval (i.e. $B =$

5), which underlines the strong coupling between evolving crystallization and effective shear viscosity. This correlates well with the results of the rheological test in Figure 4, where the viscosity increases abruptly after the onset of crystallization. On the other hand, the identified value of the FIC parameter $\beta$ is not straightforward to assess but $\beta = 11.37$ translates into a lower contribution of FIC to the total relative crystallinity. The optimized value for $\ln(K_0)$ corresponds to the upper bound of the considered searching interval, which corroborates the renowned high crytstallization kinetics of POM, even at quiescent conditions. The identified value is however greater than the one determined by Plummer and Kausch [83, 85] equal to 53. Finally, the heat transfer coefficient during packing was found equal to 500 W m$^{-2}$ °C$^{-1}$, which is a result prone to critic and discussion. In fact, this low value is typically used to describe the heat transfer condition under detachment between mold and part, which hardly occurs during the entire packaging phase. The heat transfer coefficient is a simulation parameter controlling the evolution of the temperature field in the polymer domain and therefore can be seen also as numerical artifact for tuning the crystallization kinetics and solidification in the simulation. For reference, in Moldflow® the default value for the global heat transfer coefficient during packing is 2500 W m$^{-2}$ °C$^{-1}$. This result points out the need for implementing a local definition of the heat transfer coefficient, probably depending on the instantaneous local pressure, in the commercial codes of injection molding simulation.

# 4 Calibrated High-fidelity Simulation with Crystallization Model

In the following, the high-fidelity simulation including the calibrated flow-induced crystallization model is compared with the default Moldflow® simulation that does not take crystallization explicitly into account.

## 4.1 Pressure Results

As shown in Table 10, each surrogate model was generated using a fixing set of processing conditions and mimics one particular level of the
ed crystallization model for injection molding simulation

<figure><img src="image_22.png" /><figcaption>**Fig. 17** The pressure signals for the V6 processing condition obtained using the simulation with the calibrated models in comparison to the ones obtained by the default Moldflow® simulation along with the corresponding experimental results for this condition</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>**Fig. 18** The pressure signals for the V14 processing condition obtained using the simulation with the calibrated models in comparison to the ones obtained by the default Moldflow® simulation along with the corresponding experimental results for this condition</figcaption></figure>
<figure><img src="image_24.png" /><figcaption>**Fig. 19** The pressure signals for the V27 processing condition obtained using the simulation with the calibrated models in comparison to the ones obtained by the default Moldflow® simulation along with the corresponding experimental results for this condition</figcaption></figure>

experimental design of experiments. On the other hand, the pressure signal at the sensor located after the gate (P2) was used to identify the modeling parameters given in Table 13. Therefore, a natural way to assess the performance of the implemented models is to compare the predicted pressure signals by the standard Moldflow® simulation (Default simulation) with the high-fidelity simulation including the calibrated crystallization model (Simulation with calibrated crystallization model). Figures 17, 18 and 19 present the experimental and simulated pressure signals at the sensor locations P1, P2 and P3 (refer to Figure 5) for the processing conditions V6, V14 and V27, respectively.

The high-fidelity simulation including the flow-induced crystallization model outperforms the default high-fidelity simulation for all processing conditions and improves particularly the pressure estimation at sensor locations P2 and P3. Interestingly the proposed calibrated simulation is able to recreate the pressure plateau just after the switch-over as observed in Figures 18a and 18b at P1 and P2, respectively, which leads to a pressure signal delay at P3 as shown in Figure 18c for the processing condition V14.

To quantify the improvement in pressure predictions, the normalized RMSE are calculated for the different simulations approaches and plotted in Figures 20a, 20b and 20c for the processing conditions V6, V14 and V27, respectively. As expected, the prediction error decreases significantly by including the calibrated crystallization model along with the crystallization-dependent viscosity, PVT and solidification models. Major improvements are noticeable at sensor locations P2 and P3 for the processing conditions V6 and V14, where the error decreases three folds at P2 and four times at P3. A more accurate estimation of the pressure field should in principle lead to higher accuracy when predicting part shrinkage and warpage. For the processing condition V27, the improvement of accuracy is less pronounced in comparison to the other conditions, but the error still decreased by approximately 5 %.

The previous comparison appears promising, but is biased as far as the implemented models were calibrated using the measured pressure signals for those processing conditions. Therefore, Figure 21 presents the averaged normalized RMSE at P2 for all experimental processing conditions
ed crystallization model for injection molding simulation

<figure><img src="image_25.png" /><figcaption>**Fig. 20** 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</figcaption></figure>

(all levels of the design of experiments) while using standard Moldflow® models and Moldflow® with the calibrated crystallization model.

According to Figure 21, for the majority of the different DoE levels, the estimation error is decreased by more than half with respect to the default high-fidelity simulation. Some exceptions are for example the DoE levels V7, V8 and V9. To analyze the probable reasons behind this discrepancy, Figure 22 presents the experimental pressure signal at sensor location P2 in comparison to the simulated ones (default and with calibrated crystallization models) for the mentioned DoE levels. Both simulations fail to recreate the experimental pressure signals in terms of maximum value and signal decay. In fact, all experimental signals exhibit an unusual inflection point around 10 s. This could be due to a partial remelting of the solidified material at the gate that would lead to an increased melt flux entering the cavity. Even though our proposed modeling includes the effect of remelting, our current calibrated simulation is not able to recreate such bump in pressure signal. This fact reveals that some model assumptions need to be revisited or that the models are not fully calibrated yet. Please notice however that in Figure 22b, our proposed simulation shows an inflection point around 7 s. On the other hand, the high-fidelity simulation with calibrated crystallization models is able to predict with higher accuracy the maximal pressure before switch-over and also the pressure signal response just after commuting to pressure control.

## 4.2 Predictions of short shots

As mentioned in Section 3.1.2, some experimental DoE levels did not produce a fully filled part, displaying what in the milieu is known as short shot. Figure 23a summarizes the filling status of the 27 experimental DoE levels. One of the claimed advantages of implementing crystallization-dependent material models in injection molding simulation is the ability to predict if a mold cavity can be filled completely or not. The default Moldflow® simulation is not able to predict the short shots as seen in Figure 23b. However, our high-fidelity simulation with calibrated crystallization models makes possible to predict most of these short shots as shown in Figure 23c. In total, four out of the five short
<figure><img src="image_26.png" /><figcaption>**Fig. 21** Averaged normalized RMSE of the pressure estimations at sensor location P2 using the high-fidelity simulation with the calibrated crystallization models and the default Moldflow® simulation for all 27 experimental processing conditions used for injecting the half-length 3 mm thick part</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>**Fig. 22** Comparison of pressure signals at position P2 for processing condition V7, V8 and V9: Experimental, standard Moldflow® simulation and Moldflow® simulation with calibrated crystallization models</figcaption></figure>

shots are identified along with one wrong prediction for V3 (≈ 4 mm of flow length not filled). For the V19 condition, the simulation showed a complete cavity filling, but experimentally it shorted 7 mm before end of cavity. With our proposed models, the estimated flow length is accurate up to ±15 mm. For example, the actual parts produced with the processing conditions V10 shorted 10 mm before end of cavity whereas in the simulation the leftover unfilled length was around 5 mm.

# 5 Conclusion

We successfully implemented a flow-induced crystallization model coupled with viscosity, PVT and solidification models in a commercial software for injection molding simulation (Autodesk® Moldflow® Insight). This modeling framework is applicable to any semi-crystalline thermoplastic polymer. For a given material, the identification of the material dependent parameters requires initially standard experimental characterization of the quiescent crystallization (fundamentally based on DSC), which reduces the number of unknown
ed crystallization model for injection molding simulation

<figure><img src="image_28.png" /><figcaption>**Fig. 23** The experimental fill results (a) in comparison to the ones predicted using the (b) default simulation, (b) simulation with the calibrated crystallization models</figcaption></figure>

parameters to 5. The lacking 5 unknown material-
dependent parameters are calibrated using a sur-
rogate model-based optimization taking as ref-
erence in-mold pressure sensor signals from an
actual injection molding process. In this work,
we applied the identification methodology to an
industrial grade of an unreinforced polyoxymethy-
lene (POM) homopolymer. This novel data-driven
calibration method required the generation of a
surrogate model (trained with 132 high-fidelity
simulations) and a multi-objective optimization
of the respective surrogate model. This last
step needed only 15 iterations (less than 1 s)
to solve the minimization problem for identify-
ing the 5 material-dependent parameters of the
flow-induced crystallization model under actual
injection-molding conditions. In sum, the effort
of the proposed method seems quite lower in
comparison with the costs (and time) associated
to the extensive characterization used for iden-
tifying phenomenological crystallization models
under laboratory equipment conditions.

The high-fidelity simulation with calibrated
material-dependent model parameters produced
a remarkable increase in accuracy of the pres-
sure field during injection molding. Concretely,
the error in the estimation of the pressure signals
at the sensor locations for a variety of process-
ing conditions decreased by more than half when
comparing the simulation including flow-induced
crystallization modeling in front of the default

simulation (without explicit crystallization model-
ing). This result highlights particularly the impor-
tant role of the crystallization-viscosity coupling,
which influences the pressure evolution during a
mold filling. In addition, the simulation with crys-
tallization modeling was able to predict 80% of the
short shots observed in an actual injection molding
DoE (5 from 27), whereas the default simula-
tion predicted not one. This outcome underlines
the importance of implementing a crystallization
model in injection molding simulation as well as
the use of a crystallization-dependent solidifica-
tion criterion instead of the today conventional
no-flow temperature.

As an outlook into the topic of surrogate mod-
eling, the approach presented in this work could be
extended to emulate a full time and space model
of an injection molding cycle. Such emulator could
be used to visualize and assess virtually an injec-
tion molding process in real time while varying
input parameters, also paving the way to perform
efficient uncertainty quantification. On the other
hand, those models constitute an opening to the
industrial deployment of digital twins in plastics
processing.

References

[1] Kennedy, P. & Zheng, R. *Flow analysis of injection molds* 2nd edition edn (Hanser Publishers, Munich ; Cincinnati, 2013).
[2] Poitou, A., Ammar, A., Marco, Y., Chevalier, L. & Chaouche, M. Crystallization of polymers under strain: from molecular properties to macroscopic models. *Computer Methods in Applied Mechanics and Engineering* **192** (28-30), 3245–3264 (2003). URL https://linkinghub.elsevier.com/retrieve/pii/S0045782503003499. https://doi.org/10.1016/S0045-7825(03)00349-9.

[3] Saad, S., Sinha, A., Cruz, C., Rgnier, G. & Ammar, A. Towards an accurate pressure estimation in injection molding simulation using surrogate modeling. *International Journal of Material Forming* **15** (6) (2022). URL https://link.springer.com/10.1007/s12289-022-01717-0. https://doi.org/10.1007/s12289-022-01717-0.

[4] Hu, W. *Polymer Physics* (Springer Vienna, Vienna, 2013). URL http://link.springer.com/10.1007/978-3-7091-0670-9. DOI: 10.1007/978-3-7091-0670-9.

[5] Liang, G., Bao, S. & Zhu, F. in *Theoretical Aspects of Polymer Crystallization in Multiphase Systems* 17–48 (Elsevier, 2018). URL https://linkinghub.elsevier.com/retrieve/pii/B9780128094532000025. DOI: 10.1016/B978-0-12-809453-2.00002-5.

[6] Zuidema, H. *Flow induced crystallization of polymers: application to injection moulding*. PhD thesis, Technische Universiteit Eindhoven, Eindhoven (2000). OCLC: 248341950.

[7] Boutaous, M., Bourgin, P. & Zinet, M. Thermally and flow induced crystallization of polymers at low shear rate. *Journal of Non-Newtonian Fluid Mechanics* **165** (5-6), 227–237 (2010). URL https://linkinghub.elsevier.com/retrieve/pii/S0377025709002353. https://doi.org/10.1016/j.jnnfm.2009.12.005.

[8] Zhang, M., Guo, B.-H. & Xu, J. A Review on Polymer Crystallization Theories. *Crystals* **7** (1), 4 (2016). URL http://www.mdpi.com/2073-4352/7/1/4. https://doi.org/10.3390/cryst7010004.

[9] Pantani, R., Coccorullo, I., Speranza, V. & Titomanlio, G. Modeling of morphology

evolution in the injection molding process of thermoplastic polymers. *Progress in Polymer Science* **30** (12), 1185–1222 (2005). URL https://linkinghub.elsevier.com/retrieve/pii/S0079670005001176. https://doi.org/10.1016/j.progpolymsci.2005.09.001.

[10] Pantani, R., De Santis, F., Speranza, V. & Titomanlio, G. Analysis of flow induced crystallization through molecular stretch. *Polymer* **105**, 187–194 (2016). URL https://linkinghub.elsevier.com/retrieve/pii/S0032386116309375. https://doi.org/10.1016/j.polymer.2016.10.026.

[11] Kolmogotov, A. N. A Statistical Theory for the Recrystallisation of Metals. *Izvestiya Akademii Nauk SSSR, Neorganicheskie Materialy* **3**, 355–359 (1937).

[12] Avrami, M. Kinetics of Phase Change. I General Theory. *The Journal of Chemical Physics* **7** (12), 1103–1112 (1939). URL http://aip.scitation.org/doi/10.1063/1.1750380. https://doi.org/10.1063/1.1750380.

[13] Evans, U. R. The laws of expanding circles and spheres in relation to the lateral growth of surface films and the grain-size of metals. *Transactions of the Faraday Society* **41**, 365 (1945). URL http://xlink.rsc.org/?DOI=tf9454100365. https://doi.org/10.1039/tf9454100365.

[14] Lauritzen, J. I. & Hoffman, J. D. Theory of formation of polymer crystals with folded chains in dilute solution. *Journal of Research of the National Bureau of Standards Section A: Physics and Chemistry* **64A** (1), 73 (1960). URL https://nvlpubs.nist.gov/nistpubs/jres/64A/jresv64An1p73_A1b.pdf. https://doi.org/10.6028/jres.064A.007.

[15] Ozawa, T. Kinetics of non-isothermal crystallization. *Polymer* **12** (3), 150–158 (1971). URL https://linkinghub.elsevier.com/retrieve/pii/0032386171900413. https://doi.org/10.1016/0032-3861(71)90041-3.

[16] Nakamura, K., Watanabe, T., Katayama, K. & Amano, T. Some aspects of nonisothermal crystallization of polymers. I. Relationship
ed crystallization model for injection molding simulation

between crystallization temperature, crys-
tallinity, and cooling conditions. *Journal
of Applied Polymer Science* **16** (5), 1077–
1091 (1972). URL http://doi.wiley.com/10.
1002/app.1972.070160503. https://doi.org/
10.1002/app.1972.070160503.

[17] Ziabicki, A. Crystallization of polymers in variable external conditions: 1. General equations. *Colloid & Polymer Science* **274** (3), 209–217 (1996). URL http://link.springer.com/10.1007/BF00665637. https://doi.org/10.1007/BF00665637.

[18] Koscher, E. & Fulchiron, R. Influence of shear on polypropylene crystallization: morphology development and kinetics. *Polymer* **43** (25), 6931–6942 (2002). URL https://www.sciencedirect.com/science/article/pii/S0032386102006286.

[19] Schneider, W., Kppl, A. & Berger, J. Non Isothermal Crystallization Crystallization of Polymers: System of Rate Equations. *International Polymer Processing* **2** (3-4), 151–154 (1988). URL http://www.hanser-elibrary.com/doi/10.3139/217.880150. https://doi.org/10.3139/217.880150.

[20] Eder, G., Janeschitz-Kriegl, H. & Liedauer, S. Crystallization processes in quiescent and moving polymer melts under heat transfer conditions. *Progress in Polymer Science* **15** (4), 629–714 (1990). URL https://linkinghub.elsevier.com/retrieve/pii/007967009090008O. https://doi.org/10.1016/0079-6700(90)90008-O.

[21] Doufas, A. K., McHugh, A. J. & Miller, C. Simulation of melt spinning including flow-induced crystallization: Part I. Model development and predictions. *Journal of Non-Newtonian Fluid Mechanics* **92** (1), 27–66 (2000). URL https://www.sciencedirect.com/science/article/pii/S0377025700000884.

[22] Zuidema, H., Peters, G. W. & Meijer, H. E. Development and validation of a recoverable strain-based model for flow-induced crystallization of polymers. *Macromolecular theory and simulations* **10** (5), 447–460 (2001). URL

https://onlinelibrary.wiley.com/doi/abs/
10.1002/1521-3919(20010601)10:5%3C447::
AID-MATS447%3E3.0.CO;2-C.

[23] Guo, J. & Narh, K. A. Simplified model of stress-induced crystallization kinetics of polymers. *Advances in Polymer Technology* **21** (3), 214–222 (2002). URL http://doi.wiley.com/10.1002/adv.10022. https://doi.org/10.1002/adv.10022.

[24] Kim, K. H., Isayev, A., Kwon, K. & van Sweden, C. Modeling and experimental study of birefringence in injection molding of semicrystalline polymers. *Polymer* **46** (12), 4183–4203 (2005). URL https://linkinghub.elsevier.com/retrieve/pii/S0032386105002077. https://doi.org/10.1016/j.polymer.2005.02.057.

[25] Coppola, S., Balzano, L., Gioffredi, E., Maffettone, P. L. & Grizzuti, N. Effects of the degree of undercooling on flow induced crystallization in polymer melts. *Polymer* **45** (10), 3249–3256 (2004). URL https://linkinghub.elsevier.com/retrieve/pii/S0032386104002733. https://doi.org/10.1016/j.polymer.2004.03.049.

[26] Titomanlio, G. & Lamberti, G. Modeling flow induced crystallization in film casting of polypropylene. *Rheologica Acta* **43** (2), 146–158 (2004). URL http://link.springer.com/10.1007/s00397-003-0329-4. https://doi.org/10.1007/s00397-003-0329-4.

[27] Zinet, M., El Otmani, R., Boutaous, M. & Chantrenne, P. Numerical modeling of nonisothermal polymer crystallization kinetics: Flow and thermal effects. *Polymer Engineering & Science* **50** (10), 2044–2059 (2010). URL http://doi.wiley.com/10.1002/pen.21733. https://doi.org/10.1002/pen.21733.

[28] Roozemond, P. C., van Drongelen, M., Ma, Z., Hulsen, M. A. & Peters, G. W. Modeling flow-induced crystallization in isotactic polypropylene at high shear rates. *Journal of rheology* **59** (3), 613–642 (2015). URL https://sor.scitation.org/doi/abs/10.1122/1.4913696.
[29] Spina, R., Spekowius, M. & Hopmann, C. Multiphysics simulation of thermoplastic polymer crystallization. *Materials & Design* **95**, 455-469 (2016). URL https://linkinghub.elsevier.com/retrieve/pii/S026412751630123X. https://doi.org/10.1016/j.matdes.2016.01.123.

[30] Spina, R., Spekowius, M. & Hopmann, C. Simulation of crystallization of iso-tactic polypropylene with different shear regimes. *Thermochimica Acta* **659**, 44-54 (2018). URL https://linkinghub.elsevier.com/retrieve/pii/S0040603117302782. https://doi.org/10.1016/j.tca.2017.10.023.

[31] Kulkarni, J. A. & Beris, A. N. Lattice-based simulations of chain conformations in semi-crystalline polymers with application to flow-induced crystallization. *Journal of Non-Newtonian Fluid Mechanics* **82** (2-3), 331-366 (1999). URL https://linkinghub.elsevier.com/retrieve/pii/S0377025798001712. https://doi.org/10.1016/S0377-0257(98)00171-2.

[32] Tanner, R. I. & Qi, F. A comparison of some models for describing polymer crystallization at low deformation rates. *Journal of Non-Newtonian Fluid Mechanics* **127** (2-3), 131-141 (2005). URL https://linkinghub.elsevier.com/retrieve/pii/S0377025705000583. https://doi.org/10.1016/j.jnnfm.2005.02.005.

[33] Brahmia, N. *Contribution la modification de la cristallisation des polymères sous cisaillement: application linjection des polymères semi-cristallins*. Ph.D. thesis, PhD Thesis INSA de Lyon, France (2007). URL http://theses.insa-lyon.fr/publication/2007ISAL0069/these.pdf.

[34] Mu, Y., Zhao, G., Chen, A. & Wu, X. Numerical investigation of the thermally and flow induced crystallization behavior of semi-crystalline polymers by using finite elementfinite difference method. *Computers & Chemical Engineering* **46**, 190-204 (2012). URL https://linkinghub.elsevier.com/retrieve/pii/S0098135412002104. https://doi.org/10.1016/j.compchemeng.2012.06.026.

[35] Mu, Y., Zhao, G., Chen, A., Dong, G. & Li, S. Numerical investigation of the crystallization and orientation behavior in polymer processing with a two-phase model. *Computers & Chemical Engineering* **63**, 91-107 (2014). URL https://linkinghub.elsevier.com/retrieve/pii/S0098135413003918. https://doi.org/10.1016/j.compchemeng.2013.12.010.

[36] Haas, T. W. & Maxwell, B. Effects of shear stress on the crystallization of linear polyethylene and polybutene-1. *Polymer Engineering and Science* **9** (4), 225-241 (1969). URL https://onlinelibrary.wiley.com/doi/10.1002/pen.760090402. https://doi.org/10.1002/pen.760090402.

[37] Titomanlio, G., Speranza, V. & Brucato, V. On the Simulation of Thermoplastic Injection Moulding Process: II Relevance of Interaction Between Flow and Crystallization. *International Polymer Processing* **12** (1), 45-53 (1997). URL https://www.degruyter.com/document/doi/10.3139/217.970045/html. https://doi.org/10.3139/217.970045.

[38] Pantani, R., Speranza, V. & Titomanlio, G. Effect of flow-induced crystallization on the distribution of spherulite dimensions along cross section of injection molded parts. *European Polymer Journal* **97**, 220-229 (2017). URL https://linkinghub.elsevier.com/retrieve/pii/S0014305717311618. https://doi.org/10.1016/j.eurpolymj.2017.10.012.

[39] Poitou, A. & Ammar, A. Cristallisation induite par coulement ou par deformation d'un polymère une approche thermodynamique. *Comptes Rendus de l'Académie des Sciences - Series IIB - Mechanics* **329** (1), 5-11 (2001). URL https://linkinghub.elsevier.com/retrieve/pii/S1620774200012903. https://doi.org/10.1016/S1620-7742(00)01290-3.

[40] Kim, S. K. & Jeong, A. Numerical Simulation of Crystal Growth in Injection Molded Thermoplastics based on Monte Carlo Method with Shear Rate Tracking. *International Journal of Precision Engineering and Manufacturing* **20** (4), 641-650.
ed crystallization model for injection molding simulation

(2019). URL http://link.springer.com/10.1007/s12541-019-00089-x. https://doi.org/10.1007/s12541-019-00089-x.

[41] van Meerveld, J., Htter, M. & Peters, G. W. Continuum model for the simulation of fiber spinning, with quiescent and flow-induced crystallization. *Journal of Non-Newtonian Fluid Mechanics* **150** (2-3), 177–195 (2008). URL https://linkinghub.elsevier.com/retrieve/pii/S0377025707002340. https://doi.org/10.1016/j.jnnfm.2007.10.020.

[42] Zinet, M., Refaa, Z., Boutaous, M., Xin, S. & Bourgin, P. Thermophysical Characterization and Crystallization Kinetics of Semi-Crystalline Polymers. *Journal of Modern Physics* **04** (07), 28–37 (2013). URL http://www.scirp.org/journal/doi.aspx?DOI=10.4236/jmp.2013.47A2005. https://doi.org/10.4236/jmp.2013.47A2005.

[43] Roozemond, P. C., van Erp, T. B. & Peters, G. W. Flow-induced crystallization of isotactic polypropylene: Modeling formation of multiple crystal phases and morphologies. *Polymer* **89**, 69–80 (2016). URL https://linkinghub.elsevier.com/retrieve/pii/S0032386116300325. https://doi.org/10.1016/j.polymer.2016.01.032.

[44] Troisi, E., Arntz, S., Roozemond, P., Tsou, A. & Peters, G. Application of a multi-phase multi-morphology crystallization model to isotactic polypropylenes with different molecular weight distributions. *European Polymer Journal* **97**, 397–408 (2017). URL https://linkinghub.elsevier.com/retrieve/pii/S0014305717309722. https://doi.org/10.1016/j.eurpolymj.2017.09.042.

[45] Acierno, S., Coppola, S., Grizzuti, N. & Maffettone, P. L. Coupling between kinetics and rheological parameters in the flow-induced crystallization of thermoplastic polymers. *Macromolecular Symposia* **185** (1), 233–241 (2002). URL https://onlinelibrary.wiley.com/doi/10.1002/1521-3900(200208)185:1(233::AID-MASY233)3.0.CO;2-W. https://doi.org/10.1002/1521-3900(200208)185:1(233::AID-MASY233)3.0.CO;2-W.

[46] Zheng, R. & Kennedy, P. K. A model for post-flow induced crystallization: General equations and predictions. *Journal of Rheology* **48** (4), 823–842 (2004). URL http://sor.scitation.org/doi/10.1122/1.1763944. https://doi.org/10.1122/1.1763944.

[47] Zheng, R. et al. Modeling of flow-induced crystallization of colored polypropylene in injection molding. *Korea-Australia Rheology Journal* **22** (3), 151–162 (2010). URL http://images.autodesk.com/apac_japan_main/files/kr22-03-0001.pdf.

[48] Laschet, G., Spekowius, M., Spina, R. & Hopmann, C. Multiscale simulation to predict microstructure dependent effective elastic properties of an injection molded polypropylene component. *Mechanics of Materials* **105**, 123–137 (2017). URL https://linkinghub.elsevier.com/retrieve/pii/S0167663616304240. https://doi.org/10.1016/j.mechmat.2016.10.009.

[49] Schrank, T. et al. Injection Molding Simulation of Polyoxymethylene Using Crystallization Kinetics Data and Comparison with the Experimental Process. *Polymer Crystallization* **2022**, 1–15 (2022). URL https://www.hindawi.com/journals/pcrys/2022/2387752/. https://doi.org/10.1155/2022/2387752.

[50] Lamberti, G., Peters, G. W. M. & Tito-manlio, G. Crystallinity and Linear Rheological Properties of Polymers. *International Polymer Processing* **22** (3), 303–310 (2007). URL http://www.hanser-elibrary.com/doi/abs/10.3139/217.2006. https://doi.org/10.3139/217.2006.

[51] Luye, J.-F., Regnier, G., Le Bot, P., Delaunay, D. & Fulchiron, R. PVT measurement methodology for semicrystalline polymers to simulate injection-molding process. *Journal of Applied Polymer Science* **79** (2), 302–311 (2001). URL https://onlinelibrary.wiley.com/doi/10.1002/1097-4628(20010110)79:2(302::AID-APP120)3.0.CO;2-I. https://doi.org/10.1002/1097-4628(20010110)79:2(302::AID-APP120)3.0.CO;2-I.
[52] Fulchiron, R., Koscher, E., Poutot, G., Delaunay, D. & Rgnier, G. Analysis of the pressure effect on the crystallization kinetics of polypropylene: dilatometric measurements and thermal gradient modeling. *Journal of Macromolecular Science, Part B* **40** (3-4), 297-314 (2001). URL https://www.tandfonline.com/doi/abs/10.1081/MB-100106159.

[53] Zhao, P. et al. A novel method for predicting degrees of crystallinity in injection molding during packing stage. *Proceedings of the Institution of Mechanical Engineers, Part B: Journal of Engineering Manufacture* **233** (1), 204-214 (2017). URL http://journals.sagepub.com/doi/10.1177/0954405417718593. https://doi.org/10.1177/0954405417718593.

[54] Cook, P. S., Zuo, Z. & Costa, F. S. *A research framework for cooling rate-dependent PVT models* (Society of Plastics Engineers (SPE), Indianapolis, Indiana, USA, 2016). URL https://www.researchgate.net/profile/Zhihao-Zuo/publication/303703576_A_RESEARCH_FRAMEWORK_FOR_COOLING_RATE-DEPENDENT_PVT_MODELS/links/574ed9bb08aefc38ba120d6a/A-RESEARCH-FRAMEWORK-FOR-COOLING-RATE-DEPENDENT-PVT-FRM-00000000000000000000000000000000.pdf.

[55] Wang, J. et al. Continuous Two-Domain Equations of State for the Description of the Pressure-Specific Volume-Temperature Behavior of Polymers. *Polymers* **12** (2), 409 (2020). URL https://www.mdpi.com/2073-4360/12/2/409. https://doi.org/10.3390/polym12020409.

[56] Hopmann, C., Xiao, C., Kahve, C. E. & Fellerhoff, J. Prediction and validation of the specific volume for inline warpage control in injection molding. *Polymer Testing* **104**, 107393 (2021). URL https://linkinghub.elsevier.com/retrieve/pii/S014294182100338X. https://doi.org/10.1016/j.polymertesting.2021.107393.

[57] Hijazi, S., Stabile, G., Mola, A. & Rozza, G. Data-driven POD-Galerkin reduced order model for turbulent flows. *Journal of Computational Physics* **416** (2020). URL https://linkinghub.elsevier.com/retrieve/pii/S0021999120302874. https://doi.org/10.1016/j.jcp.2020.109513.

[58] Chinesta, F. et al. PGD-Based Computational Vademecum for Efficient Design, Optimization and Control. *Archives of Computational Methods in Engineering* **20** (1), 31-59 (2013). URL http://link.springer.com/10.1007/s11831-013-9080-x. https://doi.org/10.1007/s11831-013-9080-x.

[59] Wang, G. G. & Shan, S. Review of Metamodeling Techniques in Support of Engineering Design Optimization. *Journal of Mechanical Design* **129** (4), 370-380 (2007). URL https://asmedigitalcollection.asme.org/mechanicaldesign/article/129/4/370/466824/Review-of-Metamodeling-Techniques-in-Support-of. https://doi.org/10.1115/1.2429697.

[60] Forrester, A. I. & Keane, A. J. Recent advances in surrogate-based optimization. *Progress in Aerospace Sciences* **45** (1-3), 50-79 (2009). URL https://linkinghub.elsevier.com/retrieve/pii/S003704210000566. https://doi.org/10.1016/j.paerosci.2008.11.001.

[61] Iuliano, E. & Quagliarella, D. Proper Orthogonal Decomposition, surrogate modelling and evolutionary optimization in aerodynamic design. *Computers & Fluids* **84**, 327-350 (2013). URL https://linkinghub.elsevier.com/retrieve/pii/S0045793013002223. https://doi.org/10.1016/j.compfluid.2013.06.007.

[62] Simpson, T., Poplinski, J., Koch, P. N. & Allen, J. Metamodels for Computer-based Engineering Design: Survey and recommendations. *Engineering with Computers* **17** (2), 129-150 (2001). URL http://link.springer.com/10.1007/PL00007198. https://doi.org/10.1007/PL00007198.

[63] Chen, W., Zhou, X.-h., Wang, H.-f. & Wang, W. Multi-objective optimal approach for injection molding based on surrogate
ed crystallization model for injection molding simulation

model and particle swarm optimization
algorithm. *Journal of Shanghai Jiaotong University (Science)* **15** (1), 88-93
(2010). URL http://link.springer.com/10.1007/s12204-010-9517-4. https://doi.org/10.1007/s12204-010-9517-4.

[64] Gao, Y. & Wang, X. Surrogate-based process optimization for reducing warpage in injection molding. *Journal of Materials Processing Technology* **209** (3), 1302-1309 (2009). URL https://linkinghub.elsevier.com/retrieve/pii/S0924013608002586. https://doi.org/10.1016/j.jmatprotec.2008.03.048.

[65] Wang, X., Gu, J., Shen, C. & Wang, X. Warpage optimization with dynamic injection molding technology and sequential optimization method. *The International Journal of Advanced Manufacturing Technology* **78** (1-4), 177-187 (2015). URL http://link.springer.com/10.1007/s00170-014-6621-x. https://doi.org/10.1007/s00170-014-6621-x.

[66] Kang, G.-J., Park, C.-H. & Choi, D.-H. Metamodel-based design optimization of injection molding process variables and gates of an automotive glove box for enhancing its quality. *Journal of Mechanical Science and Technology* **30** (4), 1723-1732 (2016). URL http://link.springer.com/10.1007/s12206-016-0328-x. https://doi.org/10.1007/s12206-016-0328-x.

[67] Li, C., Wang, F.-L., Chang, Y.-Q. & Liu, Y. A modified global optimization method based on surrogate model and its application in packing profile optimization of injection molding process. *The International Journal of Advanced Manufacturing Technology* **48** (5-8), 505-511 (2010). URL http://link.springer.com/10.1007/s00170-009-2302-6. https://doi.org/10.1007/s00170-009-2302-6.

[68] Kitayama, S., Miyakawa, H., Takano, M. & Aiba, S. Multi-objective optimization of injection molding process parameters for short cycle time and warpage reduction using conformal cooling channel. *The

*International Journal of Advanced Manufacturing Technology* **88** (5-8), 1735-1744 (2017). URL http://link.springer.com/10.1007/s00170-016-8904-x. https://doi.org/10.1007/s00170-016-8904-x.

[69] Feng, Q. & Zhou, X. Automated and robust multi-objective optimal design of thin-walled product injection process based on hybrid RBF-MOGA. *The International Journal of Advanced Manufacturing Technology* **101** (9-12), 2217-2231 (2019). URL http://link.springer.com/10.1007/s00170-018-3084-5. https://doi.org/10.1007/s00170-018-3084-5.

[70] Shi, H., Xie, S. & Wang, X. A warpage optimization method for injection molding using artificial neural network with parametric sampling evaluation strategy. *The International Journal of Advanced Manufacturing Technology* **65** (1-4), 343-353 (2013). URL http://link.springer.com/10.1007/s00170-012-4173-5. https://doi.org/10.1007/s00170-012-4173-5.

[71] Cheng, J., Liu, Z. & Tan, J. Multi-objective optimization of injection molding parameters based on soft computing and variable complexity method. *The International Journal of Advanced Manufacturing Technology* **66** (5-8), 907-916 (2013). URL http://link.springer.com/10.1007/s00170-012-4376-9. https://doi.org/10.1007/s00170-012-4376-9.

[72] Xia, W., Luo, B. & Liao, X.-p. An enhanced optimization approach based on Gaussian process surrogate model for process control in injection molding. *The International Journal of Advanced Manufacturing Technology* **56** (9-12), 929-942 (2011). URL http://link.springer.com/10.1007/s00170-011-3227-4. https://doi.org/10.1007/s00170-011-3227-4.

[73] Villarreal-Marroqun, M. G. et al. A comparison of two metamodel-based methodologies for multiple criteria simulation optimization using an injection molding
case study. *Journal of Polymer Engineering* **33** (3), 193–209 (2013). URL https://www.degruyter.com/document/doi/10.1515/polyeng-2013-0022/html. https://doi.org/10.1515/polyeng-2013-0022.

[74] Ivan, R., Sorgato, M., Zanini, F. & Lucchetta, G. Improving Numerical Modeling Accuracy for Fiber Orientation and Mechanical Properties of Injection Molded Glass Fiber Reinforced Thermoplastics. *Materials* **15** (13), 4720 (2022). URL https://www.mdpi.com/1996-1944/15/13/4720. https://doi.org/10.3390/ma15134720.

[75] Gao, H., Zhang, Y., Zhou, X. & Li, D. Intelligent methods for the process parameter determination of plastic injection molding. *Frontiers of Mechanical Engineering* **13** (1), 85–95 (2018). URL http://link.springer.com/10.1007/s11465-018-0491-0. https://doi.org/10.1007/s11465-018-0491-0.

[76] Xu, G. & Yang, Z. Multiobjective optimization of process parameters for plastic injection molding via soft computing and grey correlation analysis. *The International Journal of Advanced Manufacturing Technology* **78** (1-4), 525–536 (2015). URL http://link.springer.com/10.1007/s00170-014-6643-4. https://doi.org/10.1007/s00170-014-6643-4.

[77] Mukras, S. M. S., Omar, H. M. & al Mufadi, F. A. Experimental-Based Multi-objective Optimization of Injection Molding Process Parameters. *Arabian Journal for Science and Engineering* **44** (9), 7653–7665 (2019). URL http://link.springer.com/10.1007/s13369-019-03855-1. https://doi.org/10.1007/s13369-019-03855-1.

[78] Hieber, C. A. Correlations for the quiescent crystallization kinetics of isotactic polypropylene and poly (ethylene terephthalate). *Polymer* **36** (7), 1455–1467 (1995). URL https://www.sciencedirect.com/science/article/pii/003238619595925Q.

[79] Wissbrun, K. F. The equilibrium melting point of polyoxymethylene. *Journal of Polymer Science Part A-2: Polymer*

*Physics* **4** (5), 827–829 (1966). URL https://onlinelibrary.wiley.com/doi/abs/10.1002/pol.1966.160040513.

[80] Li, J., Wang, Y., Wang, X. & Wu, D. Development of Polyoxymethylene/Polylactic Blends for a Potentially Biodegradable Material: Crystallization Kinetics, Lifespan Prediction, and Enzymatic Degradation Behavior. *Polymers* **11** (9), 1516 (2019). URL https://www.mdpi.com/2073-4360/11/9/1516. https://doi.org/10.3390/polym11091516.

[81] Kitano, T., Kataoka, T. & Shirota, T. An empirical equation of the relative viscosity of polymer melts filled with various inorganic fillers. *Rheologica Acta* **20** (2), 207–209 (1981). URL http://link.springer.com/10.1007/BF01513064. https://doi.org/10.1007/BF01513064.

[82] Rudolph, N. M., Agudelo, A. C., Granada, J. C., Park, H. E. & Osswald, T. A. WLF model for the pressure dependence of zero shear viscosity of polycarbonate. *Rheologica Acta* **55** (8), 673–681 (2016). URL http://link.springer.com/10.1007/s00397-016-0945-4. https://doi.org/10.1007/s00397-016-0945-4.

[83] Plummer, C. J. G. & Kausch, H.-H. A DSC investigation of the crystallization kinetics of polyoxymethylene. *Polymer Bulletin* **32** (1), 117–124 (1994). URL https://link.springer.com/article/10.1007/BF00297423.

[84] Metzner, A. B. Rheology of Suspensions in Polymeric Liquids. *Journal of Rheology* **29** (6), 739–775 (1985). URL http://sor.scitation.org/doi/10.1122/1.549808. https://doi.org/10.1122/1.549808.

[85] Plummer, C. J. G. & Kausch, H.-H. DSC non-isothermal crystallization curves in polyoxymethylene. *Colloid and Polymer Science* **273** (3), 227–232 (1995). URL https://link.springer.com/article/10.1007/BF00657828.
ed crystallization model for injection molding simulation

**Statements and Declarations.** The authors declare that no funds, grants, or other support were received during the preparation of this manuscript.

The authors have no relevant financial or non-financial interests to disclose.

All authors contributed to the study conception and design. Material preparation, data collection and analysis were performed by Sandra Saad and Camilo Cruz. The first draft of the manuscript was written by Sandra Saad and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.