A Novel Multiscale Methodology for Simulating Droplet Morphology Evolution during Injection Molding of Polymer Blends

The morphology of polymer blends plays a critical role in determining the properties of the blends and performance of resulting injection-molded parts. However, it is currently impossible to predict the morphology evolution during injection molding and the final micro-structure of the molded parts, as the existing models for the morphology evolution of polymer blends are still limited to a few simple flow fields. To fill this gap, this paper proposed a novel model for droplet morphology evolution during the mold filling process of polymer blends by coupling the models on macro- and meso-scales. The proposed model was verified by the injection molding experiment of PP/POE blends. The predicted curve of mold cavity pressure during filling process agreed precisely with the data of the corresponding pressure sensors. On the other hand, the model successfully tracked the moving trajectory and simulated morphology evolution of the droplets during the mold-filling process. After mold-filling ended, the simulation results of the final morphology of the droplets were consistent with the observations of the scanning electron microscope (SEM) experiment. Moreover, this study revealed the underlying mechanism of the droplet morphology evolution through the force analysis on the droplet. It is validated that the present model is a qualified tool for simulating the morphology evolution of polymer blends during injection molding and predicting the final microstructure of the products.


Introduction
Polymer blends are ubiquitously used in machinery, electronics, transportation and our everyday life due to their enhanced properties such as light weight, high strength, and chemical resistance [1][2][3]. For the large unfavorable enthalpy of mixing, most physical blends of different polymers prove to be immiscible. The resultant small-scale arrangement of the phases constitutes the microstructure inside the blends. Microstructure is determined by comments properties, flow history and fractions. If the fractions of the components are sufficiently different, the sea-island-like microstructure will emerge in the blends where numerous droplets are immersed in the matrix [4]. Injection molding is a widely-used processing method in manufacturing polymer blends products. With the blend melt filling the mold cavity on macroscale, at the same time the droplets morphologies undergo drastic change on the mesoscale. The microstructure of the polymer blend has a significant effect on the mechanical, thermal, optical and physical properties of polymer blends and their injection-molded parts [5]. The latest researches showed that the orientation, polymer blends on macroscale. Coupling the models on different scales of polymer melt flow and droplet morphology evolution has been one of the long-standing challenges in multiscale modeling of polymer dynamics for years [28].
It is feasible to analyze the morphology evolution process of polymer blends during injection molding using a macro-and meso-scale coupled model. The finite volume method (FVM) was used in the proposed model to simulate the polymer melt flow on the macroscopic scale [29,30]. The lattice Boltzmann method (LBM) was employed to describe the complex morphology evolution of droplets on the mesoscopic scale. In recent years, the LBM has gained much success in studying multiphase flows, i.e., emulsion flow, bubble flow and binary flow [31,32]. Its fundamental idea is to construct simplified kinetic models that incorporate the essential physics of microscopic processes [33]. Macroscopic hydrodynamic behaviors, such as interface dynamics, naturally emerge as a result of this kinetics. This paper proposed a novel coupling framework, consisting of the droplet trajectory tracking method and the phenomenological constitutive model of polymer blends. Given the macro-and meso-scales coupled, a multiscale model for simulating the droplet morphology evolution during the mold filling process has been established.
This paper is organized as follows. Section 2 detailed the coupled macro-and mesoscale model including the FVM-based model for mold filling flow of polymer blends melt, the LBM-based model for droplet morphology evolution and a bi-directional coupling framework between them. To validate the proposed model, the simulation results were compared with two injection molding experiments in Section 3. One is the mold cavity pressure monitoring experiment and the other is the SEM imaging of the droplet morphologies of the PP/POE part. The conclusions drawn from this paper and the outlook for future research were briefed in Section 4.

Multiscale Model
The proposed multiscale approach was the integration of macroscopic LBM and mesoscopic FVM. The macro-and meso-scale models were coupled in a novel way that the macroscopic flow simulation defined the exterior conditions surrounding the droplets while the mesoscopic droplets morphologies contributed to the rheology of the polymer blends. This coupling approach was achieved through tracking the droplet trajectory following the polymer melt flow and the phenomenological constitutive model for molten polymer blends.

