長玻纖增強復合材料注塑成型構(gòu)件強度分析

，，，，

(1 湖南工業(yè)大學 機械工程學院，湖南 株洲 412007；2 株洲時代新材料科技股份有限公司，湖南 株洲 412007)

隨著生態(tài)環(huán)境日益惡化，為實現(xiàn)節(jié)能減排、低碳環(huán)保，各國都提出了相應的措施，輕量化技術(shù)被認為是目前行之有效的方法和手段[1-2]。纖維增強復合材料是一種纖維增強的高分子材料，具有高比強度、高比模量等特點，是替代金屬的輕量化材料[2-3]。

纖維增強復合材料的力學性能除了受基體材料和纖維材料的物理屬性影響外，還與纖維的體積分數(shù)、纖維長徑比及纖維排列方式等有關(guān)[4-7]。Eshelby的等效夾雜理論[8]，自洽法[9], Mori-Tanaka方法[10]等能對纖維增強復合材料的模量進行預測。基于Jeffery方程[11]，F(xiàn)olgar等[12]針對濃縮的纖維懸浮液，考慮流體和纖維的相互作用，給出了纖維角速度模型。在此基礎(chǔ)上Wang 等[13]基于保持特征值的旋轉(zhuǎn)速率不變的情況下,通過標量閉合因子降低取向張量特征值的增長速率的理念，提出了RSC取向模型。隨后，Phelps等[14]用各向異性旋轉(zhuǎn)擴散(ARD) 張量項取代 RSC模型中的各向同性標量項，提出了ARD-RSC取向模型，它能有效預測較長纖維間相互作用增強后的纖維取向分布。

纖維增強復合材料由于其注塑結(jié)構(gòu)件的可設計性，應用日益廣泛，但因其在注塑成型過程中纖維取向的分散性而影響注塑件的力學性能和使用性能。如何準確預測塑件力學性能，對纖維增強復合材料實際應用十分重要。目前，在工程應用中，一般是將長玻纖增強復合材料作為各向同性材料處理，鮮有考慮其各向異性特征，即在不同的方向抵抗變形及失效的能力不同的特點的分析。針對長玻纖復合材料注塑結(jié)構(gòu)件成型過程中，纖維含量是可控因素，在纖維含量一定的情況下，纖維取向狀態(tài)是影響長玻纖增強復合材料力學性能的主要影響因素。本工作基于廣義牛頓流體本構(gòu)方程，采用ARD-RSC纖維取向模型，考慮纖維間相互作用，應用數(shù)值模擬預測注塑件的纖維取向分布；基于復合材料細觀力學的Eshelby夾雜理論，應用Mean Field均勻化方法，建立長玻纖增強復合材料均質(zhì)化RVE模型；綜合運用材料建模、離散均質(zhì)化RVE模型場、注塑成型和結(jié)構(gòu)有限元分析技術(shù)，對構(gòu)件強度進行分析，并在此基礎(chǔ)上針對纖維取向分布不理想對產(chǎn)品力學性能產(chǎn)生重大影響進行結(jié)構(gòu)改進。

1 本構(gòu)關(guān)系和纖維取向模型

1.1 正交各向異性材料彈性本構(gòu)關(guān)系

將長玻纖增強復合材料簡化為正交各向異性材料，則存在以下應力-應變關(guān)系[15]：

(1)

用矩陣符號表示為：

ε=Sσ

(2)

其中，Sij為柔度系數(shù)，S為柔度矩陣，柔度矩陣是剛度矩陣的逆矩陣；同樣，Sij具有對稱性，只有9個獨立柔度系數(shù)，剛度和柔度系數(shù)對均質(zhì)材料都可以認為是彈性常數(shù)。

正交各向異性材料的柔度矩陣S與9個彈性常數(shù)的關(guān)系如式(3)所示：

(3)

1.2 纖維取向理論模型

Bay等對纖維取向理論模型[16]進行了研究。在長玻纖增強PA66復合材料中，纖維一般假設為不可彎曲的棒狀體，則可用單位矢量p來描述單根纖維的取向，如圖1所示:

(4)

式中：p為纖維的方向；θ表示p與坐標軸3的夾角；φ表示p在1-2平面上的投影與坐標軸1的夾角。

圖1 纖維取向矢量Fig.1 Fiber orientation vector

