基于離散元-有限元聯合算法的纖維流場變形行為

董關震,姜 茜,錢曉明,吳利偉

(1.天津工業大學 紡織科學與工程學院,天津 300387;2.天津工業大學 先進紡織復合材料教育部重點實驗室,天津 300387)

在紡織過程中,纖維的動態特性往往會影響織物結構性能,因此通過設計流場引導纖維運動是一種有效的質量控制手段。利用氣流可實現對纖維的加捻、運輸等過程。利用氣流控制纖維具有纖維磨損率小、纖維分布均勻、纖維排列可調控以及生產效率提升等優勢,其在諸多工藝環節如氣流紡紗[1-2]、氣流成網[3]、梳理成網[4-6]等均有應用。然而纖維在氣流力作用下伴隨著高度柔性變形與纖維間的碰撞,且流場分布差異對纖維運動過程與最終取向存在密切聯系,因此對纖維在流場作用下的變形轉移機制進行研究對于指導紡織生產至關重要。

利用計算流體動力學(CFD)獲得流場計算結果,再通過構建的纖維動力學方程模擬纖維在流場中的運動已成為一種重要研究手段。利用CFD可以解除環境條件的限制,且其具備模型改進便捷、時間成本低等優勢,以其針對紡織流場的研究正在不斷發展。由于二維模型始終無法考慮纖維在真實三維空間中的運動與變形,相關流場研究逐漸從二維流場擴展至三維[7-9]。基于紡織流場的特性,相關流場模型的應用也在逐漸發展,Jin等[10]采用大渦模擬研究了轉杯紡紗中轉杯氣流的渦特性,相比于k-e模型等傳統RANS模型,LES模型捕捉到了更多非定常非平衡大尺度效應與擬序結構,且與直接數值模擬相比減少了計算成本。不同纖維模型在氣流作用下會呈現不同形態,纖維模型主要分為2類,一類是剛性模型,如圓球[5,11]和圓柱模型[12],由于剛性模型無法改變其形狀,因此只能近似表征纖維聚集形態,沒有足夠的自由度來模擬纖維的變形;另一類是柔性模型[13-15],可以實現彎曲、拉伸和扭轉等情況。纖維與流體的耦合形式也會影響纖維的變形結果。在顆粒為稀相時使用流固單向耦合方法[16]計算纖維運動可以減少計算量,但隨著顆粒濃度的增加,顆粒運動對流場的反作用逐漸增強,此時使用雙向耦合[17-18]才能保證計算結果的準確性。目前大量紡織流場中的研究忽略了纖維柔性問題,使得對纖維動態轉移與彎曲變形的行為認識不足。此外柔性纖維多針對單纖維,而忽略了纖維流中多纖維的交互,不利于預測紡織品結構與性能。

本文擬構建三維流場與柔性多球鏈纖維模型的流固耦合動力學模型,基于歐拉-拉格朗日算法Euler-Lagrange建立流場與柔性多球鏈纖維模型間的耦合作用關系,揭示單根柔性纖維在氣流場中的瞬態受力變化與彎曲變形機制以及纖維簇在差異流場中的轉移演化差異,以實現高變形多數量纖維的動力學計算模型構建。

1 仿真建模

采取數值模擬中歐拉-拉格朗日算法Euler-Lagrange方法求解,將空氣視為連續相,顆粒視為離散相。耦合執行如下:①根據實際問題,可以依照時間步長、時間步數計算1個時間步內的流態,或直接進行穩態計算得到穩定流態;②通過耦合接口將流場作用在顆粒上的力傳遞給顆粒受力點;③通過引入的接觸模型與黏結模型,計算柔性纖維各點在流場作用下的運動;④更新顆粒位置,若為雙向耦合計算,則將纖維對流場的反作用反饋至流場,更新流場再次循環直至收斂,若為單向耦合計算,則纖維根據空間上更新的流場力進行迭代計算至收斂。⑤耦合計算結束,獲得纖維受力過程與形貌變形演化。整體流程如圖1所示。

圖1 CFD-DEM耦合過程

