基于非牛頓流體本構(gòu)方程的熔噴纖維直徑預(yù)測(cè)

孫光武， 李杰聰， 辛三法， 王新厚

(1. 上海工程技術(shù)大學(xué) 服裝學(xué)院， 上海 201620； 2. 東華大學(xué) 紡織學(xué)院， 上海 201620)

熔噴技術(shù)是由高溫高速氣流牽伸熔融聚合物進(jìn)而一步成纖的紡絲技術(shù)，制備出的纖維(一般為1～5 μm)在成網(wǎng)簾上逐漸聚集而形成熔噴非織造布。由于纖維較細(xì)，具有較大的比表面積，熔噴非織造布往往被用于制造空氣和液體過(guò)濾材料、吸油材料等。熔噴非織造布在應(yīng)用中的表現(xiàn)直接受到其纖維直徑的影響，具有較細(xì)纖維的熔噴非織造布的過(guò)濾性能往往較好。熔噴非織造生產(chǎn)線的工藝參數(shù)對(duì)熔噴非織造纖維直徑具有很大的影響，需要反復(fù)調(diào)整才可獲得最佳的纖維細(xì)度，因此，熔噴拉伸理論模型的研究特別是涉及到纖維直徑的預(yù)測(cè)則顯得愈發(fā)重要。因?yàn)槠淠茉诶碚撋项A(yù)測(cè)直徑，從而指導(dǎo)實(shí)際生產(chǎn)，并能從機(jī)制上闡述各工藝參數(shù)之間的關(guān)系，為商業(yè)化生產(chǎn)提供理論依據(jù)。

早期的熔噴拉伸理論由Uyttendaele等[1]通過(guò)改進(jìn)熔融紡絲模型得到，這種模型被稱(chēng)為“自由端”纖維-氣流模型。因在熔噴技術(shù)中氣流牽伸纖維，纖維的尾端并無(wú)握持力，這與熔融紡絲“非自由端”不同。文獻(xiàn)[2-3]將熔噴拉伸模型分別推廣至二維和三維空間。各類(lèi)方法如蒙特卡洛法、龍格-庫(kù)塔法、有限元法和冪律模型[4]也被廣泛用于求解熔噴拉伸模型方程。Sinha等[5]提出一種簡(jiǎn)單的準(zhǔn)一維模型，該模型計(jì)算量較小，同時(shí)也能精確預(yù)測(cè)纖維的運(yùn)動(dòng)和直徑變化。Shambaugh等[6]以靜態(tài)結(jié)晶動(dòng)力學(xué)方程為基礎(chǔ)，模擬了熔噴纖維在成纖過(guò)程中的結(jié)晶變化，該方法成功將熔噴纖維的宏觀結(jié)構(gòu)預(yù)測(cè)與微觀結(jié)構(gòu)預(yù)測(cè)相結(jié)合。不同于上述模型，Sun等[7-8]提出了一種基于拉格朗日方法的拉伸模型，通過(guò)將纖維假設(shè)成由珠子、彈簧、黏壺串聯(lián)的模型，預(yù)測(cè)每個(gè)珠子的位置和物理性能的變化，進(jìn)而預(yù)測(cè)纖維的實(shí)時(shí)位置和各項(xiàng)性能參數(shù)。

熔噴纖維拉伸模型經(jīng)過(guò)多年的發(fā)展，從簡(jiǎn)單到復(fù)雜，從一維到三維，目前文獻(xiàn)中已報(bào)道的模型可預(yù)測(cè)任意時(shí)刻纖維在整個(gè)氣流場(chǎng)中的位置、溫度、速度、直徑，甚至結(jié)晶度。然而為了簡(jiǎn)化計(jì)算，在模型中采用牛頓流體本構(gòu)方程或者經(jīng)典的Maxwell本構(gòu)方程，這些方程由于難以準(zhǔn)確描述聚合物熔體復(fù)雜的非牛頓黏彈行為，導(dǎo)致模擬的纖維直徑與實(shí)際結(jié)果差異較大。為此，本文引入4種非牛頓流體本構(gòu)方程：PTT(Phan-Thien-Tanner)方程、UCM(Upper-Convected Maxwell)方程、Giesekus方程和Rouse-Zimm方程，結(jié)合已建立的拉格朗日熔噴拉伸模型[9-11]預(yù)測(cè)纖維直徑，以定量化分析非牛頓流體本構(gòu)方程對(duì)最終纖維直徑預(yù)測(cè)結(jié)果的影響[12]，從機(jī)制上揭示非牛頓流體的黏彈行為對(duì)熔噴纖維的成纖作用。