一組纖維的取向狀態(tài)可以用概率分布函數(shù)ψ(p)表示。假定纖維在流體中流動，該分布函數(shù)的守恒方程[17]為：

(5)

對于長玻纖增強PA66復合材料，纖維分布具有隨機性。而對于某一特定區(qū)域內(nèi)分纖維取向分布情況，可用單根特定長度的纖維在空間隨機取向狀態(tài)來描述，即用纖維取向張量描述纖維取向狀態(tài)。定義單位矢量p的張量積在所有方向上的積分并點乘纖維分布函數(shù)ψ(p)：

(6)

其中，

(7)

式中：aijk…l為纖維的取向張量；pi,pj,pk…pl為第i,j,k…l根纖維的方向；分布函數(shù)Ω(p)為纖維出現(xiàn)在θ與θ+dθ，φ與φ+dφ之間的概率。

由于纖維隨流場一起運動的過程中沒有頭與尾的區(qū)別，則分布函數(shù)ψ(p)是偶函數(shù)，即ψ(p)=ψ(-p)，則奇數(shù)個p相乘的積分結(jié)果均為0。所以可以只用偶數(shù)個p取向矢量的張量積來定義取向張量。現(xiàn)以二階張量為例，aij=〈pi,pj〉是一個二階對稱張量[16]。 二階張量aij描述纖維取向狀態(tài)[16]由纖維分布函數(shù)定義為：

(8)

根據(jù)取向張量滿足對稱性及常規(guī)條件αij=αji，二階取向張量可寫成：

(9)

經(jīng)對角化后，描述纖維取向狀態(tài)采用如下二階張量的形式：

(10)

其中[16]，

λ1+λ2+λ3=1

(11)

經(jīng)對角化后的二階取向張量狀態(tài)描述，如圖2所示。圖2(a)表示單軸取向，即纖維沿單一方向排列；圖2(b)表示雙軸取向，既纖維沿2個主軸方向平均分布，同時也可表示纖維在平面中呈隨機取向分布；圖2(c)表示沿3個主軸方向的分布狀態(tài)或是空間中的隨機取向。

圖2 纖維取向張量(a)單軸取向;(b)二維隨機取向；(c)3D隨機取向Fig.2 Fiber orientation tensor(a)uniaxial aligned;(b)planar random biaxial;(c)3D random

進行注塑成型流動模擬，假設忽略熔體彈性，可采用廣義牛頓流體本構(gòu)方程[18]:

(12)

在注塑成型加工數(shù)值模擬中，利用取向分布函數(shù)能完整描述纖維取向狀態(tài)，但需要巨大的計算資源，目前多使用取向張量來描述纖維取向。在Folgar-Tucker模型[12]和RSC模型[13]的基礎(chǔ)上，Phelps等[14]用各向異性旋轉(zhuǎn)擴散(ARD) 張量項取代 RSC模型中的各向同性標量項，提出了ARD-RSC取向模型，它能有效預測較長纖維間相互作用增強后的纖維取向分布，演化而成的ARD-RSC模型方程如下所示：

2[A+(1-K)(L-M∶A)]∶D)+

5(C·a+a·C)+

10[A+(1-K)(L-M∶A)]∶C)

(13)

其中，

(14)

(15)

1.3 纖維均質(zhì)化模型

基于Eshelby等效夾雜理論[8]，Mori-Tanaka方法[10]假設復合材料的體積平均應力應等于其遠場作用的均勻應力，其等效彈性模量為：

(16)

其中，

(17)

(18)

為了有效預測纖維增強復合材料宏觀力學響應，建立細觀纖維均質(zhì)化RVE模型如圖3所示。首先根據(jù)Mori-Tanaka模型對每一個偽谷粒均質(zhì)化，然后使用Voigt模型[7]對所有偽谷粒進行集體均質(zhì)化計算。

圖3 纖維均質(zhì)化RVE模型Fig.3 Fiber homogenized RVE model

2 長玻纖增強復合材料注塑成型推力桿強度分析與結(jié)構(gòu)改進

2.1 推力桿的材料及性能