1.1 連續相模型

兩相耦合模型中,氣體設為連續相。假設氣體不可壓縮,控制方程為:

div(ρu)=0

(1)

(2)

式中:ρ為氣體的密度,kg/m3;u為氣流速度矢量;ui為u在坐標分量xi方向的分量(i=1,2,3),m/s;t為時間,s;μ為動力學黏度,N·s/m2;p為壓力,Pa;Si為xi方向的廣義源項。

采用穩態流場求解器,對流項采用二階迎風格式離散,使用基于壓力基的SIMPLE算法求解流場。

1.2 離散相模型

本文構建的纖維可被視為由n個球形顆粒與n-1根彈性桿組成的柔性多球鏈,柔性多球鏈纖維模型如圖2所示。在纖維鏈中,第a個顆粒只與第a-1與第a+1個顆粒相鄰,第1個和最后1個顆粒分別組成纖維頭端與尾端。

圖2 柔性多球鏈纖維模型

在初始階段,顆粒之間通過黏結鍵連接形成與纖維類似的柔性結構,粒子所受的力(Fn,t)扭矩(Tn,t)設為0,并在每一時間步遞增調整,如式(3)～(8)所示。

δFn=-vnSnAδt

(3)

δFt=-vtStAδt

(4)

δMn=-ωnStJδt

(5)

(6)

(7)

(8)

式中:RB是黏結的半徑,m;Sn、St分別為法向剛度和剪切剛度,N/m3;δt為時間步長,s;vn、vt為粒子的法向和切向速度,m/s;wnwt是法向和切向角速度,rad/s;Fn、Ft是法向和切向接觸力,N;Mn、Mt是法向和切向扭矩,N/m;A為黏結面積,m2;J為轉動慣量,kg·m2。

當法向σmax和切向τmax剪應力超過某個預定值時,黏結就會斷裂:

(9)

(10)

黏結形成后,顆粒在氣流作用下進行平移與旋轉,動力學控制方程如式(11)(12)所示。

(11)

(12)

球形顆粒與流體產生相互作用,顆粒間的碰撞使用軟球模型,彈性桿起到傳遞內力與連接形成纖維的作用,整體可以實現纖維的拉伸、彎曲與扭轉,較為完整地考慮了纖維的物理特性。

1.3 流場結構與纖維模型設置

本文在三維漸擴管道狀的流體域內進行計算模擬,并基于短纖維特征參數(見表1)設置離散元柔性鏈。流體域結構與特征參數見圖3。

表1 短纖維參數設計

圖3 流體域結構與特征參數

2 結果與討論

2.1 喂入纖維與氣流夾角對纖維變形的影響

為了比較喂入纖維與氣流夾角對纖維變形的影響,設置左端速度入口為2 m/s,得到了圖4所示的流道中心沿Y-Z剖面的流速分布與初始纖維位置。纖維在T=0時刻喂入,且與氣流方向的夾角分別為0°、10°、20°。

圖4 流道中心沿Y-Z剖面的流速分布與初始纖維位置

經過搭建的耦合計算平臺獲得了圖5所示的纖維在出口附近固定位置的最終取向變形。圖5(a)中纖維初始取向與中心流速方向平行,最終時刻纖維仍呈現伸直狀態,這意味著氣流只在沿纖維主體方向存在速度差時,纖維沒有形成彎鉤的趨勢。只有纖維初始取向與中心流速方向呈一定夾角時圖5(b)(c),纖維彎鉤才會形成,且隨著夾角的增加,最終形成彎鉤的形態越明顯。在角度為10°與20°情況下,指定的纖維頭端、中端尾端位置顆粒為典型位置點如圖6所示,由于Y方向為主要受力方向,因此對其進行Y方向的受力分析以進一步說明纖維形成彎鉤的過程。

圖5 纖維與氣流方向不同初始夾角時的最終彎曲情況

圖6 纖維典型位置點