Macroscale Model for Mold Filling Flow of Polymer Melt
On the macroscale, the FVM with SIMPLE algorithm was used to simulate the mold filling flow of polymer blends melt. The governing equations of velocity, pressure and temperature inside the mold cavity could be written as: Continuity equation: Momentum equation: Energy equation The variables and their meanings of the above equations were listed in Table 1 below: As the mold filling process was a kind of free surface flow, the volume of fraction (VOF) method was employed to re-construct the melt-air interface advancement, and the transportation equation for the filled fraction went as follows: where χ was the filled fraction of control volume: χ(x, t) = 1, filled with polymer melt 0, filled with air (5) The conservative Equations (1)-(3) together with the transportation Equation (4) constituted the governing equations of the mold filling flow of the polymer blends melt.
The boundary conditions for the flow continuity and momentum equations must be prescribed to complete the mathematical description of the whole flow domain. The boundary conditions on the inlet and mold wall went as: At the flow front, interfaces between the melt and air, the traction free boundary condition and the atmospheric pressure condition were imposed.
For the energy equation, the initial conditions for the flow domain and the boundary conditions for the injection gate and mold wall were defined as follows: Besides, the initial and boundary conditions for χ of Equation (4) were: The governing Equations (1)-(4) can be written in the form of the following generalized transport equations The physical meanings of φ, Λ, Γ and Q φ depend on the specific equation represented by Equation (13), as detailed in Table 2. Table 2. Physical meanings of the symbols in the generalized transport equation.
By integrating Equation (13) and applying the Gaussian divergence theorem, the discretized equation was obtained after some algebraic operations: Finally, the discrete equations of velocity, pressure and temperature were solved using the SIMPLE algorithm while the CICSAM method proposed by Ubbink and Issa et al. was utilized to maintain the sharpness of the melt-air interface and the boundedness of χ [34].

Mesoscale Droplet Morphology Evolution
On the mesoscale, the Shan-Chen pseudo-potential scheme of LBM [35,36] was employed to model the morphology evolution of the droplet immersed in the matrix.
Under the LBM framework, the continuous velocity space was discretized into a finite set of N vectors, and the density distribution function f i (x, t) denoted the probability of finding a particle with the discrete velocity e i at position x and time t. In this paper, the D3Q19 lattice model was used, as shown in Figure 1.

transient term
The physical meanings of φ , Λ , Γ and Q φ depend on the specific equation represented by Equation (13), as detailed in Table 2.
By integrating Equation (13) and applying the Gaussian divergence theorem, the discretized equation was obtained after some algebraic operations: Finally, the discrete equations of velocity, pressure and temperature were solved using the SIMPLE algorithm while the CICSAM method proposed by Ubbink and Issa et al. was utilized to maintain the sharpness of the melt-air interface and the boundedness of χ [34].

Mesoscale Droplet Morphology Evolution
On the mesoscale, the Shan-Chen pseudo-potential scheme of LBM [35,36] was employed to model the morphology evolution of the droplet immersed in the matrix.
Under the LBM framework, the continuous velocity space was discretized into a finite set of N vectors, and the density distribution function ( ) , i f t x denoted the probability of finding a particle with the discrete velocity i e at position x and time t . In this paper, the D3Q19 lattice model was used, as shown in Figure 1.  For the components of droplets and the matrix, they were represented by f d i (x, t) and f m i (x, t), respectively, and their governing equations were expressed as: where f σ,eq i (x, t) was the equilibrium distribution function corresponding to the component σ with velocity e i and was determined by: with the weight coefficients: Polymers 2021, 13, 133 6 of 18 τ σ was the dimensionless relaxation time and related to the kinematic viscosity of the component by ν σ = (τ σ − 0.5)/3. ρ σ and u σ were the macroscopic density and velocity of component σ and calculated from the zeroth and first moments of f σ i (x, t): It should be noted that the expression of the equilibrium velocity u σ,eq took into account the interaction from other components acted on component σ: In the Shan-Chen pseudo-potential scheme, the total force on component σ results from the interaction of component σ at neighboring lattices, as illustrated in Figure 2.
where G denotes the strength of the interaction and would determine the resultant interfacial tension as well as the interface thickness.  (19) It should be noted that the expression of the equilibrium velocity ,eq σ u took into account the interaction from other components acted on component σ : (20) In the Shan-Chen pseudo-potential scheme, the total force on component σ results from the interaction of component σ at neighboring lattices, as illustrated in Figure 2.
x e e (21) where G denotes the strength of the interaction and would determine the resultant interfacial tension as well as the interface thickness.   For the size of the droplets was far smaller than that of a control volume, when a droplet, marked as P, passed through a control volume, the droplet can be treated as a particle and its trajectory was tracked through first-order integral of its velocity. The droplet acceleration was obtained by the chain rule of derivative, i.e., in the x direction: Integrating the two sides of Equation (22) between time t n and t n+1 (t n+1 = t n + ∆t) gave: Using a similar derivation process, the coordinates of the droplet in the y-direction and the z-direction can also be updated from t n to t n+1 : The tracking operation was repeated in all control volumes where droplet P traversed one by one. During the mold filling flow, the droplets were largely subjected to the shear force of the flow field that depend on the location of the droplet. Given the droplet's trajec-tory, the history of shear force on the droplet during mold filling was decided accordingly by interpolation, and then the exterior conditions for the morphology evolution of the droplets were well prepared.

