Investigation on the Coupling Effects between Flow and Fibers on Fiber-Reinforced Plastic (FRP) Injection Parts

Glass or carbon fibers have been verified that can enhance the mechanical properties of the polymeric composite injection molding parts due to their orientation distribution. However, the interaction between flow and fiber is still not fully understood yet, especially for the flow–fiber coupling effect. In this study, we have tried to investigate the flow–fiber coupling effect on fiber reinforced plastics (FRP) injection parts utilizing a more complicated geometry system with three ASTM D638 specimens. The study methods include both numerical simulation and experimental observation. Results showed that in the presence of flow–fiber coupling effect, the melt flow front advancement presents some variation, specifically the “convex-flat-flat” pattern will change to a “convex-flat-concave” pattern. Furthermore, through the fiber orientation distribution (FOD) study, the flow–fiber coupling effect is not significant at the near gate region (RG). It might result from the strong shear force to repress the appearance of the flow–fiber interaction. However, at the end of filling region (ER), the flow–fiber coupling effect tries to diminish the flow direction orientation tensor component A11 and enhance the cross-flow orientation tensor component A22 simultaneously. It results in the dominance in the cross-flow direction at the ER. This orientation distribution behavior variation has been verified using a micro-computerized tomography (micro-CT) scan and image analysis technology.


Introduction
As the human population keeps increasing dramatically, CO 2 emissions from transportation based on the burning of fossil fuels will continue to be a major contributor to air pollution [1][2][3][4]. To decrease the CO 2 emissions from transportation, one of the most effective methods is using lightweight materials, such as fiber-reinforced plastics (FRP), to enhance the fuel efficiency [5][6][7]. To understand how fibers can enhance the mechanical strength of FRP products, many researchers have tried to study their features experimentally [8][9][10][11]. However, there are at least two major challenges which people still encounter in FRP injection molding. Specifically, one is what the fiber features are under the influence of the polymer matrix flow dynamic, and the other is how the fiber will affect the polymer matrix flow dynamic and generate the coupling between polymer flow and fibers.
One of the key fiber features in FRP injection molding is the flow-induced fiber orientation behavior. Folger and Tucker [12], Advani and Tucker [13], and Advani [14] have proposed some numerical models to predict fiber orientation in short fiber reinforced plastics injection molding. These models have given great guidance for people to understand fiber orientation behavior and

Models for Polymer Fluid Mechanics
The polymeric melt is considered as compressible, non-Newtonian fluids. The governing equations for the polymer fluid mechanics, which are regarded as 3D transient non-isothermal motion, are: where ρ is density; u is velocity vector; t is time; σ is total stress tensor; τ is extra stress tensor; g is the acceleration vector of gravity; P is pressure; C P is specific heat; T is temperature; k is thermal conductivity; D is the rate-of-deformation tensor. For the polymer melt, the extra stress tensor without flow-fiber coupling can be expressed as: where η is the shear viscosity of a polymer melt; it is a function of temperature and shear rate as described in Equations (6) and (7). Moreover, one of the most important properties to influence the flow mechanics is viscosity, which is a strong function of temperature and shear rate. To estimate the shear rate and temperature dependence of viscosity, the viscosity characterization has been performed using a rheometer from CoreTech System (Moldex3D) Co. Ltd. Furthermore, the measured data has been executed with curve-fitting using various viscosity models. One of the best fitted is based on the modified-Cross model. Hence, the modified-Cross model with Arrhenius temperature dependence is employed to describe the viscosity of a polymer melt in this study: where n is the power law index; η 0 is the zero shear viscosity; τ* is the parameter that describes the transition region between zero shear rate and the power law region of the viscosity curve; B is constant; T b is a reference temperature.

