Simulation of Three-dimensional Complex Flows in Injection Molding Using Immersed Boundary Method

In this paper, an immersed boundary method (IBM) has been developed to simulate three-dimensional (3D) complex flows in the injection molding process, in which the irregular boundary of mould is treated by a level set function. The melt front (melt-air interface) is captured and treated using the coupled level set and volume of fluid (CLSVOF) method. The finite volume method on the nonstaggered meshes is implemented to solve the governing equations, and the melt filling process is simulated in a rectangular mould with both thickand thin-wall sections. The numerical result shows good agreement with the available data. Finally, the melt filling processes in three different moulds, i.e., ring-shaped channels with concentric, eccentric and thin cylindrical insets, are investigated respectively. It is found that the so-called “race-tracking effect” and air traps phenomena could be observed successfully, which demonstrates the applicability of IBM to complex flows in injection molding.


Introduction
The Flow and heat transfer are two very complicated physical processes in plastic injection molding, especially in complex-shaped moulds. The simulation of the filling processes is of most importance to understand and design the moulds, so many researches have be done to simulate the complex flows in injection molding using different methods, among which the finite element method (FEM) is usually used for it can treat complex boundaries more easily [Polynkin et al. (2005); Ilinca and Hétu (2002)]. However, it takes too much computational time, especially in 3D space. The finite volume method (FVM) is an alternative for its lower computational cost [Chang (2001); Chau (2006); Araujo (2008) ;]. To deal with complex boundaries, in FVM, unstructured grids or body-fitting grids are usually employed to conform to complex moulds [Chang (2001); Araujo (2008)]. But the mesh generation is time consuming. In recent years, the immersed boundary method (IBM) has achieved great progress, one of its big advantages is that the governing equations can be solved on Cartesian grid easily with a body force prescribed on boundaries, and it does not require complex morphing that can deteriorate the mesh or re-meshing of the computational domain. Hereafter, the developments and applications of IBM will be briefly reviewed.
The IBM, one of the most useful computational methods in studying fluid structure interaction, is applied in various aspects of fluid dynamics. The IBM was originally proposed [Peskin (1972); Peskin (1977)] to simulate the blood flow in a human heart model. The prominent feature of this method was that the entire simulation was carried out on a Cartesian grid system, which did not conform to the geometry of the heart, and a novel procedure was imitated for imposing the effect of the immersed boundary on the flow. Since Peskin introduced this method, a number of modifications and refinements have been proposed [Roma et al. (1999); Kim et al. (2001);Ji et al. (2012); Tseng and Ferziger (2003); Seo and Mittal (2011);Mariano et al. (2010); Silva et al. (2009);Yang et al. (2009)]. The Cartesian mesh (or body-fitting mesh) was usually used in the IBM by applying a body force into the virtual boundary, the boundary can reach to a no-slip condition. Choi et al. (2007) presented an IBM for time-dependent, three-dimensional, incompressible flows. The immersed boundary objects were rendered as level sets in the computational domain, and the flow induced by realistic human walking motion was simulated as an example of a problem involving multiple moving immersed objects.
In addition, there is another class of methods, usually referred to as "Cartesian grid methods", which was originally developed for simulating inviscid flows with complex embedded solid boundaries on Cartesian grids [Aftosmis et al. (1998);De Zeeuw and Powell (1991)]. These methods had been extended to simulate unsteady viscous flows [Udaykumar et al. (1996); Ye et al. (1999)] and thus have capabilities similar to those of IBM.
Another branch of Cartesian grid method, cut-cell method, has succeeded in simulating two-phase flows with embedded solid boundaries, which has been successfully applied to the Euler equations in two and three dimensions, to flows involving both moving bodies and moving material interfaces [Ilinca, Hartmann et al.(2011);Bouchon et al. (2012); Meinke et al. (2012)]. The conceptually simple approach "cuts" solid bodies out of a background Cartesian mesh. In this review, the meaning of IBM encompasses all such methods that simulate viscous flows with immersed (or embedded) boundaries on Cartesian grid system. This paper mainly focuses on the numerical simulation of melt injection in complex mould. The two-phase flow model proposed in ] is adopted to simulate the melt filling, where the Cross-viscosity model is employed to describe the viscous behavior of polymer melt. The governing equations of two-phase flows are solved using the finite volume and immersed boundary methods (FV-IBM). To my knowledge this is the first example in which the immersed boundary technique is presented and used to treat the irregular boundary of the mould. Meanwhile, two level set functions are employed, one for treating the complex moulds, another for tacking melt front with the aid of volume-of-fluid (VOF) method, which is the so-called CLSVOF (coupled level set and volume of fluid) method [Sussman and Puckett (2000); Son (2003); ]. The content of this paper is listed as follows. First, the mathematical model is proposed in Section 2. Section 3 presents the numerical implementation of the FV-IBM. In Section 4 the proposed IBM is tested and verified by comparing the numerical results with available data. As a case study, melt filling processes are simulated and analyzed in detail in three different moulds, i.e., ring-shaped channels with concentric, eccentric and thin cylindrical insets. Some conclusions and future research directions conclude this paper.