汽車推力桿主要用于保持車橋位置的相對固定，保證車輪和車身之間有確定的運動關(guān)系，使汽車具有良好的駕駛性能。推力桿中裝有橡膠球鉸，在汽車行駛過程中，橡膠球鉸頻繁受到拉壓，為提高球鉸壽命，在安裝球鉸時，要使球鉸處于預壓狀態(tài)。推力桿的材料為熱塑性長玻纖增強PA66(GF50)復合材料，由PA66尼龍樹脂基體和玻璃纖維增強相兩部分組成，兩相均為各向同性材料，其力學性能見表1。

表1 兩組分材料性能參數(shù)Table 1 Performance parameters on two kinds of material

參照GB/T 1447-2005 纖維增強塑料拉伸性能實驗方法，使用長玻纖增強PA66(GF50)復合材料I型注塑拉伸試樣，在JYW-93 Z100高低溫電子萬能材料試驗機上進行拉伸實驗，加載速率為2mm/min，測試試樣為5根，5組數(shù)據(jù)獲取的拉伸模量和泊松比見表2。試樣的纖維長度方向與拉伸方向一致，表中的拉伸強度是表征其近似纖維長度分布方向的強度。

表2 長玻纖增強PA66(GF50)復合材料拉伸性能Table 2 Tensile properties of long glass fiber reinforced PA66 (GF50) composites

2.2 推力桿臺架實驗

桿體設計質(zhì)量為2kg，為了節(jié)省開發(fā)費用，作者選取二分之一模型進行注塑加工，如圖4所示，質(zhì)量為1kg。在實際注塑成型加工中，桿體質(zhì)量范圍0.8～0.93kg。

圖4 推力桿二分之一模型Fig.4 Model of a half-thrust rod

對桿體樣品進行臺架實驗，在球鉸壓裝實驗中，單邊過盈量為0.15mm，質(zhì)量低于0.87kg的桿體，桿體末端都出現(xiàn)了裂縫，如圖5左所示，質(zhì)量大于0.91kg以上的桿體沒有出現(xiàn)裂縫；對壓裝合格的推力桿進行拉伸和壓縮實驗，實驗樣品數(shù)量各為5個，拉伸實驗中桿體能承受的拉伸載荷集中在60～80kN，如圖5中所示，在壓縮實驗中，試樣加載到-110kN出現(xiàn)了脆響，-150kN出現(xiàn)明顯裂縫，如圖5右所示。

圖5 破壞位置Fig.5 Damage location

2.3 推力桿的各向同性強度分析

采用傳統(tǒng)的 FEA 分析方法，不考慮纖維取向?qū)Ξa(chǎn)品性能的影響，將長玻纖增強PA66(GF50)復合材料視為各向同性材料，選取表2所示的彈性模量16.8GPa和泊松比0.33。橡膠球鉸壓入量為0.15mm，構(gòu)件能承受一定的拉伸和壓縮載荷，推力桿結(jié)構(gòu)模型如圖6所示。單元類型為C3D8I，單元數(shù)為52824，節(jié)點數(shù)為66077。

圖6 推力桿結(jié)構(gòu)模型Fig.6 Thrust rod structure model

通過剛性面來模擬球鉸壓入工序，點A，B為剛性面的幾何參考中心與桿體兩孔中心重合，剛性面與桿體設置接觸連接。分兩步加載對桿體系統(tǒng)進行數(shù)值模擬，首先是對A，B兩點進行固定約束，通過兩剛性面與桿體進行擠壓來模擬球鉸壓裝；然后解除點B在Z軸方向的約束，通過設置分析步，在B點的Z方向施加集中力F，F(xiàn)值的正和負，分別表示桿件拉、壓。模擬推力桿拉伸或壓縮，有限元分析結(jié)果顯示，三種狀況下顯示的最大主應力分別為96.2，107，132MPa，與實驗結(jié)果對比如圖7所示。從圖7可知，傳統(tǒng)FEA 各向同性仿真分析得到的球鉸壓裝和推力桿拉伸危險點位置與實驗結(jié)果吻合度不高。

圖7 危險點與破壞位置(a)球鉸壓裝0.15mm；(b)拉伸69kN；(c)壓縮120kNFig.7 Danger point and damage location(a)spherical hinge pressure 0.15mm;(b)stretch 69 kN;(c)compression 120 kN

2.4 推力桿的各向異性強度分析

2.4.1 分析流程