Model Set-Up for Droplet Morphology
In multiscale fluid dynamics and modeling, the complicated flow field on the macro scale always appeared as simple flow forms on the local mesoscale [37]. Macroscopically, the mold filling flow of polymer melt was a typical laminar flow [38][39][40][41], while mesoscopically the droplets immersed in the mold filling flow were mainly subjected to shear flow forces [42][43][44]. Based on this analysis, a local model of droplet morphology evolution on the mesoscale was set up accordingly.
As shown in Figure 3, a virtual shear flow field of time-varying shear rate . γ(t) was produced by two oppositely-moving plates. An initially sphere droplet was placed in the center and the droplet morphology was about to change under the shear force.
The tracking operation was repeated in all control volumes where droplet P traversed one by one. During the mold filling flow, the droplets were largely subjected to the shear force of the flow field that depend on the location of the droplet. Given the droplet's trajectory, the history of shear force on the droplet during mold filling was decided accordingly by interpolation, and then the exterior conditions for the morphology evolution of the droplets were well prepared.

Model Set-Up for Droplet Morphology
In multiscale fluid dynamics and modeling, the complicated flow field on the macro scale always appeared as simple flow forms on the local mesoscale [37]. Macroscopically, the mold filling flow of polymer melt was a typical laminar flow [38][39][40][41], while mesoscopically the droplets immersed in the mold filling flow were mainly subjected to shear flow forces [42][43][44]. Based on this analysis, a local model of droplet morphology evolution on the mesoscale was set up accordingly.
As shown in Figure 3, a virtual shear flow field of time-varying shear rate ( ) t γ was produced by two oppositely-moving plates. An initially sphere droplet was placed in the center and the droplet morphology was about to change under the shear force.  The shear rate . γ was not constant but calculated according to the location of the droplet in the flow filed. For instance, given the locations of the droplet at time t n and t n+1 by tracking the droplet trajectory, then the values of shear rate and viscosity at time t n and t n+1 was obtained from the macroscopic simulation. During the timestep, the shear rate on the droplet was assumed linearly varying from .
where (x n , y n , z n ) and (x n+1 , y n+1 , z n+1 ) were the coordinates of the droplet at time t n and t n+1 separately. The morphology evolution of the droplet during one timestep from t n to t n+1 was simulated using the Shan-Chen pseudo-potential LBM, namely solving Equation (15).

Polymer Blends Constitutive Equation
Compared with single-component polymers, the droplets morphologies directly contributed to the rheology of polymer blends and affected the macroscopic flow of the polymer blends. Yu quantified the effect of the interface structure on the local stress of the blends based on the ellipsoidal description of droplets [45]. However, the number density of droplets in polymer blends was extremely large, the morphologies of droplets were irregular and complex so that developing a theoretical constitutive model for polymer blends was nearly impossible. From another perspective, it is valuable to propose a phenomenological constitutive model without explicitly including the interface term. In this paper, the widely-used Cross-WLF equation, Equations (28) and (29), was utilized to calculate the viscosities of the polymer blends melt for the macroscopic simulation of the mold filling flow.
where n was the power law index in the high shear rate regime, . τ was the critical stress level at the transition to shear thinning, A 1 , A 2 , D 1 , D 2 and D 3 were the material coefficients. All these coefficients in the Cross-WLF can be determined by fitting the rheological experiment data of the polymer blends.