1 理論模型

拉格朗日方法是一種描述流體運(yùn)動(dòng)的基本方法，即對(duì)流體質(zhì)點(diǎn)進(jìn)行標(biāo)記，流體質(zhì)點(diǎn)的空間坐標(biāo)與物理量就表示為質(zhì)點(diǎn)及時(shí)間的函數(shù)。建立熔噴纖維拉伸模型的目的是為了研究熔噴纖維在拉伸過(guò)程中的細(xì)化規(guī)律、纖維內(nèi)外應(yīng)力變化、纖維在氣流場(chǎng)中的運(yùn)動(dòng)情況以及纖維段之間的相互作用。模型的研究對(duì)象是纖維片段，其目的是求解出纖維片段的各物理量與空間坐標(biāo)在流場(chǎng)中隨時(shí)間變化的規(guī)律，而這與拉格朗日方法一致。

圖1示出采用珠鏈模型描述的熔噴纖維。纖維可假設(shè)為由珠子和串聯(lián)珠子的熔體單元(類(lèi)似鏈條)組成。

圖1 纖維珠鏈模型Fig.1 Bead-chain model of fiber

在不考慮纖維相互糾纏的前提下，熔體從噴絲孔逐漸細(xì)化落至成網(wǎng)簾上，纖維所受到的作用力主要有4種：氣流拉伸力、黏彈力、表面張力和重力。其中：氣流拉伸力和重力為外應(yīng)力，主要起推動(dòng)纖維運(yùn)動(dòng)、拉伸并細(xì)化纖維的作用；而表面張力和黏彈力屬內(nèi)應(yīng)力，主要起抵抗纖維拉伸和彎曲形變的作用。由于整體的動(dòng)力學(xué)模型較為復(fù)雜，且除黏彈力之外的其他作用力并非本文的研究重點(diǎn)，在此僅列出數(shù)理方程，詳細(xì)推導(dǎo)過(guò)程可參見(jiàn)文獻(xiàn)[9-11]。

氣流拉伸力是纖維在高速流場(chǎng)中受到的主要作用力，氣流力的大小直接影響纖維在氣流中的運(yùn)動(dòng)和最終成纖時(shí)的纖維直徑。纖維段(i-1,i)上所受到的氣流力Fd,i-1,i可以分解成壓差阻力Fp,i-1,i和摩擦阻力Ff,i-1,i，有

Fd,i-1,i=Ff,i-1,i+Fp,i-1,i

(1)

(2)

(3)

式中：Cp、Cf分別為壓差阻力系數(shù)、摩擦阻力系數(shù)；ρa(bǔ)為空氣密度，其值為1.29 kg/m3；vrt、vrn為纖維相對(duì)速度(纖維實(shí)際速度與所在位置氣流速度的差)沿其軸向和徑向的分量，m/s；di-1,i為纖維段(i-1,i)的直徑，m；li-1，i為纖維段(i-1,i)的長(zhǎng)度，也是珠子i-1與i的間距，m；ft、fn分別為纖維軸向單位向量、徑向單位向量。

表面張力與纖維的氣流力方向相反，纖維段(i-1,i)的表面張力Fb,i-1,i由下式確定：

(4)

式中：(xi,yi)、(xi-1,yi-1)分別為珠子i、i-1的坐標(biāo)；b、j為x和y方向的單位向量；θ為纖維表面張力系數(shù)，本文取0.7 kg/s2 [1-3]；ki-1,i為纖維段(i-1,i)的曲率。

黏彈力是表征聚合物非牛頓行為的主要作用力，也是本文的重點(diǎn)研究對(duì)象。纖維段(i-1,i)的黏彈力Fve,i-1,i[12]可表示為