為了能夠準確分析預測玻纖增強復合材料注塑結(jié)構(gòu)件的承載能力和破壞位置，采樣各向異性方法進行強度分析的核心是建立玻纖增強復合材料注塑結(jié)構(gòu)件各向異性有限元離散模型。基于細觀力學原理，綜合運用復合材料細觀建模、注塑成型和結(jié)構(gòu)有限元分析技術(shù)，本工作提出了長玻纖增強復合材料注塑構(gòu)件各向異性分析方法，分析流程如圖8所示。

圖8 長玻纖增強復合材料注塑構(gòu)件各向異性分析流程Fig.8 Anisotropic analysis workflow of long glass fiber reinforced composite injection molding engineering structure

各向異性分析方法工作流程主要由4部分組成。具體如下：

①基于復合材料細觀力學的Eshelby夾雜理論，應用Mean Field均勻化方法，建立玻纖增強復合材料纖維均質(zhì)化RVE模型，并求出其正交各向異性彈性常數(shù)；

②基于廣義牛頓流體本構(gòu)方程，采用ARD-RSC纖維取向模型，考慮纖維間相互作用，應用數(shù)值模擬纖維增強復合材料注塑加工過程，預測注塑結(jié)構(gòu)件的纖維取向分布；

③對于注塑構(gòu)件整體而言，由于不同位置點的纖維取向不盡相同，導致材料各向異性彈性常數(shù)隨位置坐標的變化而變化。本工作將推力桿構(gòu)件劃分的每一個有限單元看作是一個RVE模型，構(gòu)建離散的RVE模型場，對每一個RVE模型沿纖維取向建立局部坐標系，通過局部坐標系和全局整體坐標之間的轉(zhuǎn)化關(guān)系，可獲構(gòu)件在整體坐標系下的材料各向異性彈性常數(shù)分布，進而建立構(gòu)件有限元3D離散均質(zhì)化RVE場網(wǎng)格模型。為了保證產(chǎn)品質(zhì)量，塑件的壁厚應盡量控制在10mm以下，壁厚方向的節(jié)點設置相對較密。在建立有限單元局部坐標系時，六面體網(wǎng)格的厚度方向為2方向，距離最小，熔體流動方向為3方向，垂直于厚度方向和熔體流動方向為1方向，具體局部坐標系如圖9所示；

圖9 局部坐標系設計Fig.9 Local coordinate system design

④基于結(jié)構(gòu)有限元分析技術(shù)，設置加載和約束方式，對構(gòu)件有限元3D離散均質(zhì)化RVE網(wǎng)格模型進行應力應變分析，獲取構(gòu)件的各向異性強度分析結(jié)果。

2.4.2 塑件模流分析及纖維取向分布

推力桿體采用長玻纖增強PA66(GF50)復合材料，玻纖初始長度為11mm，隨著玻纖初始長度的增加纖維間的相互作用就會增加，從而加劇注塑成型過程中纖維的斷裂。本工作基于廣義牛頓流體本構(gòu)方程，采用ARD-RSC纖維取向模型，建立成型工藝分析網(wǎng)格模型，設置實際工藝參數(shù)，數(shù)值模擬塑件注塑過程，預測注塑件的纖維取向分布，所得結(jié)果見圖10。

圖10 纖維取向分布Fig.10 Fiber orientation distribution

數(shù)值模擬采用Moldflow材料庫中的Verton RF-700-10-EM: SABIC Innovative Plastics US, LLC(50%長玻璃纖維增強型PA66)，推薦的熔體溫度范圍280～300℃。對推力桿實體模型先用Moldflow軟件對其劃分雙面層網(wǎng)格，接著對雙面層網(wǎng)格進行修復，然后將其轉(zhuǎn)為3D網(wǎng)格，網(wǎng)格尺寸為6mm，單元數(shù)為571667。工藝路線：注塑-保壓-冷卻，模具溫度90℃，熔體溫度為290℃，注塑時間為6s，保壓時間130s，保壓壓力110MPa，冷卻時間40s。圖10中的數(shù)值表示纖維取向張量主因子λ3，即纖維在熔體流動方向的集度，數(shù)值越趨近于1，纖維取向方向越趨于一致，也越能體現(xiàn)纖維在這一方向的承載能力。

2.4.3 材料細觀模型與FEA分析

基于復合材料細觀力學的Eshelby夾雜理論，應用Mean Field均勻化方法，建立長玻纖增強PA66(GF50)復合材料均質(zhì)化RVE模型。注塑成型后纖維長度主要分布在 2～4mm，粒子直徑為 15μm， 設置纖維長徑比為 220，纖維含量為50%(質(zhì)量分數(shù))。