Algorithm Summary
To summarize, the proposed multiscale model can be implemented with the following steps for the droplet during the mold filling process: (1) Given the initial macroscopic flow field of u 0 , p 0 and T 0 , the droplet morphology and coordinate r 0 at time t 0 ; (2) Calculate the shear rate on the droplet . γ 0 according to its coordinate r 0 ; (3) Solve Equations (1)-(4) to update the velocity, pressure and temperature to u 1 , p 1 and T 1 to time t 1 = t 0 + ∆t; (4) Update the coordinate of the droplet to r 1 at time t 1 according Equations (23)-(25); (5) Calculate the shear rate on the droplet . γ 1 according to its coordinate r 1 ; (6) Solve Equations (15) and (27) to simulate the droplet morphology evolution during the timestep between t 0 and t 1 ; (7) Go back to step (2) until the mold filling ends.
To make it clearer, the bidirectional coupled framework for the droplet morphology evolution during the mold filling of polymer blends was outlined in Figure 4 below. (6) Solve Equations (15) and (27) to simulate the droplet morphology evolution during the timestep between 0 t and 1 t ; (7) Go back to step (2) until the mold filling ends.
To make it clearer, the bidirectional coupled framework for the droplet morphology evolution during the mold filling of polymer blends was outlined in Figure 4 below.

Experimental Validation
As the proposed model spanned over macro-and meso-scales, two sets of injection molding experiments of polymer blends were conducted. Both sets of experiments were carried out in the HTL-90-F5B injection molding machine (produced by Ningbo Haitai Plastic Machinery Co., Ltd., Ningbo, China). The processing parameters were defined mainly according to the experience of the injection molding in consideration of the geometry, size and material of the molded product, assisted by the parameters provided by Moldflow Synergy as well. Meanwhile, corresponding simulations were carried out and the comparison between the experiments and simulations was made in order to validate the model.

Experimental Validation
As the proposed model spanned over macro-and meso-scales, two sets of injection molding experiments of polymer blends were conducted. Both sets of experiments were carried out in the HTL-90-F5B injection molding machine (produced by Ningbo Haitai Plastic Machinery Co., Ltd., Ningbo, China). The processing parameters were defined mainly according to the experience of the injection molding in consideration of the geometry, size and material of the molded product, assisted by the parameters provided by Moldflow Synergy as well. Meanwhile, corresponding simulations were carried out and Polymers 2021, 13, 133 9 of 18 the comparison between the experiments and simulations was made in order to validate the model.

Cavity Pressure Variation during Mold Filling
Cavity pressure prediction was the benchmark test for the macroscopic simulation of mold-filling flow. As shown in Figure 5, the molded part was a thin slat with the length, width and thickness of 220, 40 and 2.5 mm, respectively. Two pressure sensors were fixed on the surface of the mold cavity, one near the injection gate, P1, and the other one near the end, P2. The processing parameters of this experiment was prescribed in Table 3. The polymers used in the experiment were the blends of polypropylene (PP, trademark T36f, supplied by Sinopec Corp., Wuhan, China) and polyolefin elastomers (POE, trademark 8150, supplied by DowDuPont Inc., Wilmington, NC, USA).The major physical properties of PP and POE were listed in Table 4. The granules of 25% POE and 75% PP were blended using a Haake (Vreden, Germany) twin-screw extruder machine. There are two reasons for the fixed fraction of POE. Firstly, this research was in the context of the manufacturing product of polymer blends. Trough mechanical testing of the products of various fraction, it is found that the molded parts of PP/POE blends with this fraction demonstrates the optimal mechanical performance and formability. Secondly, the fraction of components decides the morphology of the disperse phase and this paper was mainly focused on the sea-island microstructure. The processing parameters of this experiment was prescribed in Table 3. Table 3. Processing parameters in of the cavity pressure experiment. The polymers used in the experiment were the blends of polypropylene (PP, trademark T36f, supplied by Sinopec Corp., Wuhan, China) and polyolefin elastomers (POE, trademark 8150, supplied by DowDuPont Inc., Wilmington, NC, USA).The major physical properties of PP and POE were listed in Table 4. The granules of 25% POE and 75% PP were blended using a Haake (Vreden, Germany) twin-screw extruder machine. There are two reasons for the fixed fraction of POE. Firstly, this research was in the context of the manufacturing product of polymer blends.