(5)

式中：ri與ri-1分別為珠子i與i-1的空間位置；N2,i-1,i為纖維段(i-1,i)所受第二法向應(yīng)力。

根據(jù)張量分析原理，在理想拉伸狀態(tài)下，對(duì)于某一纖維段，有N2=τt-τn。其中τt和τn分別為軸向偏應(yīng)力張量和徑向的偏應(yīng)力張量。在本文研究中，引入4種非牛頓流體本構(gòu)方程分別計(jì)算出4種不同的偏應(yīng)力張量，進(jìn)而計(jì)算出不同的黏彈力，再結(jié)合熔噴纖維拉伸模型預(yù)測(cè)纖維直徑。這4種非牛頓流體本構(gòu)方程[13]分別為：

UCM方程

(6)

PTT方程

(7)

Giesekus方程

(8)

Rouse-Zimm方程

(9)

張量分解公式

(10)

(11)

式中：Tf,i-1,i為纖維段(i-1,i)的溫度，K。在實(shí)際計(jì)算中，利用有限差分方法計(jì)算上述非牛頓流體本構(gòu)方程，得到4個(gè)不同的第二法向應(yīng)力值，并分別代入式(5)計(jì)算黏彈力。

為使整個(gè)模型封閉，還需纖維的運(yùn)動(dòng)控制方程和能量守恒方程，分別為：

(12)

(13)

式中：g為重力加速度，取值為9.8 m/s2；Cb為黑體輻射系數(shù)，取值為5.67 W/(m2·K4)[10]；hi-1,i為纖維段(i-1,i)與周?chē)鷼饬髦g的傳熱系數(shù)，取值為0.026 W/(m2·K·s)[10]；ψ為纖維黑度，取值為0.4[10]；Ta,i-1,i為纖維段(i-1,i)周?chē)臍饬鳒囟?，K；mi-1,i為纖維段(i-1,i)的質(zhì)量，kg；Si-1,i為參與熱對(duì)流的纖維面積，m2；Cg為纖維段(i-1,i)的比熱容，J/(kg·K)，可依據(jù)下式[14]進(jìn)行計(jì)算：

Si-1,i=πdi-1,ili-1,i

(14)

(15)

Cg,i-1,i=366.9+2.42Tf,i-1,i

(16)

式中，ρf,i-1,i為纖維段(i-1,i)的密度，kg/m3。以聚丙烯為例，其密度與溫度的經(jīng)驗(yàn)方程[14]為

(17)

2 算法描述

熔噴纖維拉伸模型在第1節(jié)已經(jīng)建立，在求解時(shí)仍缺少計(jì)算氣流力所必需的流場(chǎng)速度和溫度，為此需采用Fluent有限元分析軟件模擬熔噴氣流場(chǎng)(設(shè)置流場(chǎng)入口速度為60 m/s，溫度為640 K，模頭溫度為640 K，出口為標(biāo)準(zhǔn)大氣環(huán)境)，并導(dǎo)出氣流場(chǎng)三維坐標(biāo)、速度和溫度的數(shù)值解[7-9]。模型的求解方法可概括為：首先依據(jù)纖維段的初始物理參數(shù)值，并結(jié)合式(13)與各作用力的數(shù)理方程分別計(jì)算出氣流力、表面張力和黏彈力，再依據(jù)式(12)求解纖維段的速度與位移，最后依據(jù)式(15)求解纖維段的直徑。由于采用了4種不同的非牛頓流體本構(gòu)方程，黏彈力需分別求解4次，以觀察直徑結(jié)果的差異。理論模型的求解過(guò)程如下。

步驟1：設(shè)置i=2，并賦予纖維段(i-1,i)初始長(zhǎng)度為10 mm，初始直徑為400 μm，初始溫度為640 K。依據(jù)式(11)與式(15) ～ (17)計(jì)算纖維段的初始黏度、質(zhì)量、比熱容、密度。導(dǎo)入Fluent計(jì)算的流場(chǎng)坐標(biāo)、速度與溫度數(shù)值解。