Governing Equations
In the melt filling process, since the air velocity is low, both the air and liquid phases can be regarded as incompressible flows [Yang et al. (2010);]. The continuity, momentum and energy equations of the incompressible fluids can be written as the unified equations in dimensionless form f is the body force(virtual force), which can be prescribed on a regular mesh in IBM. µ is the solid volume fraction in a computational cell, i.e. µ = 0 for the fluid, µ = 1 for the solid, and 0 < µ < 1 for the solid/fluid interface. CMES, vol.94, no.3, pp.207-224, 2013 In three-dimensional space, the governing equations could be written as

CLSVOF method and shaped level set function
In the CLSVOF method, the level set function is adopted to capture the melt front with the aid of VOF function, while the level set and VOF functions are advected using the following equation, respectively where Φ is level set function φ or VOF function F. Please see [Son and Hur (2002); Son (2003)] for important details about CLSVOF method.
The melt front is represented by the level set function φ (x,t) . And the sign of φ (x,t) is that φ > 0 in the melt phase and φ < 0 in the air phase, respectively. The level set function is employed to treat the discontinuities of density and viscosity near the interface. The density and viscosity are assumed to be constant in each phase, while the fluid properties across the interface are smoothed over a transition band using the level set function where ρ, η, C and κ are density, viscosity, thermal capacity and conductivity, respectively. The subscript a and m represent air and melt, respectively. The where ε is a parameter related to the interface thickness, herein ε = 1.5∆x. ∆x is the grid width along the x direction.
The Cross-viscosity model is adopted to describe the rheological property of the polymer melt in this paper [Chau (2008); .
where τ * is the model constant that shows the shear stress rate, n is the model constant which symbolizes the pseudoplastic behavior slope of the melt as (1 − n), η 0 is the melt viscosity under zero-shear-rate conditions. η 0 can be expressed as Arrhenius formula Eq. (14) is 5-parameter model, where B, T b and β are material parameters.
As far as the author is aware, the method based on Cartesian mesh, such as IBM or cut-cell method, has not been applied to numerical simulation of polymer melt injection molding yet. The melt front and mould shape are described using two level sets, one for the liquid/air phase, another for the solid/fluid phase, and the irregular boundary Γ of the mould can be implicitly defined as the zero-contour of a level set (signed distance) function φ f , which is stored at each cell center of the Cartesian background mesh.
An irregular solid domain Ω s is embedded in the computational domain Ω, such that Ω f = Ω\Ω s represents the fluid domain where the governing equations are to be discretized. To describe the irregular boundary of the mould, a level set function is employed such that φ f is negative in the fluid region Ω f and positive in the solid region Ω s . Meanwhile, the boundary Γ corresponds to the zero level set of this function, i.e., where dist represents the distance between x and the nearest point on the irregular boundary of the mould. In order to build the level set function φ f that represents the irregular boundary Γ of the mould, the Constructive Solid Geometry (CSG) method is used, which can construct complex domains out of basic geometries such as circles, hyperplanes and spheres [Cheny and Botella (2009); Cheny and Botella (2010)]. For example, let Ω 1 and Ω 2 be the cylinder boundaries Γ 1 and Γ 2 , respectively. And their level set functions are, where R 1 and R 2 are radii of a concentric cylindrical, and (x 0 ,y 0 ) is the center. Thus, the level set function of the concentric cylindrical is simply φ f = − min (φ 2 , −φ 1 ) (see Fig. 1). Other complex geometries could also be created by Boolean CSG operations [Osher and Fedkiw (2002)]. . The VOF function F is employed to conserve the melt mass which would be lost in the level set method. To solve the VOF function F, the flux-splitting algorithm is used. Please refer to [Son (2003)] for detailed description.

Body Force Solvers
In IBM, the body force f , including f x , f y and f z , should be calculated. For simplicity, here take the calculation of f x as an example. In each computational cell the fluid velocityu , which is firstly calculated without body force f x , is recomputed using(1 − µ) u , and then the body force f x could be obtained using Eqs. (18-19) [ Fadlun et al. (2000)] where u * is the fluid velocity recomputed using (1 − µ) u . So the following formula can be obtained, where µ is the solid volume fraction which can be calculated using following formulas from the shaped level set function φ f [ Van der Pijl et al. (2008)].
Similarly, f y and f z can be obtained using the above method, so I will not repeat them here.

The numerical simulation of complex filling process
In this section, the injection molding processes of Cross fluid are simulated using IBM and CLSVOF method, where the shaped level set function is employed to treat complex mould boundaries.

Filling process in a rectangular mould with both thick-and thin-wall sections
As a test example, the filling process is simulated in a rectangular cavity with both thick-and thin-wall sections. Fig. 2 shows the geometry and size of the rectangular mould. The selected material is Taiwan PP6733, whose parametric constants corresponding to n, τ * , D 1 , D 2 , D 3 , A 1 and A 2 of the viscosity model are 0.219, 56, 1.3×10 14 , 263.15, 0, 28.44, 51.6. The temperatures of melt and mould wall are 240˚C and 50˚C, respectively. The injection velocity is 50cm 3 /s, and the space and time intervals are 0.075 and 0.012, respectively. Fig. 3 shows the melt interfaces of filling process at different time. If the mould has both thick-and thin-wall sections, the thick-wall section has lower flow resistance than the thin-wall section. Therefore the melt flows fast in the thick-wall section, which is so-called "race-tracking effect". As shown in Fig. 3(a), the phenomenon of race-tracking effect is clearly observed when using the 3D model mentioned above. From Fig.3 the melt interfaces at different time are quite the same as those obtained by HseAE3D ( Fig. 3(b)) and Moldflow softwares (Fig. 3(c)) in [Zhou et al. (2008)], which indicates that the IBM proposed in this paper, with CLSVOF method to capture the moving interface, could simulate the complex flows accurately in 3D mould.
Above example is relatively simple. Hereafter the filling processes in three ringshaped channels will be simulated. And the typical phenomena, such as racetracking effect and air traps, will be observed clearly. Unless otherwise stated, the polymer material and other parameters for simulations remain unchanged.

Filling process in ring-shaped channel with concentric cylindrical inset
The geometry and size of the mould is shown in Fig. 4, and its thickness is 1. Fig. 5 demonstrates the melt interfaces in ring-shaped channel with concentric cylindrical  [Zhou et al. (2008)] (c)Moldflow, MPI_3D in [Zhou et al. (2008)] inset at different time. After the melt turns around the circular insert two melt branches encounter and the cavitation forms among melt and insert (Fig. 5(d)), and then the seam line or weldline begins to form (Fig. 5(d)-(e)), which is undesirable in the injection molding.

Filling process in ring-shaped channels with eccentric cylindrical inset
Since the core location is an important parameter in the filling process, in order to further investigate the complex flows in injection molding, the influence of core location on the flow is considered in this subsection.   6 gives the filling process in ring-shaped channel with eccentric cylindrical inset whose core centre is located at O (3.5, 1.2) and the radius is 0.55, while the other parameters are the same as those in Fig. 4. From Fig. 6, the difference between two melt branches is unobservable in a short time after the polymer melt touches the core (Fig. 6(a)). As time goes by, the velocity difference between the two branches is gradually enlarged since the flow resistance is lower in the right channel. Then the right branch goes faster than the left one due to the race-tracking effect (Figs. 6(b)-(c)). After the two branches contact, the seam line is formed. And air traps also occur among inset and melt. Compared Fig. 6 with Fig. 5, it can be found that big eccentricity ratio results in great difference of the two melt branches, the core location plays an important role in the molding filling process.

Filling process in ring-shaped channels with thin cylindrical inset
From above subsections, it is known that the race-tracking effect can be also observed in the filling processes in ring-shaped cannels with eccentric cylindrical inset, which could also be simulated successfully in 2D space. However, the filling flows in more complex shaped mould would only be simulated successfully in full 3D space. Hereafter, as a typical example, the filling process will be simulated in a ring-shaped channel with thin cylindrical inset.    8 shows the filling process in ring-shaped channel with thin cylindrical inset at different time. When polymer melt flows around a cylindrical inset which is thinner than other segments of the mould, in a short time, the difference between two melt branches is unobservable. At t=0.432, the two branches emerge (Fig. 8(d)). Then the difference enlarges gradually due to race-tracking effect. After flowing around the inset, the two branches meet each other (Fig. 8(f)) and the seam line begins to form which gradually disappears over time ( Fig. 8(g)). Fig. 8

Conclusions
In this work, the immersed boundary and CLSVOF methods are proposed to simulate the melt filling processes in complex moulds, which is tested by the filling process in a rectangular mould with both thick-and thin-wall sections. Then the filling processes are simulated in ring-shaped channels with concentric, eccentric and thin cylindrical insets, respectively. And the conclusions can be drawn as follows.
1. With a shaped level set function to treat the complex boundaries of the mould, the immersed boundary and CLSVOF methods are demonstrated the capability in handling complex flow problems, which could be applied to simulate melt filling process in arbitrary shaped moulds.
2. The race-tracking effect is a phenomenon that polymer melt front along the thicker edge moves faster than that in the thinner area, which could cause weld line to form or air to be trapped in the polymer melt. Due to race-tracking effect, the melt flow is unbalanced in the thick and thin segments in the mould, which will influence weightily the quality and performance of the final plastic products. So in injection molding process the injection port should be added in the thin sections or the thickness of each part should be consistent in the mould.