Processing Conditions Parameters
Trough mechanical testing of the products of various fraction, it is found that the molded parts of PP/POE blends with this fraction demonstrates the optimal mechanical performance and formability. Secondly, the fraction of components decides the morphology of the disperse phase and this paper was mainly focused on the sea-island microstructure.
Measured by the capillary rheometer, the shear viscosities versus shear rate of PP, POE and PP/POE blends at 210 • C were plotted in Figure 6. Fitting the rheological data to the Cross-WLF viscosity model, seven model parameters of PP/POE blends were calculated and listed in Table 5.  Figure 7 compared the cavity pressure variation over time at P1 and P2 that were obtained numerically and experimentally. It can be seen from the comparison that the pressure transmission and rise at P1 was more complete and rapid than that of P2 because of different distances to the injection gate. The agreement between the experimental data and simulation results indicated that the macroscopic part of the present model was definitely validated.  Fitting the rheological data to the Cross-WLF viscosity model, seven model parameters of PP/POE blends were calculated and listed in Table 5.  Figure 7 compared the cavity pressure variation over time at P1 and P2 that were obtained numerically and experimentally. It can be seen from the comparison that the pressure transmission and rise at P1 was more complete and rapid than that of P2 because of different distances to the injection gate. The agreement between the experimental data and simulation results indicated that the macroscopic part of the present model was definitely validated. Fitting the rheological data to the Cross-WLF viscosity model, seven model parameters of PP/POE blends were calculated and listed in Table 5.  Figure 7 compared the cavity pressure variation over time at P1 and P2 that were obtained numerically and experimentally. It can be seen from the comparison that the pressure transmission and rise at P1 was more complete and rapid than that of P2 because of different distances to the injection gate. The agreement between the experimental data and simulation results indicated that the macroscopic part of the present model was definitely validated.   Figure 8 showed the part for the observation and simulation of the droplet morphology evolution. For convenience of making brittle sections at the designated locations on the part, a smaller-sized (80 × 10 × 4 mm 3 ) product was used, instead of the slat in the previous experiment. After molded, the parts would be frozen in liquid nitrogen for 24 h and then made brittle sections, and the solidified droplets morphologies were observed with the SEM.    Besides validating the droplet morphology simulation, the effect of injection rate on droplet morphology was investigated numerically. For this purpose, three cases (I, II and III) were implemented and injection rate was the single variable among them, as was marked in Table 6 below.  Figure 10 showed the simulated trajectories and morphologies of the droplets of case II. Two droplets, marked A and B, were involved in the simulation. To better present the simulated process of droplets morphologies evolution, the droplets morphologies at location (2) were supplemented in addition to (1) and (3).