步驟2：計(jì)算開(kāi)始。讀取珠子i-1和i的空間坐標(biāo)及其周?chē)牧鲌?chǎng)溫度和速度。根據(jù)式(1) ～ (3)計(jì)算纖維段(i-1,i)受到的氣流力；根據(jù)式(4)計(jì)算其表面張力；從式(6) ～ (9)中選擇一種非牛頓流體本構(gòu)方程，并根據(jù)式(5)計(jì)算纖維段(i-1,i)的黏彈力；依據(jù)式(13)計(jì)算纖維段(i-1,i)的溫度；依據(jù)式(15)計(jì)算纖維段的直徑；依據(jù)式(12)中纖維段受到的合力的大小和方向，計(jì)算珠子i-1和i在Δt時(shí)間內(nèi)(可自行設(shè)置Δt值，數(shù)值越小則輸出數(shù)據(jù)越多，結(jié)果越精確，程序運(yùn)行越慢；反之則結(jié)果精確性較低，程序運(yùn)行較快)的位移。輸出各作用力的值、珠子空間坐標(biāo)和纖維段直徑的計(jì)算結(jié)果。

步驟3：根據(jù)步驟2計(jì)算的位移，更新珠子i-1和i到新的位置。

步驟4：判斷珠子i與噴絲孔(原點(diǎn))的距離是否達(dá)到給定值(該值可自由設(shè)定，越大則纖維越長(zhǎng)，計(jì)算耗時(shí)越久)。如果是，則執(zhí)行步驟5；反之，返回至步驟2，以計(jì)算珠子i-1和i在下一個(gè)Δt后的位置。

步驟5：Uyttendaele等[1]提出內(nèi)外應(yīng)力一致的“凝固點(diǎn)”用以描述纖維直徑停止變化的邊界條件。在模型中，本文也采用同樣的方法，即判斷纖維段受到的氣流力與重力之和是否等于黏彈力與表面張力之和。如果是，停止計(jì)算纖維段的直徑變化；反之，在噴絲孔處加入珠子i+1。賦予系統(tǒng)中新纖維段(i,i+1)相應(yīng)的初始物理屬性。返回步驟2。

由于需要研究4種不同的非牛頓流體本構(gòu)方程對(duì)計(jì)算結(jié)果的影響，所以模型需進(jìn)行多次計(jì)算，每次計(jì)算采用一種不同的本構(gòu)方程。在綜合考慮計(jì)算時(shí)間與計(jì)算機(jī)性能的基礎(chǔ)上，本文采用400個(gè)珠子模擬纖維在流場(chǎng)中的運(yùn)動(dòng)和細(xì)化程度，程序采用MatLab編程求解。當(dāng)?shù)?00個(gè)珠子出現(xiàn)在噴絲孔處時(shí)，所有運(yùn)算停止。

本文研究采用聚丙烯的物理屬性作為原料參數(shù)，但該模型可模擬幾乎所有的熔噴原料及工藝過(guò)程，僅需進(jìn)行相應(yīng)的修改。如更換其他材料，需獲得該材料的黏度、比熱容和密度方程，即修改式(11)、(16)和(17)；還需修改物理參數(shù)值至合適的值，包括應(yīng)力松弛時(shí)間、傳熱系數(shù)和表面張力系數(shù)。若更改其他熔噴工藝條件，可在Fluent模擬氣流場(chǎng)時(shí)，設(shè)置不同的壓縮氣流速度、氣流溫度、模頭溫度和外部氣流環(huán)境等；同時(shí)需依據(jù)噴絲孔直徑設(shè)置纖維的初始直徑，依據(jù)模頭溫度設(shè)置纖維的初始溫度。

3 結(jié)果與討論

3.1 模擬熔噴過(guò)程中纖維應(yīng)力的變化

圖2示出纖維與噴絲孔距離z對(duì)模擬獲得的纖維內(nèi)外應(yīng)力的影響。

圖2 不同非牛頓流體本構(gòu)方程得到的纖維應(yīng)力隨z坐標(biāo)的變化Fig.2 Fiber stress calculated by different non-Newtonian fluid constitutive equations and its variation with z coordinate