對塑件進行有限元結(jié)構(gòu)分析的加載方式與2.3節(jié)相同。由表1可知，玻纖強度遠高于PA66基體的強度，因此玻纖增強復合材料注塑結(jié)構(gòu)件沿著纖維分布方向的增強效果最明顯，承載能力也最強。根據(jù)這一特性，根據(jù)注塑成型過程中形成的纖維取向分布，推力桿的桿頭部分纖維呈現(xiàn)同心圓取向分布，應用柱坐標下顯示的應力來判斷其是否破壞，桿體中間部分近似看成長方體，應用直角坐標下顯示的應力結(jié)果來判斷其是否破壞。有限元結(jié)構(gòu)分析結(jié)果與實驗結(jié)果對比見圖11所示。

圖11 危險點與破壞位置(a)球鉸壓裝0.15mm；(b)拉伸69kN；(c)壓縮120kNFig.11 Danger point and damage location(a)spherical hinge pressure 0.15mm;(b)stretch 69kN;(c)compression 120kN

最大主應力分別為64.7,103,146MPa。考慮到澆注口纖維分布較為混亂降低了該處的結(jié)構(gòu)強度，以及熔接縫的強度較低，對比桿體臺架實驗破壞位置圖，顯示仿真危險位置與實際破壞位置較為吻合，證明了各向異性分析方法的精準性與可行性。在推力桿壓縮實驗中，試樣加載到-110kN出現(xiàn)了脆響，-150kN出現(xiàn)明顯裂縫，說明在-110kN時，推力桿基體局部已發(fā)生破壞。在壓縮載荷為120kN時，各向異性分析與各向同性仿真分析應力結(jié)果分別為146，132MPa，表明在同樣載荷下，能表征纖維增強效果的各向異性分析能提前發(fā)現(xiàn)危險點，相對差值為10.6%。

2.5 推力桿結(jié)構(gòu)改進

由于桿體中部圓孔對桿體強度影響很大，同時考慮加強澆注口處與熔接縫強度，對推力桿結(jié)構(gòu)進行改進設計，新的桿體結(jié)構(gòu)如圖12所示。單元類型為C3D8I，單元數(shù)為32192，節(jié)點數(shù)為43305。參照2.4節(jié)各向異性分析方法，纖維取向分布如圖13所示，單元類型為3D，網(wǎng)格尺寸為6mm，單元數(shù)為700533。有限元分析結(jié)果見圖14。

圖12 推力桿新結(jié)構(gòu)模型Fig.12 New structure model of thrust rod

圖13 纖維取向分布Fig.13 Fiber orientation distribution

對推力桿結(jié)構(gòu)進行改進后，桿體中間部分強度有了顯著提升，結(jié)構(gòu)改進前后性能比對如表3所示。

3 結(jié)論

(1)基于廣義牛頓流體本構(gòu)方程，采用ARD-RSC纖維取向模型，考慮纖維間相互作用，仿真預測長玻纖增強復合材料注塑構(gòu)件的纖維取向分布；應用復合材料細觀力學Eshelby夾雜理論和Mean Field均勻化方法，建立長玻纖增強復合材料均質(zhì)化RVE模型。

(2)綜合運用復合材料細觀建模、離散RVE模型場、注塑成型和結(jié)構(gòu)有限元分析技術(shù)，提出了長玻纖增強復合材料注塑構(gòu)件各向異性強度分析方法。

圖14 危險點主應力(a)球鉸壓裝0.15mm；(b)拉伸69kN；(c)壓縮120kNFig.14 Dangerous point principal stress(a)spherical hinge pressure 0.15mm;(b)stretch 69kN;(c)compression 120kN

LoadPrincipalstress/MPaOriginalNewOptimizationeffect/%164．778．2+20．872103．044．1-57．183146．030．3-79．25

(3)對推力桿注塑構(gòu)件進行各向異性強度分析，顯示仿真危險位置與實際破壞位置較為吻合。在此基礎(chǔ)上對推力桿進行結(jié)構(gòu)改進，結(jié)果表明桿體中間部分在拉伸載荷下的最大主應力降低了57.18%，在壓縮載荷下的最大主應力降低了79.25%。