Experiments and Simulation Set-Up
On droplets trajectories, it was obvious that the two droplets diverged almost shortly after they entered the mold cavity through the injection gate and subsequently flowed     Besides validating the droplet morphology simulation, the effect of injection rate on droplet morphology was investigated numerically. For this purpose, three cases (I, II and III) were implemented and injection rate was the single variable among them, as was marked in Table 6 below.  Figure 10 showed the simulated trajectories and morphologies of the droplets of case II. Two droplets, marked A and B, were involved in the simulation. To better present the simulated process of droplets morphologies evolution, the droplets morphologies at location (2) were supplemented in addition to (1) and (3).
On droplets trajectories, it was obvious that the two droplets diverged almost shortly after they entered the mold cavity through the injection gate and subsequently flowed parallel to the injection direction. Droplet A drifted mainly in the shear layer, roughly 1.0 Besides validating the droplet morphology simulation, the effect of injection rate on droplet morphology was investigated numerically. For this purpose, three cases (I, II and III) were implemented and injection rate was the single variable among them, as was marked in Table 6 below.  Figure 10 showed the simulated trajectories and morphologies of the droplets of case II. Two droplets, marked A and B, were involved in the simulation. To better present the simulated process of droplets morphologies evolution, the droplets morphologies at location (2) were supplemented in addition to (1) and (3). from the part surface. This trajectory pattern coincided with the laminar streamline of the mold-filling of polymer melt.
On droplets morphologies, the deformation, breakup and retraction of the droplets were figured out on their trajectories. Droplet A underwent drastic deformation and even broke up before the polymer flow ceased. In sharp contrast, droplet B was only slightly deformed throughout its trajectory and nearly retracted to a sphere finally.  (1) and (3) inside the injection molded part of case II respectively. Both sets of SEM images were focused on the shear layer, precisely 1.0 mm from the part surface, so they were compared with the simulation results of droplet A. From Figures 11b and 12b, the cross-sectional profile of the droplets on the perpendicular section was approximately ellipses, suggesting that the droplets orientated at some an angle to the injection direction. By comparing Figures 11a  and 12a, it could be seen that the droplets were substantially deformed at location (1) and there were abundant evidence of droplet breakup and retraction at location (2). Both phenomenon on the SEM images agreed very well with the above simulations.   On droplets trajectories, it was obvious that the two droplets diverged almost shortly after they entered the mold cavity through the injection gate and subsequently flowed parallel to the injection direction. Droplet A drifted mainly in the shear layer, roughly 1.0 mm from the part surface, while droplet B mainly in the core layer, approximately 2.0 mm from the part surface. This trajectory pattern coincided with the laminar streamline of the mold-filling of polymer melt.
On droplets morphologies, the deformation, breakup and retraction of the droplets were figured out on their trajectories. Droplet A underwent drastic deformation and even broke up before the polymer flow ceased. In sharp contrast, droplet B was only slightly deformed throughout its trajectory and nearly retracted to a sphere finally. Figures 11 and 12 demonstrated the SEM images taken at locations (1) and (3) inside the injection molded part of case II respectively. Both sets of SEM images were focused on the shear layer, precisely 1.0 mm from the part surface, so they were compared with the simulation results of droplet A. From Figures 11b and 12b, the cross-sectional profile of the droplets on the perpendicular section was approximately ellipses, suggesting that the droplets orientated at some an angle to the injection direction. By comparing Figures 11a and 12a, it could be seen that the droplets were substantially deformed at location (1) and there were abundant evidence of droplet breakup and retraction at location (2). Both phenomenon on the SEM images agreed very well with the above simulations. from the part surface. This trajectory pattern coincided with the laminar streamline of the mold-filling of polymer melt. On droplets morphologies, the deformation, breakup and retraction of the droplets were figured out on their trajectories. Droplet A underwent drastic deformation and even broke up before the polymer flow ceased. In sharp contrast, droplet B was only slightly deformed throughout its trajectory and nearly retracted to a sphere finally.  (1) and (3) inside the injection molded part of case II respectively. Both sets of SEM images were focused on the shear layer, precisely 1.0 mm from the part surface, so they were compared with the simulation results of droplet A. From Figures 11b and 12b, the cross-sectional profile of the droplets on the perpendicular section was approximately ellipses, suggesting that the droplets orientated at some an angle to the injection direction. By comparing Figures 11a  and 12a, it could be seen that the droplets were substantially deformed at location (1) and there were abundant evidence of droplet breakup and retraction at location (2). Both phenomenon on the SEM images agreed very well with the above simulations.
(a) parallel section (b) vertical section Figure 11. Morphologies of the droplets in shear layer at position (1) of case II part. Figure 11. Morphologies of the droplets in shear layer at position (1) of case II part.   Figure 13 showed the curve of shear stress and interfacial tension of the two droplets over time during the mold filling process. The combined action of shear stress and interfacial tension were key factors determining the droplets morphologies. Droplet A was subjected to greater shear stress than droplet B throughout mold filling. This led to larger deformation of droplet A than droplet B. On the other hand, according to Young-Laplace law, larger droplet deformation meant that the interfacial tension experienced by droplet A was also greater. With the weakening of shearing stress and the strengthening of interfacial tension, the deformation of the droplets gradually reached the maximum degree when shear stress and interfacial tension stroke a balance. Whether the droplets broke up or not, both droplets retracted more or less. The different interfacial tension and shear stress of the droplets, which was caused by their distinct trajectories in the flow field, explained the various morphology evolution of droplets.