Models for Fiber Orientation Kinetics
During the injection molding processing, the orientation of the fiber inside the FRP composite will be affected by the polymer fluid mechanics. In general, the orientation of a single fiber can be at an arbitrary direction, or at one of the three principal directions, as shown in Figure 1. However, since there are hundreds of thousands of fibers inside the polymer matrix in injection molding, it is too difficult to observe manually. Instead, to predict the fiber orientation kinetics, each single fiber is theoretically regarded as an axisymmetric bond with rigidness. The bond's unit vector p along its axis direction can be described as the fiber orientation. The orientation state of a group of fibers is given by the second-order orientation tensor, where ψ(p) is the probability density distribution function over the orientation space. Tensor A 4 is a fourth order orientation tensor, defined as: where this tensor is also symmetric. Among the nine components of A ij , only five are independent. Three of the remaining components are determined by symmetry: The last component is determined by The acceptable calculation is obtained through the eigenvalue-based optimal fitting approximation of the orthotropic closure family.
Polymers 2020, 12, x FOR PEER REVIEW 4 of 17 where ) (p  is the probability density distribution function over the orientation space. Tensor A4 is a fourth order orientation tensor, defined as: where this tensor is also symmetric. Among the nine components of Aij, only five are independent.
Three of the remaining components are determined by symmetry: The last component is determined by The acceptable calculation is obtained through the eigenvalue-based optimal fitting approximation of the orthotropic closure family. To handle this complicated tensor system, Tseng et al. [21][22][23][24] extended the ARD-RSC models [18][19][20] to develop a fiber orientation model to couple with Jeffery's hydrodynamic (HD) model, namely, the iARD-RPR model (known as the Improved Anisotropic Rotary Diffusion model combined with the Retarding Principal Rate model), where ̇ represents the material derivative of A. Parameters CI and CM describe the fiber-fiber interaction and fiber-matrix interaction, while parameter α can slow down the response of fiber orientation. The Jeffery's hydrodynamic (HD) term can be written as To handle this complicated tensor system, Tseng et al. [21][22][23][24] extended the ARD-RSC models [18][19][20] to develop a fiber orientation model to couple with Jeffery's hydrodynamic (HD) model, namely, the iARD-RPR model (known as the Improved Anisotropic Rotary Diffusion model combined with the Retarding Principal Rate model), A represents the material derivative of A. Parameters C I and C M describe the fiber-fiber interaction and fiber-matrix interaction, while parameter α can slow down the response of fiber orientation. The Jeffery's hydrodynamic (HD) term can be written as where W and D are the vorticity tensor and the rate-of-deformation tensor, respectively. ξ is a shape factor of a particle. The rest of the details for the RPR model and the iARD model are available elsewhere [29].

Models for Evolution of Viscosity by Flow-Fiber Coupling
It was mentioned earlier that the fibers will affect the polymer fluid mechanics, and flow will further influence fibers back and forth; it introduces the coupling effect between the flow and fibers. To predict this flow-fiber coupling, the revised IISO constitutive equation has developed to handle the complex flow-fiber coupling by Tseng and Favaloro [29]. This revised IISO constitutive equation has been implemented into the commercial injection molding software, Moldex3D. During the injection molding simulation, the governing equations of flow field and fiber orientation are solved by the 3D finite volume method based on 3D geometry. After the velocity field and the fiber orientation filed are obtained at the beginning, they are used to determine the revised IISO viscosity. Then the iteration keeps going till they are converged. The details are available elsewhere [29]. Specifically, the revised IISO viscosity is presented as below.
where η S the nonlinear Newtonian viscosity for the fiber-filled polymer fluids and is described by the modified-Cross model; R T is the dimensionless Trouton ratio parameter as a function of the strain rate; R 0 T is the initial value of R T ; K S is a stretching kernel which is related to the flow fields and to the fiber orientation state; . γ C is the critical strain rate (1/s).

Simulation Model and Related Information
The geometrical model is shown as in Figure 2. Specifically, the dimension of the whole model is 400 mm × 165 mm × 3 mm. There are three standard specimens based on ASTM D638, marked as Model I, II, and III, respectively. The dimension for each standard specimen is 172 mm × 20 mm × 3 mm. To study the flow field behavior, we have designed those three models with different gate types. In particular, one has an edge gate, another has a sprue gate, and the other has a double gate, as presented. To fill the polymer melt into those three specimens, there is one plunger type melt entrance with 40 mm diameter that is located in the center of the cavity. Moreover, the moldbase and cooling channel layout are presented in Figure 3.
Furthermore, several measuring nodes in different models are specified to study the fiber orientation distribution (FOD), as in Figure 4. In particular, they are specified as three regions for observation for each model: the gate region (GR), the center region (CR), and the end of filling region (ER), respectively. Moreover, to perform the injection molding simulation, the associated process condition setting is as follows. Specifically, the filling time is 1.49 s, the packing time is 5 s, the cooling time is 15 s, and other parameters are in Table 1. The materials used in this study are pure polypropylene (later called PP) and polypropylene with short fiber of 3 mm length (later called PP + SF). Both PP and PP + SF materials are commercially available and supplied by LCY Chemical. Moreover, the associated material properties including viscosity, the specific volume against temperature and pressure (pvT), heat capacity (Cp), and thermal conductivity (K) have been measured and saved as material data by Moldex3D. In addition, to study the fiber orientation effect during the injection molding, the iARD-RPR model is adopted. In general, those three parameters are suggested in 0 < C I < 1, 0 < C M < 1, and 0 < α < 1 [29]. Here, the related parameters for the iARD-RPR model in this study are suggested from Moldex3D, as listed in Table 2. Furthermore, to conduct the flow-fiber effect, the revised IISO model is utilized. The associated parameters are also suggested from Moldex3D, which are presented in Table 3.

Experimental Equipment and Related Information
Figure 5a displays the injection machine (Chuan Lih Fa Machinery Works Co., Tainan, Taiwan Model: CLF-180TXL). Figure 5b shows the real mold structure and cooling channel layout. The PP and PP + SF materials are supplied by LCY Chemical as mentioned earlier. Specifically, the grade name: PP is Globalene ST868M and PP + SF is Globalene SF7351. Moreover, to perform the real injection molding testing, the associated process condition setting is the same as that of the simulation, as listed in Table 1. In order to validate the flow-fiber coupling effect calculated from numerical simulation, the mechanical properties via tensile test for the real injected parts have been executed using the universal tensile machine from Shimadzu, Kyoto, Japan (AGS-J mode). In addition, to validate the fiber orientation behavior, micro-computerized tomography (μ-CT) technology has been performed using Bruker Skyscan 2211 with 40-190 kV and a resolution of 5 μm, supported by Material and Chemical Laboratories (MCL) Multiscale X-ray CT laboratory, Industrial   Table 2. The parameters for iARD-RPR model in injection molding simulation.  PP is Globalene ST868M and PP + SF is Globalene SF7351. Moreover, to perform the real injection molding testing, the associated process condition setting is the same as that of the simulation, as listed in Table 1. In order to validate the flow-fiber coupling effect calculated from numerical simulation, the mechanical properties via tensile test for the real injected parts have been executed using the universal tensile machine from Shimadzu, Kyoto, Japan (AGS-J mode). In addition, to validate the fiber orientation behavior, micro-computerized tomography (µ-CT) technology has been performed using Bruker Skyscan 2211 with 40-190 kV and a resolution of 5 µm, supported by Material and Chemical Laboratories (MCL) Multiscale X-ray CT laboratory, Industrial Technology Research Institute, Hsinchu, Taiwan.

Results and Discussion
To understand what the coupling effect is between flow and fiber in the injection molding process, we have tried to observe the mechanical property (macroscopic property) in the presence of fibers for FRP first. Furthermore, the fiber orientation behavior (microscopic feature) has been discovered. Through the relationship between mechanical properties and fiber orientation behavior, the flow-fiber coupling effect can be further discussed. The details are as follows.

Mechanical Property Test
To conduct the fiber reinforced capability, tensile testing has been performed using the universal tensile testing machine of Shimadzu (AGS-J mode). The testing method and procedures are based on the Ref. [31]. During the testing, the deformation and the force were recorded. The deformation and the force can be converted to the stress and strain using the following equations: where σσσ is tensile stress; F is the pulling force on the original cross-section area (A0) of the narrow portion of the standard specimen; εis the elongation ratio; (L0) is the original length of the narrow portion; L is deformation; E is the tensile modulus. For each tensile testing, those three models of ASTM D638 standard specimens have been installed into the tensile machine under a constant strain at 20 mm/min. For each model, five specimens have been executed for the same testing. After the five tests, the average stress-strain relation is presented in Figure 6. From this stress-strain relation in Figure 6a, the tensile modulus is calculated around the strain range of 2.5-3.5% to avoid the influence from the testing apparatus. The

Results and Discussion
To understand what the coupling effect is between flow and fiber in the injection molding process, we have tried to observe the mechanical property (macroscopic property) in the presence of fibers for FRP first. Furthermore, the fiber orientation behavior (microscopic feature) has been discovered. Through the relationship between mechanical properties and fiber orientation behavior, the flow-fiber coupling effect can be further discussed. The details are as follows.

Mechanical Property Test
To conduct the fiber reinforced capability, tensile testing has been performed using the universal tensile testing machine of Shimadzu (AGS-J mode). The testing method and procedures are based on the Ref. [31]. During the testing, the deformation and the force were recorded. The deformation and the force can be converted to the stress and strain using the following equations: where σ σσ is tensile stress; F is the pulling force on the original cross-section area (A 0 ) of the narrow portion of the standard specimen; ε is the elongation ratio; (L 0 ) is the original length of the narrow portion; L is deformation; E is the tensile modulus. For each tensile testing, those three models of ASTM D638 standard specimens have been installed into the tensile machine under a constant strain at 20 mm/min. For each model, five specimens have been executed for the same testing. After the five tests, the average stress-strain relation is presented in Figure 6. From this stress-strain relation in Figure 6a, the tensile modulus is calculated around the strain range of 2.5-3.5% to avoid the influence from the testing apparatus. The tensile modulus are obtained as 1463.67 N/mm 2 for Model I, 1427.33 N/mm 2 for Model II, 1299.25 N/mm 2 for Model III, respectively. Based on the tensile modulus, Model I is stronger than Model II, while Model III is weakest due to the double gate design. Moreover, the tensile strength is also conducted for those three models as listed in Figure 6b. Additionally, the tensile strength for each testing has been recorded as in Table 4. Specifically, when material is changed from pure PP to PP + SF, the tensile strength of Model I is increased from 23.35 to 75.71 N/mm 2 . It is increased from 23.04 to 73.34 N/mm 2 for Model II, and from 22.65 to 34.38 N/mm 2 for Model III. Although the tensile modulus and tensile strength of the Model I and Model II are close, the actual measured data (including stress-strain curve, stress modulus, and tensile strength) of Model I are greater than that of Model II, as shown in Figure 6 and Table 4. Overall, the mechanical properties of Model I are stronger than those of Model II, while Model III's are weakest due to the double gate design.
Polymers 2020, 12, x FOR PEER REVIEW 9 of 17 modulus and tensile strength of the Model I and Model II are close, the actual measured data (including stress-strain curve, stress modulus, and tensile strength) of Model I are greater than that of Model II, as shown in Figure 6 and Table 4. Overall, the mechanical properties of Model I are stronger than those of Model II, while Model III's are weakest due to the double gate design.

Fiber Orientation Effects With and Without Flow-Fiber Coupling
To realize why the mechanical properties of Model I are stronger than those of Model II and then to discover the flow-fiber coupling effect, the fiber orientation features have been utilized. The fiber orientation tensor components are predicted numerically by using Moldex3D R16 ® . They are also validated experimentally by using a micro-computerized tomography (micro-CT) scan of the specimens and utilizing Avizo ® for image analyses. The details are as follows. Figure 7a shows the fiber orientation distribution (FOD) predicted numerically without considering flow-fiber coupling at Point B for Model I. The flow direction orientation tensor component A 11 is increased from 0.7 to 0.78 at the skin-layer then decreased to 0.6 at the core-layer. Overall, the fiber orientation is dominant along the flow direction at gate region (GR). Moreover, Figure 7c shows the FOD observation experimentally. Except at little skin region, the flow direction orientation tensor component A 11 is around 0.8 to 0.9 at the skin-layer then decreased to 0.8 at the core-layer. It also shows that the dominant fiber orientation is at the flow direction in experiments. The simulation prediction of the FOD without considering flow-fiber coupling is in reasonable agreement with that of the experimental result for Model I at GR. Moreover, the simulation results of the fiber orientation tensor components are very close for both without and with flow-fiber coupling for Model I at GR, as shown in Figure 7a,b. Their trends are all consistent with the experimental observation. Since there is a strong entrance effect around GR [30], the fibers are influenced significantly by the higher shear rate (shear force) that will repress the flow-fiber coupling effect.

Fiber Orientation Effects With and Without Flow-Fiber Coupling
To realize why the mechanical properties of Model I are stronger than those of Model II and then to discover the flow-fiber coupling effect, the fiber orientation features have been utilized. The fiber orientation tensor components are predicted numerically by using Moldex3D R16 ® . They are also validated experimentally by using a micro-computerized tomography (micro-CT) scan of the specimens and utilizing Avizo ® for image analyses. The details are as follows. Figure 7a shows the fiber orientation distribution (FOD) predicted numerically without considering flow-fiber coupling at Point B for Model I. The flow direction orientation tensor component A11 is increased from 0.7 to 0.78 at the skin-layer then decreased to 0.6 at the core-layer. Overall, the fiber orientation is dominant along the flow direction at gate region (GR). Moreover, Figure 7c shows the FOD observation experimentally. Except at little skin region, the flow direction orientation tensor component A11 is around 0.8 to 0.9 at the skin-layer then decreased to 0.8 at the core-layer. It also shows that the dominant fiber orientation is at the flow direction in experiments. The simulation prediction of the FOD without considering flow-fiber coupling is in reasonable agreement with that of the experimental result for Model I at GR. Moreover, the simulation results of the fiber orientation tensor components are very close for both without and with flow-fiber coupling for Model I at GR, as shown in Figure 7a,b. Their trends are all consistent with the experimental observation. Since there is a strong entrance effect around GR [30], the fibers are influenced significantly by the higher shear rate (shear force) that will repress the flow-fiber coupling effect. Moreover, considering the center region (CR), Figure 8a shows the FOD predicted numerically at Point E for Model I in the absence of flow-fiber coupling. The flow direction orientation tensor component A11 is increased from 0.7 to 0.85 at the skin-layer then decreased to 0.75 at the core-layer. Overall, the fiber Moreover, considering the center region (CR), Figure 8a shows the FOD predicted numerically at Point E for Model I in the absence of flow-fiber coupling. The flow direction orientation tensor component A 11 is increased from 0.7 to 0.85 at the skin-layer then decreased to 0.75 at the core-layer. Overall, the fiber orientation is still dominant along the flow direction at CR. In addition, Figure 8c shows the FOD observation experimentally. The flow direction orientation tensor component A 11 is around 0.6 to 0.95 at the skin-layer then decreased to 0.85 at the core-layer. It also shows that the dominant fiber orientation is at the flow direction in experiments. Once again, the simulation prediction of the FOD in the absence of flow-fiber coupling is in a reasonable agreement with that of the experimental result for Model I at CR. Moreover, the simulation results of the fiber orientation tensor components for both without and with flow-fiber coupling for Model I at CR are very similar, as shown in Figure 8a,b. Due to the stronger shear rate (shear force) by the convergent section, the flow-fiber coupling is not significant at this region. Furthermore, in the end of filling region (ER), Figure 9a shows the FOD predicted numerically without considering flow-fiber coupling at Point H for Model I. The flow direction orientation tensor component A11 is increased from 0.65 to 0.85 at the skin-layer then decreased to 0.5 at the core-layer. Overall, the fiber orientation is still dominant along the flow direction at ER. Moreover, Figure 9c shows the FOD observation experimentally. The flow direction orientation tensor component A11 is around 0.65 to 0.8 at the skin-layer then decreased dramatically to 0.4 at the core-layer. The fiber orientation is dominant at the flow direction A11 at the skin-layer and then changes to the cross flow direction A22 at the core-layer. From the experimental observation, the fiber orientation is strong in the cross-flow direction. It is noted that the simulation prediction without flow-fiber coupling over-predicted A11 and under-predicted A22 for Model I at ER. In addition, the width of the core-layer from the experiment is wider than that of the numerical prediction. Moreover, Figure 9b shows the simulation result under the influence of the flow-fiber coupling. Compared to the experimental observation, the simulation result with the flow-fiber coupling is much closer to the experimental observation. It is noted that during the melt flowing into the ER, the influence of the transport driving force is weaker and then the flow-fiber coupling becomes more significant. It could be the reason why the A11 is over-predicted and A22 is under-predicted by numerical simulation without considering the flow-fiber coupling effect. Furthermore, in the end of filling region (ER), Figure 9a shows the FOD predicted numerically without considering flow-fiber coupling at Point H for Model I. The flow direction orientation tensor component A 11 is increased from 0.65 to 0.85 at the skin-layer then decreased to 0.5 at the core-layer. Overall, the fiber orientation is still dominant along the flow direction at ER. Moreover, Figure 9c shows the FOD observation experimentally. The flow direction orientation tensor component A 11 is around 0.65 to 0.8 at the skin-layer then decreased dramatically to 0.4 at the core-layer. The fiber orientation is dominant at the flow direction A 11 at the skin-layer and then changes to the cross flow direction A 22 at the core-layer. From the experimental observation, the fiber orientation is strong in the cross-flow direction. It is noted that the simulation prediction without flow-fiber coupling over-predicted A 11 and under-predicted A 22 for Model I at ER. In addition, the width of the core-layer from the experiment is wider than that of the numerical prediction. Moreover, Figure 9b shows the simulation result under the influence of the flow-fiber coupling. Compared to the experimental observation, the simulation result with the flow-fiber coupling is much closer to the experimental observation. It is noted that during the melt flowing into the ER, the influence of the transport driving force is weaker and then the flow-fiber coupling becomes more significant. It could be the reason why the A 11 is over-predicted and A 22 is under-predicted by numerical simulation without considering the flow-fiber coupling effect.
influence of the flow-fiber coupling. Compared to the experimental observation, the simulation result with the flow-fiber coupling is much closer to the experimental observation. It is noted that during the melt flowing into the ER, the influence of the transport driving force is weaker and then the flow-fiber coupling becomes more significant. It could be the reason why the A11 is over-predicted and A22 is under-predicted by numerical simulation without considering the flow-fiber coupling effect. Moreover, Figure 10a presents the FOD predicted numerically without consideration of flowfiber coupling at Point B for Model II at GR. The flow direction orientation tensor component A11 is Moreover, Figure 10a presents the FOD predicted numerically without consideration of flow-fiber coupling at Point B for Model II at GR. The flow direction orientation tensor component A 11 is 0.82 at the skin-layer then decreased to 0.4 at the core-layer. The fiber orientation is dominant at the flow direction A 11 at the skin-layer and then changes to the cross flow direction A 22 at the core-layer. Similarly, Figure 10c exhibits the FOD observation experimentally at Point B for Model II. It shows that the simulation prediction of the FOD without consideration of flow-fiber coupling is consistent with that of the experimental result for Model II at GR in the trend. However, the width of the core-layer by experimental observation is wider than that of simulation prediction. Furthermore, in Figure 10a,b, there is some difference on the simulation results of the fiber orientation tensor components between both without and with flow-fiber coupling cases for Model II at GR. Specifically, under the influence of flow-fiber coupling, more fibers will stay at the cross flow direction A 22 at the core-layer, as shown in Figure 10b. The trend of the coupled case is more consistent with the experimental observation, as shown in Figure 10c. Similarly, Figure 10c exhibits the FOD observation experimentally at Point B for Model II. It shows that the simulation prediction of the FOD without consideration of flow-fiber coupling is consistent with that of the experimental result for Model II at GR in the trend. However, the width of the corelayer by experimental observation is wider than that of simulation prediction. Furthermore, in Figure  10a,b, there is some difference on the simulation results of the fiber orientation tensor components between both without and with flow-fiber coupling cases for Model II at GR. Specifically, under the influence of flow-fiber coupling, more fibers will stay at the cross flow direction A22 at the core-layer, as shown in Figure 10b. The trend of the coupled case is more consistent with the experimental observation, as shown in Figure 10c. Furthermore, Figure 11a presents the FOD predicted numerically in absence of the flow-fiber coupling at Point E for Model II. The flow direction orientation of tensor component A11 is 0.83 at the skin-layer then decreased to 0.75 at the core-layer. The fiber orientation is dominant at the flow direction A11. Additionally, Figure 11c exhibits the FOD observation experimentally at Point E for Furthermore, Figure 11a presents the FOD predicted numerically in absence of the flow-fiber coupling at Point E for Model II. The flow direction orientation of tensor component A 11 is 0.83 at the skin-layer then decreased to 0.75 at the core-layer. The fiber orientation is dominant at the flow direction A 11 . Additionally, Figure 11c exhibits the FOD observation experimentally at Point E for Model II. The flow direction orientation tensor component A 11 is 0.7 to 0.95 at the skin-layer, then decreases to 0.7 at the core-layer. Although the fiber orientation is still dominant at the flow direction A 11 , the fiber orientation prediction at the core-layer, predicted numerically without considering flow-fiber coupling, is over-predicted. Moreover, Figure 11a,b show the simulation results of the fiber orientation tensor components for both without and with flow-fiber coupling cases for Model II at CR. Due to the stronger shear rate (shear force) by the convergent section, flow-fiber coupling is not significant at this region except for the core-layer area. However, in the presence of flow-fiber coupling, the numerical prediction is closer to the experimental observation. Furthermore, Figure 11a presents the FOD predicted numerically in absence of the flow-fiber coupling at Point E for Model II. The flow direction orientation of tensor component A11 is 0.83 at the skin-layer then decreased to 0.75 at the core-layer. The fiber orientation is dominant at the flow direction A11. Additionally, Figure 11c exhibits the FOD observation experimentally at Point E for Model II. The flow direction orientation tensor component A11 is 0.7 to 0.95 at the skin-layer, then decreases to 0.7 at the core-layer. Although the fiber orientation is still dominant at the flow direction A11, the fiber orientation prediction at the core-layer, predicted numerically without considering flow-fiber coupling, is over-predicted. Moreover, Figure 11a,b show the simulation results of the fiber orientation tensor components for both without and with flow-fiber coupling cases for Model II at CR. Due to the stronger shear rate (shear force) by the convergent section, flow-fiber coupling is not significant at this region except for the core-layer area. However, in the presence of flow-fiber coupling, the numerical prediction is closer to the experimental observation. Moreover, Figure 12a shows the FOD predicted numerically without flow-fiber coupling at Point H for Model II. The flow direction orientation of tensor component A 11 is increased from 0.7 to 0.8 at the skin-layer, then decreased to 0.5 at the core-layer. Overall, the fiber orientation is still dominant along the flow direction at ER. In addition, Figure 12c shows the FOD observation experimentally. The flow direction orientation tensor component A 11 is around 0.6 to 0.8 at the skin-layer, then decreased dramatically to 0.3 at the core-layer. The fiber orientation is dominant at the flow direction A 11 at the skin-layer and then changes to the cross flow direction A 22 at the core-layer. From the experimental observation, the fiber orientation is strong in the cross-flow direction at the core-layer. The simulation prediction in the absence of flow-fiber coupling over-predicted A 11 and under-predicted A 22 for Model II at ER. Additionally, the width of the core-layer from the experiment is significantly wider than that of numerical prediction. Moreover, Figure 12b shows the fiber orientation tensor component at the ER for Model II under the influence of the flow-fiber coupling effect. Compared to the experimental observation, the simulation with flow-fiber coupling is much closer to the experimental result. It is noted that during the melt flowing into the ER, the influence of the transport driving force of Model II is much weaker and then the flow-fiber coupling becomes much more significant compared to Model I cases. Overall, from Model I and Model II, it can be concluded that when the melt flows through the higher shear rate (shear force) regions, such as GR and CR of Model I, the flow-fiber coupling will be repressed. However, when the shear rate (shear force) of the melt is reduced gradually (such as at the ER of Model I and at the CR and ER of Model II), the flow-fiber coupling effect will become more significant. significant compared to Model I cases. Overall, from Model I and Model II, it can be concluded that when the melt flows through the higher shear rate (shear force) regions, such as GR and CR of Model I, the flowfiber coupling will be repressed. However, when the shear rate (shear force) of the melt is reduced gradually (such as at the ER of Model I and at the CR and ER of Model II), the flow-fiber coupling effect will become more significant.

Correlation between Mechanical Property and Fiber Micro-Feature
Moreover, comparing the orientation tensor component A11 of the experimental observation from GR via CR to ER for both Model I and Model II (from Figures 7c-12c), the A11 component at the central line is increased by the convergent structure and then gradually reduced by the divergent structure. The variation of the FOD is quite different for Model I and II. Specifically, from Figures 7c and 10c, the A11 of the Model I is much higher than that of Model II at GR. In addition, the A11 of the Model I is also higher than that of Model II at CR and at ER. Overall, the stronger alignment in flow direction presented in A11 of Model I is expected to provide stronger mechanical properties (including tensile modulus and tensile stress), which is consistent with that in mechanical property testing. The reason for this difference is due to the entrance effect that occurs when melt flows through the edge gate design. That entrance effect will further provide some flow-induced fiber orientation to the melt,

Correlation between Mechanical Property and Fiber Micro-Feature
Moreover, comparing the orientation tensor component A 11 of the experimental observation from GR via CR to ER for both Model I and Model II (from Figures 7c, 8c, 9c, 10c, 11c and 12c), the A 11 component at the central line is increased by the convergent structure and then gradually reduced by the divergent structure. The variation of the FOD is quite different for Model I and II. Specifically, from Figures 7c and 10c, the A 11 of the Model I is much higher than that of Model II at GR. In addition, the A 11 of the Model I is also higher than that of Model II at CR and at ER. Overall, the stronger alignment in flow direction presented in A 11 of Model I is expected to provide stronger mechanical properties (including tensile modulus and tensile stress), which is consistent with that in mechanical property testing. The reason for this difference is due to the entrance effect that occurs when melt flows through the edge gate design. That entrance effect will further provide some flow-induced fiber orientation to the melt, which will further enhance the tensile properties of Model I, as described in Ref. [30]. Moreover, from the comparison of the fiber orientation tensor components between the numerical simulation and experimental results, the numerical simulation, with consideration of flow-fiber coupling, is much closer to the experimental observation.

Melt Flow Behavior under the Influence of Flow-Fiber Coupling
Finally, during the injection molding processing, the molten plastics are delivered from the sprue via the runner and gate and flow into the cavity. The melt flow front behavior is presented in Figure 13. Clearly, when the flow-fiber coupling is considered, the melt flow front exhibits "convex-flat-concave" behavior, while it shows "convex-flat-flat" behavior without flow-fiber coupling. This phenomenon is due to the free surface advancing faster along the side walls of the cavity due to the flow-fiber coupling effect. This flow-fiber coupling-induced melt flow behavior is similar to that described in Ref. [29]. sprue via the runner and gate and flow into the cavity. The melt flow front behavior is presented in Figure 13. Clearly, when the flow-fiber coupling is considered, the melt flow front exhibits "convexflat-concave" behavior, while it shows "convex-flat-flat" behavior without flow-fiber coupling. This phenomenon is due to the free surface advancing faster along the side walls of the cavity due to the flow-fiber coupling effect. This flow-fiber coupling-induced melt flow behavior is similar to that described in Ref. [29].

Conclusions
In this study, we have tried to investigate the flow-fiber coupling effect on FRP injection parts using a more complicated geometry system with three ASTM D638 specimens through both simulation prediction and experimental observation. Results showed that in the presence of the flow-fiber coupling effect, the melt flow front advancement presents some variation, specifically the original "convex-flat-flat" pattern will change to a "convex-flat-concave" pattern. This observation is consistent with that of Tseng and Favaloro [30], in that it is a simple end-gated plate system. Furthermore, through the fiber orientation study, the flow-fiber coupling effect is not significant at the gate region (RG). It might result from the strong shear force to hold down the appearance of the flow-fiber coupling. However, at the end region (ER), since the shear force becomes lower, the flow-fiber coupling effect tries to diminish the flow direction orientation tensor component A 11 and enhance the cross-flow orientation tensor component A 22 simultaneously. It ends up with the cross-flow direction dominance at the ER. This fiber orientation distribution (FOD) behavior variation to cause the flow-fiber coupling behavior has been verified using micro-computerized tomography (µ-CT) scan and images analyses. This overall A 11 dominance of the FOD arrangement for Model I has provided stronger tensile properties that are consistent with the experimental observation via tensile test.