本文采用4種非牛頓流體本構(gòu)方程計(jì)算得到4種不同的內(nèi)應(yīng)力。由圖可知，所有曲線的應(yīng)力值都是在出現(xiàn)1個(gè)峰值后，隨z值的增大而逐漸下降。即隨著纖維遠(yuǎn)離噴絲孔，應(yīng)力值逐漸下降至穩(wěn)定。在纖維應(yīng)力值達(dá)到最大峰值之前，纖維距離噴絲孔較近，由于噴射的2股高速射流在該區(qū)域相互融合形成湍流，導(dǎo)致該區(qū)域內(nèi)纖維受到的各項(xiàng)應(yīng)力不穩(wěn)定，此時(shí)的曲線波動(dòng)性較大。同時(shí)由于此處的高速射流速度最大，所以產(chǎn)生了最大的氣流拉伸應(yīng)力，纖維受到這種強(qiáng)烈外應(yīng)力作用需要抵抗拉伸力，故而內(nèi)應(yīng)力也達(dá)到峰值。隨著纖維遠(yuǎn)離噴絲孔，纖維周?chē)鷼饬魉俣戎饾u衰減，其受到的外應(yīng)力逐漸下降，內(nèi)應(yīng)力也隨之下降，然而由于熔融纖維的非牛頓流體黏彈特性，黏彈力會(huì)有一定的遲滯，這產(chǎn)生了內(nèi)外應(yīng)力的應(yīng)力差，在圖2中體現(xiàn)為外應(yīng)力曲線與4種內(nèi)應(yīng)力曲線在相同z值處的間距。在纖維距離噴絲孔的間距小于約3 cm的區(qū)域內(nèi)，內(nèi)外應(yīng)力差較大，這說(shuō)明該區(qū)域內(nèi)纖維細(xì)化程度較大。在間距超過(guò)3 cm后，內(nèi)外應(yīng)力差逐漸縮小，這說(shuō)明纖維細(xì)化程度逐漸減小。當(dāng)間距超過(guò)5 cm后，內(nèi)外應(yīng)力逐漸趨于穩(wěn)定，這意味著纖維直徑也逐漸停止變化。

3.2 纖維直徑

為與Uyttendaele等的在線測(cè)量實(shí)驗(yàn)進(jìn)行對(duì)比以展示模擬結(jié)果的準(zhǔn)確程度，在模型中采用了與該在線測(cè)量實(shí)驗(yàn)相同的工藝條件(具體工藝條件已在第2節(jié)說(shuō)明)，實(shí)驗(yàn)與模擬的纖維直徑變化示于圖3中。為進(jìn)一步闡明圖3中纖維直徑的變化機(jī)制，表1示出采用4種非牛頓流體本構(gòu)方程計(jì)算獲得的指標(biāo)結(jié)果。纖維最大內(nèi)外應(yīng)力差是指纖維受到的內(nèi)外應(yīng)力差的最大值，是圖2中外應(yīng)力曲線分別與4種內(nèi)應(yīng)力曲線之間的最大差值；纖維平均細(xì)化速度是圖3中纖維直徑變化的平均斜率；凝固點(diǎn)位置是指纖維固化時(shí)與噴絲孔的間距，在圖3中表示為曲線停止變化時(shí)對(duì)應(yīng)的z值；最終纖維直徑是指纖維經(jīng)過(guò)牽伸后的直徑，在圖3中表示為曲線尾端z為10 cm處對(duì)應(yīng)的纖維直徑。

圖3 不同非牛頓流體本構(gòu)方程得到的纖維直徑隨z坐標(biāo)的變化Fig.3 Fiber diameter calculated by different non-Newtonian fluid constitutive equations and its variation with z coordinate

非牛頓流體方程名稱(chēng)最大內(nèi)外應(yīng)力差/kPa平均細(xì)化速度/(μm·cm-1)凝固點(diǎn)位置/cm最終纖維直徑/μmGiesekus18.72128.006.4574.12PTT12.54126.865.2779.18UCM14.36127.025.3683.42Rouse-Zimm9.21125.807.6274.76