Force Analysis of Droplets
(a) droplet A ( b) droplet B Figure 13. Shear stress and interfacial tension during mold filling. Figures 14 and 15 showed the simulation results of the morphology evolution of droplet A along its trajectory under different injection rates (case I and III). Three snapshots of the droplets when they were at locations (1), (2) and (3) were presented. Compared with the droplet morphology evolution of case II, the droplet was more deformed  Figure 13 showed the curve of shear stress and interfacial tension of the two droplets over time during the mold filling process. The combined action of shear stress and interfacial tension were key factors determining the droplets morphologies. Droplet A was subjected to greater shear stress than droplet B throughout mold filling. This led to larger deformation of droplet A than droplet B. On the other hand, according to Young-Laplace law, larger droplet deformation meant that the interfacial tension experienced by droplet A was also greater. With the weakening of shearing stress and the strengthening of interfacial tension, the deformation of the droplets gradually reached the maximum degree when shear stress and interfacial tension stroke a balance. Whether the droplets broke up or not, both droplets retracted more or less. The different interfacial tension and shear stress of the droplets, which was caused by their distinct trajectories in the flow field, explained the various morphology evolution of droplets.   Figure 13 showed the curve of shear stress and interfacial tension of the two droplets over time during the mold filling process. The combined action of shear stress and interfacial tension were key factors determining the droplets morphologies. Droplet A was subjected to greater shear stress than droplet B throughout mold filling. This led to larger deformation of droplet A than droplet B. On the other hand, according to Young-Laplace law, larger droplet deformation meant that the interfacial tension experienced by droplet A was also greater. With the weakening of shearing stress and the strengthening of interfacial tension, the deformation of the droplets gradually reached the maximum degree when shear stress and interfacial tension stroke a balance. Whether the droplets broke up or not, both droplets retracted more or less. The different interfacial tension and shear stress of the droplets, which was caused by their distinct trajectories in the flow field, explained the various morphology evolution of droplets.  Figures 14 and 15 showed the simulation results of the morphology evolution of droplet A along its trajectory under different injection rates (case I and III). Three snapshots of the droplets when they were at locations (1), (2) and (3) were presented. Compared with the droplet morphology evolution of case II, the droplet was more deformed  Rate   Figures 14 and 15 showed the simulation results of the morphology evolution of droplet A along its trajectory under different injection rates (case I and III). Three snapshots of the droplets when they were at locations (1), (2) and (3) were presented. Compared with the droplet morphology evolution of case II, the droplet was more deformed and stretched into a worm-like shape under higher injection rate. Similarly, the droplet broke up more intensely and split into more than two sub-droplets of case II. It was noteworthy that at location (2) the middle sub-droplet even had a tendency of secondary breakup under higher injection rates. However, under lower injection rate, the droplet only deformed to a certain degree without any breakup. and stretched into a worm-like shape under higher injection rate. Similarly, the droplet broke up more intensely and split into more than two sub-droplets of case II. It was noteworthy that at location (2) the middle sub-droplet even had a tendency of secondary breakup under higher injection rates. However, under lower injection rate, the droplet only deformed to a certain degree without any breakup.  (1), the droplets were more deformed under higher injection rates than that under lower injection rates. The biggest difference between the droplets morphologies in Figures 16 and 17 was at location (3), where larger amount of subdroplets than that under lower injection rates, a sign of droplet breakup, were spotted in both SEM images. It suggested that higher injection rates significantly promoted deformation and breakup of droplets and agreed well with the simulation results (Figures 14c and 15c). and stretched into a worm-like shape under higher injection rate. Similarly, the droplet broke up more intensely and split into more than two sub-droplets of case II. It was noteworthy that at location (2) the middle sub-droplet even had a tendency of secondary breakup under higher injection rates. However, under lower injection rate, the droplet only deformed to a certain degree without any breakup.  (1), the droplets were more deformed under higher injection rates than that under lower injection rates. The biggest difference between the droplets morphologies in Figures 16 and 17 was at location (3), where larger amount of subdroplets than that under lower injection rates, a sign of droplet breakup, were spotted in both SEM images. It suggested that higher injection rates significantly promoted deformation and breakup of droplets and agreed well with the simulation results (Figures 14c and 15c).  (1) and (3), specifically 1.0 mm from the surface, corresponding to droplet A of the simulation. The SEM experiments results revealed the impact of injection rate on droplet morphology which were consistent with the simulation. At location (1), the droplets were more deformed under higher injection rates than that under lower injection rates. The biggest difference between the droplets morphologies in Figures 16 and 17 was at location (3), where larger amount of sub-droplets than that under lower injection rates, a sign of droplet breakup, were spotted in both SEM images. It suggested that higher injection rates significantly promoted deformation and breakup of droplets and agreed well with the simulation results (Figures 14c and 15c).    Here the droplets were assumed axisymmetric and ellipsoidal while the semi-major and semi-minor axes of the droplets, L and W, were measured on the SEM images using the image analysis program Nano measurer. Then the equivalent radius R and deformability D of the droplets were given by:

Effect of the Injection
With the equivalent radius and deformability of all the droplets calculated, the average radius and deformability and their standard deviation of the droplets at position (1) and (3) of the parts molded under high and low injection rates were obtained using the following formula:  Here the droplets were assumed axisymmetric and ellipsoidal while the semi-major and semi-minor axes of the droplets, L and W, were measured on the SEM images using the image analysis program Nano measurer. Then the equivalent radius R and deformability D of the droplets were given by: With the equivalent radius and deformability of all the droplets calculated, the average radius and deformability and their standard deviation of the droplets at position (1) and (3) of the parts molded under high and low injection rates were obtained using the following formula: From the statistics it is straightforward that higher injection rate resulted in smaller size and larger deformability of droplets on average, while along the trajectory from location (1) to (3) the averaged deformability decreased but the averaged radius increased. It is the trade-off result of droplet retraction and deformation. It also suggested that coalescence of droplets was promoted when the shear action weakened. The statistical results of the simulation further validated the proposed model.
From the statistics it is straightforward that higher injection rate resulted in smaller size and larger deformability of droplets on average, while along the trajectory from location (1) to (3) the averaged deformability decreased but the averaged radius increased. It is the trade-off result of droplet retraction and deformation. It also suggested that coalescence of droplets was promoted when the shear action weakened. The statistical results of the simulation further validated the proposed model.

Conclusions and Outlook
In this study, a novel model for the droplet morphology evolution during the injection molding of polymer blends was proposed using the multiscale approach of the FVM, LBM and local scale-coupling method. The macro-and meso-scales were successfully bridged by tracking droplet trajectory and developing the constitutive model of polymer blends. As the cornerstone of the multiscale model, the macroscopic simulation of the mold filling flow was verified by the mold cavity pressure experiment and the simulation results precisely agreed with the pressure sensors' data. On the mesoscale, the model accurately tracked the trajectory of the droplet and simulated the subtle morphology evolution of the droplet on its trajectory. The simulation results and SEM experiments agreed very well so that the mesoscopic simulation of droplet morphology was validated. Through the numerical force analysis of the droplet, it is found that the shear stress sharply increased in the vicinity of the injection gate, resulting in abrupt deformation of the droplet. As the droplet moved away from the gate, the shear stress fell drastically while the interfacial tension gradually gained dominance, making the droplet tend to re-

Conclusions and Outlook
In this study, a novel model for the droplet morphology evolution during the injection molding of polymer blends was proposed using the multiscale approach of the FVM, LBM and local scale-coupling method. The macro-and meso-scales were successfully bridged by tracking droplet trajectory and developing the constitutive model of polymer blends. As the cornerstone of the multiscale model, the macroscopic simulation of the mold filling flow was verified by the mold cavity pressure experiment and the simulation results precisely agreed with the pressure sensors' data. On the mesoscale, the model accurately tracked the trajectory of the droplet and simulated the subtle morphology evolution of the droplet on its trajectory. The simulation results and SEM experiments agreed very well so that the mesoscopic simulation of droplet morphology was validated. Through the numerical force analysis of the droplet, it is found that the shear stress sharply increased in the vicinity of the injection gate, resulting in abrupt deformation of the droplet. As the droplet moved away from the gate, the shear stress fell drastically while the interfacial tension gradually gained dominance, making the droplet tend to retract. Finally, the impact of injection rate on the droplet morphology evolution was investigated and the simulation results qualitatively agreed with the SEM observations, which further validated the proposed model.
As for the outlook for future research, there are two issue worthy of attention. Firstly, the POE fraction of the blends in this paper is fixed 25% for the sake of subsequent manufacturing quality. The impact of blends fraction on the morphology evolution could be further studied. Secondly, fountain flow is typical of injection molding of polymers.
However, the morphology evolution of the droplets in the flow front is not shown in this paper due to the difficulty of experimental observation. To overcome this thorny problem, it requires novel experiment design.
In summary, this paper has achieved some pioneering success in multiscale simulation of injection molding of polymer blends, much work needs doing and it is firmly believed that this field is interesting and meaningful. Data Availability Statement: Data available on request due to restrictions eg privacy or ethical. The data presented in this study are available on request from the corresponding author. The data are not publicly available due to research conditions.

Conflicts of Interest:
The authors declare no conflict of interest.