圖7(a)(b)示出纖維在3個典型位置出現明顯受力差大小與時間的差異。當初始角度為10°時,纖維頭端、中端和尾端在0.015 s 開始出現了1.84×10-10N 的受力差異,而當初始角度為20°時則在0.001 s時出現了1.2×10-9N的受力差異,說明20°時纖維彎鉤效果最好,這是因為較早出現受力差異使纖維開始變形的時間提前,且較大的受力差異使纖維變形過程更為迅速。

圖7 纖維典型位置點的受力、速度與時間的關系

纖維端受力出現差異到產生變形差異也需要時間,力對變形速度和位移的影響效果仍需確認。如圖7(c)(d)所示,在纖維初始角度為20°時,典型位置點出現高速度差的時間更早(0.009 s<0.016 s),數值更大(0.30 m/s>0.15 m/s)。在尾端粒子與中端粒子的組成的A部分,由于尾端粒子處的氣體流速較慢,而中端流速較快,因此尾端與中端粒子在氣流力與沿纖維方向的拉力合力作用下,尾端粒子與中端粒子部分運動受限形成一個整體,且受力與速度在各時刻幾乎保持一致;而在頭端粒子與中端粒子的B部分,由于中端粒子處的氣流流速比頭端粒子處更快,因此頭端受氣流力方向與受纖維拉力方向一致,此時對應頭端速度曲線的降低階段。纖維為了達到受力穩定狀態,纖維頭端在合力作用下逐漸翻轉,該過程中端對頭端的拉力方向逐漸相反,當纖維中端拉力開始主導纖維頭端運動時,頭端速度曲線達到峰值。翻轉后纖維各位置受力與速度逐漸趨于一致,纖維頭端的受力與速度曲線向中端靠近,最終趨于穩定。速度曲線并未出現受力曲線因纖維內應力導致的震蕩,因此纖維因氣流作用導致的內應力差異對纖維變形的影響可以略不計。

通過圖8纖維典型位置隨時間在Y、Z方向位移的變化對纖維變形過程進行了探討。從Z方向上的位置變化可以看出,纖維頭端在Z方向上先上升后減小,這對應纖維頭端的先上揚后下落形成彎鉤的過程。隨著纖維初始喂入角度從10°增加至20°,上揚峰值從0.5 mm增加至1.0 mm,且上揚至峰值的時間減少,符合上述的受力與速度分析結果,纖維中端與尾端位移無明顯變化趨于平穩。從Y方向的位置變化可以看出,隨著時間推移,纖維頭端位置曲線逐漸沿Y方向靠近尾端,且由于纖維初始喂入角度的增加,頭端纖維靠近尾端的速度越快,20°時頭端向尾端靠近的程度最高,形成了近似中心彎鉤的效果。而在纖維中端與尾端,二者在Y軸方向上隨時間的位移變化始終保持一致,呈現較好的整體性。

圖8 纖維典型位置點的位移與時間的關系

2.2 喂入纖維位置對纖維變形的影響

纖維在流場中的變形程度取決于氣流力在纖維上的累積。上述研究已知,相同位置下纖維的喂入角度越大,氣流累積作用越強使其具有更好的彎曲變形效果。為深入探究更多影響氣流力累積的因素,進一步比較了在高纖維喂入角度下不同喂入位置對氣流的累計效果。在纖維初始喂入角度為20°時,將纖維分別放置于入口速度為2 m/s的差異流場計算域的左端入口位置 I與中心位置 Ⅱ 處(見圖3),圖9為不同纖維喂入位置下纖維在位置 Ⅲ 的最終取向變形結果。相比與從位置Ⅱ喂入的纖維,將纖維從位置I喂入流場中時在位置Ⅲ呈現出更強的彎鉤效應,說明纖維喂入位置也會顯著影響纖維變形。從圖10(a)(b)纖維頭端、中端與尾端典型位置點的受力與速度變化曲線可知,從位置I喂入的纖維由于初始氣流場較大的流速差異,在T=0.001 s時刻纖維頭端與中端的受力差到達峰值為1.2×10-9N,遠大于位置Ⅱ 喂入后在T=0.006 s時產生的2.57×10-10N峰值受力差,并導致纖維頭端與中端的流速峰值差異分別到達0.34 m/s與0.26 m/s(如圖10(c)(d)所示),最終使纖維上揚量從0.49 mm增加至1.00 mm。盡管流速差異會影響形態,但在此流場中不同初始取向、不同喂入位置的纖維各典型位置點會趨近于穩定速度0.7 m/s。