由圖3可知，Uyttendaele實(shí)驗(yàn)結(jié)果[1]與模擬結(jié)果變化趨勢(shì)較為相符。纖維自噴絲孔擠出后，其直徑變化并不穩(wěn)定，這是由于噴絲孔附近的氣流場(chǎng)形成的湍流所致。隨后纖維直徑快速下降，這是由于纖維較大的內(nèi)外應(yīng)力差所致。從表1可明顯發(fā)現(xiàn)，纖維細(xì)化速度(即曲線在z小于3 cm區(qū)域的平均斜率)與內(nèi)外應(yīng)力差高度相關(guān)，最大內(nèi)外應(yīng)力差從大到小排序?yàn)椋篏iesekus > UCM > PTT >Rouse-Zimm，而纖維平均細(xì)化速度也有同樣的順序。當(dāng)纖維與噴絲孔間距z大于3～4 cm時(shí)，由于內(nèi)外應(yīng)力差縮小，纖維細(xì)化速度減緩。間距z超過(guò)5～6 cm后，內(nèi)外應(yīng)力差穩(wěn)定導(dǎo)致纖維直徑也逐漸趨于穩(wěn)定。這期間，當(dāng)纖維達(dá)到凝固點(diǎn)時(shí)，纖維直徑停止變化。采用Rouse-Zimm方程計(jì)算得到纖維達(dá)到凝固點(diǎn)時(shí)距離噴絲孔最遠(yuǎn)為7.62 cm，而采用PTT方程計(jì)算得到纖維達(dá)到凝固點(diǎn)時(shí)距離噴絲孔最近為5.27 cm。凝固點(diǎn)出現(xiàn)的位置較遠(yuǎn)可使纖維有充分的空間拉伸，更容易產(chǎn)生較細(xì)的纖維。

最終預(yù)測(cè)的纖維直徑從大到小排序?yàn)椋篣CM > PTT > Rouse-Zimm > Giesekus。采用Giesekus方程計(jì)算得到的纖維直徑更符合Uyttendaele實(shí)驗(yàn)結(jié)果。所有預(yù)測(cè)的最終纖維直徑均大于實(shí)驗(yàn)測(cè)量的結(jié)果(約60 μm)。這是由于模型中僅考慮了單個(gè)纖維段的受力與形變而沒(méi)有考慮纖維段之間的相互纏結(jié)影響。在Bresee等[17-18]的實(shí)驗(yàn)觀察中，纖維之間明顯的接觸和纏繞現(xiàn)象會(huì)導(dǎo)致纖維直徑減小。

4 結(jié) 論

本文以熔噴纖維三維拉伸模型為基礎(chǔ)，采用拉格朗日方法對(duì)纖維在空氣中的動(dòng)力學(xué)進(jìn)行分析，同時(shí)結(jié)合不同的非牛頓流體本構(gòu)方程，預(yù)測(cè)纖維直徑的變化情況。雖然模擬結(jié)果與實(shí)驗(yàn)結(jié)果基本相符，但是不同的非牛頓流體本構(gòu)方程的應(yīng)用使計(jì)算出的黏彈力有差異，從而導(dǎo)致纖維在拉伸形變時(shí)內(nèi)外應(yīng)力差異的變化。纖維直徑的模擬結(jié)果顯示，纖維的最終直徑不僅受到纖維在氣流場(chǎng)中的內(nèi)外應(yīng)力差的影響，還受到凝固點(diǎn)位置的影響，即：纖維的內(nèi)外應(yīng)力差越大，纖維細(xì)化速度越快；纖維的凝固點(diǎn)距離噴絲孔越遠(yuǎn)，則纖維有更加充分的空間進(jìn)行拉伸，其最終直徑也越細(xì)。本文的數(shù)學(xué)模型研究揭示了熔噴纖維的細(xì)化機(jī)制，同時(shí)非牛頓流體本構(gòu)方程的引入，也成功定量對(duì)比分析了熔噴纖維黏彈力學(xué)行為的非牛頓性，為熔噴纖維拉伸模型非牛頓方程的選擇提供了計(jì)算依據(jù)。

FZXB