[1] 馬鳴圖，易紅亮，路洪洲，等.論汽車輕量化[J].中國工程科學，2009，11(9):20-27.

MA M T,YI H L,LU H Z,et al. On the light-weighting of automobile [J]. Engineering Sciences,2009,11(9):20-27.

[2] 馮美斌. 汽車輕量化技術(shù)中新材料的發(fā)展及應用[J].汽車工程，2006，28(3):213-220.

FENG M B. Development and applications of new materials in automotive light-weighting technologies [J]. Automotive Engineering,2006,28(3):213-220.

[3] 熊愛華，柳和生，黃益賓，等.長玻纖增強聚合物注塑充填行為和纖維取向模擬[J].高分子材料科學與工程， 2015，31(10):110-114.

XIONG A H,LIU H S,HUANG Y B,et al. 3D simulation of filling behavior and fiber orientation for long glass fiber reinforced polymer composites during injection molding process [J]. Polymeric Materials Science and Engineering, 2015,31(10):110-114.

[4] CHUN H H. Young’s modulus of unidirectional discontinuous-fibre composites [J]. Composites Science & Technology, 2000, 60(14):2671-2680.

[5] MONDALI M, ABEDIAN A, GHAVAMI A. A new analytical shear-lag based model for prediction of the steady state creep deformations of some short fiber composites [J]. Materials & Design, 2009, 30(4):1075-1084.

[6] SHARMA M, RAO I M, BIJWE J. Influence of fiber orientation on abrasive wear of unidirectionally reinforced carbon fiber-polyetherimide composites [J]. Tribology International, 2010, 43(5/6):959-964.

[7] 任超，陳建鈞，潘紅良. 隨機短纖維增強復合材料彈性模量預測模型[J].復合材料學報， 2012， 29(4):191-194.

REN C,CHEN J J,PAN H L. Prediction model for elastic modulus of random short fiber reinforced composite [J]. Acta Materiac Compositae Sinica,2012, 29(4):191-194.

[8] ESHELBY J D. The determination of the elastic filed of ellipsoidal inclusion and related problems [J]. Proceedings of the Royal Society London Series A, 1957,241(1226):367-396.

[9] HILL R. A self-consistent mechanics of composite materials. [J].Journal of the Mechanics & Physics of Solids,1965,13(4):213-222.

[10] MORI T,TANAKA K. Average stress in matrix and average elastic energy of materials with misfitting inclusions [J].Acta Metallurgica, 1973, 21(5):571-574.

[11] BRETHERTON F P. The motion of rigid particles in a shear flow at low reynolds number[J]. Journal of Fluid Mechanics, 1962, 14(2):284-304.

[12] FOLGAR F,TUCKER C L. Orientation behavior of fibers in concentrated suspensions [J]. Journal of Reinforced Plastics & Composites,1984,3(2): 98-119．

[13] WANG J,OGARA J F,TUCKER C L. An objective model for slow orientation kinetics in concentrated fiber suspensions: theory and rheological evidence [J]. Journal of Rheology,2008, 52(5):1179-1200.

[14] PHELPS J H, TUCKER C L. An anisotropic rotary diffusion model for fiber orientation in short-and long-fiber thermoplastics [J]. Journal of Non-Newtonian Fluid Mechanics,2009, 156(3):165-176.

[15] 沈觀林，胡更開，劉彬.復合材料力學[M].北京：清華大學出版社，2006：41-52.

SHEN G L, HU G K, LIU B. Mechanics of composite materials [M]. Beijing: Tsinghua University Press,2006: 41-52.

[16] BAY R S, TUCKER C L. Stereological measurement and error estimates for three-dimensional fiber orientation [J]. Polymer Engineering &Hyperlink Science, 1992, 32(32): 240-253.

[17] ADVANI S G, TUCKER C L. The use of tensors to describe and predict fiber orientation in short fiber composites [J].Journal of Rheology, 1987, 31(8):751-784.

[18] 江青松，柳和生，熊愛華，等. 長纖維增強聚合物注塑流動纖維取向分布數(shù)值模擬[J].高分子材料科學與工程，2015，31(3):123-127.

JIANG Q S,LIU H S,XIONG A H,et al. Numerical simulation on the fiber orientation distribution during flow in injection molding of long fiber reinforced polymer [J]. Polymeric Materials Science and Engineering, 2015,31(3):123-127.