圖9 固定位置的纖維變形情況

圖10 纖維典型位置點的受力、速度隨時間變化規律

2.3 纖維交互作用對纖維變形的影響

為了研究單根纖維與纖維簇在流場運動過程中的差異,以10 m/s的入口速度條件對漸擴管流場特征進行計算。纖維簇模型包含同上的25根相同特征參數的單纖維。纖維之間的接觸力通過Hertz-Mindlin接觸模型計算,計算結果如圖11所示。單根纖維與纖維簇隨時間變化的運動位置基本吻合,但變形演化情況存在較大差異。T=0 s時,單根纖維與纖維簇被分別放入模型左端相同位置;T=0.002 s時,單根纖維整體發生屈曲變形,而纖維簇中只有上層纖維存在屈曲變形趨勢,上層纖維阻礙了下層纖維的運動與氣流力作用,導致下層纖維仍保持原有取向;T=0.004 s時,由于流場中心流速高,纖維中端逐漸引導下端伸直,使曲折部分逐漸轉移至纖維頭端,使頭端有上揚趨勢,在纖維簇中,下層纖維開始阻礙氣流運動導致上層纖維受氣流作用加劇,上層纖維變形大于單纖,纖維簇開始發生分層。T=0.006 s時刻,單根纖維頭端上揚形成彎鉤,纖維簇中9根上層纖維呈現大彎鉤效果,變形優于單根纖維,16根下層纖維呈現小變形,上揚趨勢不明顯。

圖11 單根纖維與纖維簇在流場作用下的運動演化過程

對周圍纖維對單根纖維的變形影響進行量化探究得到圖12在纖維相互作用下纖維典型位置點受力隨時間變化規律。從12(a)可以看出,只喂入1根纖維時,其在0.000 4 s時頭端與中端形成了5.01×10-9N的峰值受力差,而將這根纖維放在纖維簇中心再放入流場時,0.000 4 s時其頭端與中端只形成了4.44×10-9N的峰值受力差,相鄰纖維的交互作用使峰值受力差降低了11%。另外還比較了纖維簇層間受力的差異,從高變形量的上層區域與低變形量的下層區域分別選取2根典型纖維,從圖12(b)中得知,高變形量纖維在0.000 4 s時分別形成了7.3×10-9N與6.79×10-9N的峰值受力差,而低變形量纖維在分別在0.000 4 s時和0.000 8 s時形成了2.2×10-9與2.64×10-9N的峰值受力差,層間變形的差異是由于纖維交互作用與氣流分布差異的共同結果。

圖12 纖維相互作用下纖維典型位置點受力隨時間變化規律

3 結 論

通過構建柔性多球鏈纖維模型并結合CFD-DEM耦合模型進行三維仿真計算的研究方法,從纖維初始取向與氣流方向夾角、纖維喂入位置和纖維間相互作用方面對柔性纖維與纖維簇在差異流場中的彎曲運動進行了研究,得出以下結論。

①當纖維與氣流對齊對稱地放置在通道中心線時,纖維沿水平路徑直線前進,纖維的傾斜打破了對稱性,導致纖維傳輸路徑的偏差和波動,初始纖維的傾斜程度越高時,纖維越早發生屈曲。

②纖維的變形程度與氣流力在纖維上的累積密切相關,增加纖維在流場中的運動距離有利于氣流力的累積,纖維到達相同位置處的變形加劇,直至各纖維位置點趨于相同速度。

③纖維簇中纖維在流速差和層間相互作用下呈現分層現象,上層少數纖維會比單根纖維更快彎曲,下層多數纖維則更慢,纖維層間作用使峰值受力差降低了11%,纖維簇阻礙了單根纖維變形,不利于纖維的均勻混雜。