黏彈性表面活性劑溶液中顆粒沉降特性研究

蘇曉輝，張弛，徐志鋒，金輝，王治國

（1 西安石油大學機械工程學院，陜西西安 710065； 2 西安交通大學動力工程及多相流國家重點實驗室，陜西西安 710049）

引 言

黏彈性流體懸浮顆粒流廣泛存在于自然界和工業生產制造過程，如生物體液中細胞運動[1-2]、聚合物和日化用品加工、油氣藏水力壓裂[3-4]、微流控裝置中顆粒的連續操控與分離[5-8]。表面活性劑溶液達到一定濃度時表現出很強的黏彈性，具有獨特的流變學性質，在工業生產和制造領域受到廣泛關注。研究人員[9-10]發現多元聚合物-表面活性劑中表面活性劑的用量對于原油乳狀液的穩定性具有重要影響，溶液彈性能夠有效提高驅油效率。此外，表面活性劑溶液具有湍流減阻、強化傳熱作用，在流體管道輸送和換熱領域也得到了廣泛應用[11-12]。由于其強非線性應力應變率特性和應力松弛效應，表現出十分復雜的流動現象[11]。當黏彈性表面活性劑溶液中有顆粒存在時，體系的流動特性更為復雜，例如顆粒在其中沉降時顆粒速度的“過沖”現象[13]，顆粒沉降速度的持續波動[14-16]等。目前對于具有復雜流變性質的黏彈性表面活性劑溶液中顆粒沉降規律仍缺乏清晰的認識，顆粒所呈現的復雜沉降現象背后的機理尚不明確。

黏彈性流體中顆粒沉降特性取決于流體的流變性質以及操作條件。國內外學者[17-23]針對黏彈性流體中顆粒運動特性開展了相關研究，主要集中于剪切條件下流體中顆粒橫向運移以及成鏈特性，探索了顆粒宏觀運動與流體流變性質之間的關系。顆粒群在流體中運動主要取決于流體作用于顆粒的曳力，為此學者們[24-25]對黏彈性溶液中顆粒曳力進行了研究，發現彈性(用Deborah 數De=λu/a來表征，其中λ為松弛時間；u為流體特征速度；a為特征尺寸，取顆粒半徑)對顆粒曳力具有顯著影響，在常黏度的Boger 流體中，曳力隨著De增加先減小至某一最小值后又增加，而在剪切稀化流體中，曳力隨著De單調減小，并且減小速度和幅度與剪切稀化的強度相關。基于此，提出了黏彈性流體中顆粒曳力理論公式，用于描述兩相流動相間動量交換[24,26]。此外，對黏彈性流體中顆粒沉降特性開展了研究[2,27]。有學者[13,28-30]采用虛擬區域法研究黏彈性流體中單顆粒和多顆粒的沉降特性，發現流體彈性和壁面效應均會影響顆粒沉降，他們均采用只考慮黏彈性的Oldroyd-B 本構模型。Oldroyd-B 本構模型由于可以簡潔明確地描述流體的黏彈性，在涉及描述表面活性劑溶液體系的研究中得到了廣泛應用[10,13]。

然而，實際工程中的黏彈性表面活性劑溶液，不僅具有時間松弛效應，其應力應變率響應也具有非線性特性，如剪切稀化、拉伸硬化和剪切帶等，其對顆粒沉降也有顯著影響。學者們[14-15,31]通過實驗對一定濃度下的表面活性劑溶于水后形成的蠕蟲狀膠束溶液中顆粒沉降進行研究，發現顆粒在這種黏彈性溶液中沉降速度呈現隨機波動，不能達到穩定的終端速度，并且認為這一現象與顆粒尾跡區的復雜拉伸應變特性相關。Wang 等[16]也在實驗中發現了這一效應，并且發現在某些條件下，不僅顆粒速度在沉降過程中隨機波動，顆粒沉降方向也發生了變化，即顆粒在流體中出現了沉降-反彈的彈跳現象，認為這一現象與此類蠕蟲狀膠束流體的流變性質密切相關[3-4,32]，推測顆粒沉降過程中下游的“負尾跡”可能是引起顆粒沉降速度波動的原因，因此有必要對其機理展開研究。

本文采用同時考慮表面活性劑溶液黏彈性和剪切稀化及拉伸硬化效應的FENE-P 模型以及考慮了大分子之間相互作用的用于濃表面活性劑溶液的Giesekus模型，研究顆粒在其中的沉降特性，討論溶液的流變性質對顆粒沉降特性的影響規律以及顆粒沉降過程中尾跡區的特征，研究結果有助于進一步認識和解釋黏彈性表面活性劑溶液中顆粒沉降速度的波動機理。

1 數學模型

考慮顆粒沉降運動過程，為了簡化模擬，選取顆粒為參考坐標系，顆粒在沉降過程中受流體曳力以及重力和浮力作用，其運動由牛頓運動定律描述。同時通過非慣性系的體積力將其與描述流體流動的偏微分方程耦合。對于瞬態、不可壓縮黏彈性表面活性劑溶液的絕熱流動，其控制方程包括連續方程、動量守恒方程以及描述附加彈性應力的黏彈性本構方程。

連續方程

動量方程

其中，ρ為流體密度；u為流體速度；p為流體壓力；μS為溶劑黏度；τP為附加黏彈性應力，由黏彈性流體的本構方程確定；F為流體在非慣性系中運動時的慣性力。對于小球參考系，慣性力F只有垂直分量Fz

黏彈性本構方程

其中，λ為松弛時間；μP為溶質黏度；fr為松弛函數；fp為黏度系數。對于不同的本構方程，fr和fp具有不同的形式。

上隨體導數

Oldroyd-B本構方程

Oldroyd-B 模型把大分子溶液看作是牛頓溶劑流體中懸浮的胡克彈簧體系構成，由分子動力學理論推導而來，可用來預測基本的黏彈性特征，但是其假定黏度為常數，并且可以無限拉伸，因而在拉伸流動中常常得到不符合物理實際的結果，是比較簡單的一種黏彈性模型。

FENE-P本構方程

其中，L為伸長率。

FENE（finitely extensible nonlinear elastic）模型也是基于分子動力學理論，但是將大分子看作是非線性彈簧，并且考慮了大分子的有限拉伸率。基于FENE模型的Peterlin封閉模型同時考慮了大分子有限拉伸以及流體剪切稀化效應。

Giesekus本構方程

其中，α為遷移因子。

Giesekus 本構模型用于描述半稀和濃溶液，考慮了大分子之間相互作用產生的流動阻力，因此方程中增加了二次非線性項。

顆粒運動由牛頓運動定律確定

其中，m為顆粒質量；z為顆粒位移；Fz為流體對顆粒的作用力；Fg為顆粒重力。

流體對顆粒曳力為

其中，I是單位張量；n為球表面的單位法向量；ez為流體流動方向單位矢量；?Ωs是球形顆粒的外表面積。

顆粒所受重力為

其中，ρs為顆粒密度；Vs為顆粒體積；a為顆粒半徑。

雷諾數和維森貝格數定義為

其中，d為顆粒粒徑；零剪切黏度μ0=μS+μP；特征速度Us取小球的瞬時沉降速度。

2 數值模擬方法

模擬流動區域和計算區域如圖1所示，半徑為a的球形顆粒位于半徑為R的圓柱體容器的對稱軸上，顆粒在重力、浮力和曳力作用下從靜止開始沉降，豎直方向上瞬時沉降速度為Us，初始時刻t=0 s時沉降速度Us=0。由于取運動顆粒自身為參考系，那么對于計算域的流體運動只需考慮流體相對于顆粒的運動即可，顆粒相對于流體具有平均速度Us，因此在遠離顆粒的計算域入口處流體相對于顆粒參考系的速度u=Us，即流體從計算域下方入口?Ωin以顆粒瞬時沉降速度Us均勻垂直滿流流入，從上方出口?Ωout流出，模擬中通過改變通道阻塞比a/R表征壁面對顆粒沉降的影響，考慮到流場的對稱性，只需求解如圖1所示的二維軸對稱流動，然后通過周向旋轉即可得到三維的顆粒在流體中沉降的流場。

圖1 數值計算選取的流動區域和計算區域Fig.1 The flow domain and computational domain for numerical simulation

采用有限元軟件COMSOL 的弱形式對控制方程進行求解，采用混合Galerkin 方法對偏微分方程進行離散，對速度分量采用二次形函數，對壓力和黏彈性應力采用一次形函數。通道入口和壁面采用恒定速度u=u(U,0,0)和零壓力梯度?p=0，出口采用壓力出口，顆粒表面流體速度u=u(0,0,0)，壓力梯度?p=0。對于黏彈性應力，入口處為均勻速度，出口處充分發展，黏彈性應力分量為零。容器壁和顆粒表面處溶液中表面活性劑分子只對壁面產生切向應力，其法向應力滿足n(n·τP)=0，n是邊界的單位法矢量。

3 流體流變參數擬合及顆粒沉降模擬

采用控制應力流變儀(HAKKE-RS300)對黏彈性表面活性劑溶液[（2.8%(質量)OTAC 和0.56%(質量) NaSal）]的線性黏彈性和穩態流變進行實驗研究，結果如圖2 所示，由小幅振蕩剪切實驗結果，可以獲得溶液流變參數松弛時間λ，通過穩態剪切實驗結果可以獲得表征剪切稀化效應的參量如Giesekus 模型中遷移因子α和FENE 模型中的拉伸率L，實驗獲得溶液流變參數的方法詳見文獻[16]。通過實驗發現流體具有剪切稀化效應，根據線性黏彈性，松弛時間λ=0.27 s，μ0=1.944 Pa·s，G0=7.285 Pa，采用Giesekus 本構模型進行參數擬合，擬合結果如圖2（a）所示，采用此流變參數對顆粒沉降進行模擬，顆粒半徑a為0.003 m，容器半徑R為0.031 m，顆粒密度ρs=7800 kg·m-3，對比不同遷移因子時的結果，模擬所得顆粒沉降速度隨時間的分布如圖2（b）所示，可以看出，表面活性劑溶液中顆粒沉降速度隨時間的變化明顯不同于牛頓流體，α較小時（剪切稀化作用較弱），顆粒沉降初期速度出現了明顯的波動性，隨著α增加，剪切稀化作用增強，顆粒沉降速度不穩定性減弱。

圖2 黏彈性表面活性劑溶液（2.8%OTAC+0.56%NaSal）流變性質實驗數據及Giesekus本構模型擬合結果(a)；顆粒沉降速度隨時間的分布(b)Fig.2 The rheological properties from experiments for viscoelastic surfactant solution (2.8%OTAC+0.56%NaSal)and the data fitting using Giesekus constitutive equation(a);The variation of settling velocity of a particle with time(b)

4 模型驗證

首先進行網格無關性驗證，考慮到球表面附近速度和應力梯度較大，在球表面處采用不同密度的網格M1、M2和M3，對應球表面網格數量分別為106、160 和320，總 網 格 數 量 分 別 為71818、106238 和112613，圖3 所示為采用不同網格數量時顆粒瞬時沉降速度的計算結果，可以看出，所采用的網格均可以很好地捕捉顆粒沉降的物理特征，通過對比顆粒沉降終速發現，采用M2網格相對于M1網格的顆粒沉降終速相對變化量為0.03%，而M3網格相對于M2網格的顆粒沉降終速相對變化量小于0.01%，據此采用M2網格作為計算網格，認為進一步加密網格對于模擬結果的影響可以忽略，對于不同a/R，球表面附近網格數量始終取160，而對遠離顆粒的區域采用相對較粗的網格，既能保證計算結果的精度，又能顯著減小計算量。

圖3 網格無關性驗證Fig.3 Grid independence verification

通過本研究模擬結果與前人研究對比，驗證本文采用的模擬方法的準確性。對比結果如圖4 所示，采用的本構方程為Giesekus模型，為了與文獻已有的結果（采用Oldroyd-B 模型）對比，假定本構方程中遷移因子α為0，此時，Giesekus 模型退化為Oldroyd-B 模型，可以看出，采用三種不同通道阻塞比a/R為0.1、0.2和0.5時，模擬獲得的顆粒沉降速度隨時間的變化與文獻中結果吻合較好，說明了本研究采用的模擬方法是可靠的。

由圖4可以看出，阻塞比a/R顯著影響顆粒沉降速度，當a/R為0.1 時，顆粒沉降速度在t=0.1 s 處出現一個明顯的過沖，即沉降速度很快增加到一個最大值，而后又逐漸減小，在0.4 s 后達到穩定的終端速度0.16 m·s-1。當a/R增加到0.2時，在t=0.1 s處仍然出現了速度過沖，其速度峰值小于a/R=0.1 時的值，而后速度減小，減小到某個值后，仍然沒有趨于穩定，而是隨著時間又逐漸增加，在0.6 s 后基本達到穩定的終端速度0.124 m·s-1，其小于a/R=0.1 時的終端速度值。當a/R進一步增加到0.5 時，壁面影響最強，此時在0.04 s處就出現速度最大值0.10 m·s-1，其最大值相對于a/R=0.1時的最大值減小約一半，沉降速度減小到0.022 m·s-1，之后沉降速度又有明顯增加，再次達到一個峰值，后又有輕微減小和增加，如此，經過兩次波動在0.5 s左右才達到穩定的終端速度0.038 m·s-1，其終速值也遠小于a/R為0.1 和0.2 時的終速值。由此說明隨著阻塞比a/R增大，即壁面效應的增強，顆粒沉降終速減小，而且顆粒沉降速度不穩定性增強，這是因為壁面效應較強時，顆粒沉降過程中應變率和應力較大，其彈性應力更大，彈性效應表現更為顯著，顆粒沉降速度經過波動之后，流體中積聚的彈性能量充分耗散后，才達到穩定的終端速度。

圖4 遷移因子α=0時Giesekus流體中顆粒沉降速度模擬結果與文獻[13,33-34]對比Fig.4 Comparison of settling velocity of a particle in Giesekus fluids(mobility factor α=0)between present simulation results and literature data[13,33-34]

5 結果和討論

5.1 Giesekus 和FENE-P 本構模型的流變性質模擬

本文采用的Giesekus 模型適用于半稀和濃溶液，FENE-P 模型適用于稀溶液，兩種模型都考慮了剪切稀化效應，為了明晰其流變性質對流體流動和顆粒沉降特性的影響，首先對Giesekus 和FENE-P本構模型的流變性質進行計算，考慮到所采用本構模型的非線性，采用數值方法求解在平面均勻剪切以及單軸拉伸流動下的流體應力分布。首先研究了松弛時間λ對流體剪切和拉伸流變性質的影響，如圖5所示，可以看出，松弛時間影響流體黏度曲線關于應變率的分布，當松弛時間較大時，出現剪切稀化和拉伸硬化對應的應變率較小，當松弛時間減小時，出現剪切稀化和拉伸硬化對應的應變率增大，說明了彈性對于應力-應變率的非線性性質也有影響。

圖5 松弛時間λ對流變性質的影響Fig.5 The effects of the relaxation time λ on rheological properties

此外研究了兩種本構模型中參數對流變性質的影響，如圖6、圖7 所示，此處采用無量綱化坐標，對剪切流動，橫坐標Wi=λγ˙，γ˙為剪切速率，縱坐標τXY為無量綱總剪切應力，對拉伸流動，橫坐標Wie=λε˙，ε˙為拉伸速率，縱坐標τXX為無量綱總拉伸應力，其中λ=1.0 s，μS/μ0=0.1。可以看出，對于兩種本構模型，在剪切速率0.1 s-1以上，都呈現出蠕蟲狀膠束溶液的典型的剪切稀化效應，并且同樣剪切速率，Giesekus 剪切稀化效應更強，遷移因子α在0.8 以上時，Giesekus 流體出現剪切帶。此外，對于拉伸流動，兩種本構模型均呈現拉伸硬化的特征，如圖6(b)和圖7(b)所示。

圖6 Giesekus流體在平面均勻剪切(a)和單軸拉伸流動(b)時的流變性質（Wi=λγ˙，Wie=λε˙）Fig.6 The dimensionless stress-strain relationships for Giesekus fluids under the planar shear (a)and the uniaxial tension(b)condition (Wi=λγ˙，Wie=λε˙)

圖7 FENE-P流體在平面均勻剪切(a)和單軸拉伸流動(b)時的流變性質（Wi=λγ˙，Wie=λε˙）Fig.7 The dimensionless stress-strain relationships for FENE-P fluids under the planar shear(a)and the uniaxial tension(b)condition(Wi=λγ˙，Wie=λε˙)

5.2 顆粒沉降特性

為了闡明流體流變性質對沉降的影響，先不考慮壁面影響，即a/R趨近于0，為此取圓柱通道的半徑R=200a，計算表明，當R/a從150增加到200時，模擬所得顆粒曳力系數只變化了0.1%，據此認為繼續增加R/a對計算結果影響可以忽略不計，因此取R/a=200，此時可以看作顆粒在無界流場中沉降，無壁面效應影響。研究顆粒浸沒在表面活性劑溶液中運動特性時，認為顆粒為完全潤濕[13-16]。此外，需要指出對于實際的顆粒懸浮體系，顆粒群效應和顆粒形狀的不規則均會對顆粒沉降特性產生影響，為了簡化模型，本研究考慮簡單球形顆粒的沉降來解釋顆粒沉降的非穩態效應。

5.2.1 流變參數對顆粒沉降速度的影響 在剔除了容器壁面效應的影響后，顆粒沉降速度只受表面活性劑溶液流變性質的影響，因此研究不同De、α以及L對顆粒沉降的影響。

先不考慮剪切稀化效應，即流體本構模型退化為Oldroyd-B 方程，此時顆粒沉降速度隨時間的變化規律如圖8 所示，取橫軸為無量綱時間因子t/λ，縱軸為無量綱速度u/Ut（Ut為斯托克斯沉降速率），當松弛時間λ為0.01 s 時，顆粒整體沉降速度較小，速度的增加較為緩慢，經過較長的時間才能達到一個穩定的終端速度，隨著λ增加，顆粒初期沉降速度逐漸增大，這是因為當顆粒沉降過程中對流體施加應變時，由于應力松弛效應黏彈性應力的增加有一個滯后，此時流體應力對顆粒作用力較小，因此沉降速度加速度較大，隨著時間增加，應力增加，對顆粒阻力作用逐漸增加，顆粒加速度逐漸減小，初始階段的應力松弛效應，導致了流體中沉降時顆粒速度的過沖現象，因此彈性效應是造成顆粒沉降不穩定性的原因。λ增加至0.05 s時，顆粒沉降速度出現了輕微的過沖現象，當λ繼續增加時，顆粒沉降速度過沖幅度顯著增加，并且波動周期數也增加，這在壁面約束時的顆粒沉降速度分布（圖4）上更加明顯。由無量綱時間看出，隨著λ增加，過沖量達到最大值的無量綱時間越小，達到終速的無量綱時間越小，這是由于隨著λ增加，流體彈性應力松弛效應增強，顆粒沉降過程中，流體發生變形后相應的變形應力增加緩慢，作用于顆粒的曳力增加緩慢，導致顆粒初始加速度和速度較大，因此產生較為明顯的速度過沖。這說明顆粒在沉降初期的不穩定性主要是由流體自身的彈性效應引起，彈性效應越強，顆粒沉降速度隨時間不穩定性越強。

圖8 無剪切稀化時松弛時間λ對顆粒沉降瞬時速度的影響Fig.8 The effects of relaxation time λ on transient settling velocity of a particle without shear thinning

除了彈性效應，對剪切稀化對顆粒沉降的影響進行了研究，隨著Giesekus 本構方程遷移因子α的增加，剪切稀化效應增強，由圖9 可以看出，隨著α增加，剪切稀化增強會削弱顆粒沉降速度的波動，當α大于0.3 時，顆粒沉降速度不再出現過沖現象，而是緩慢增加最后達到穩定的終端速度，這是由于剪切稀化效應導致隨著顆粒沉降速度增加，顆粒周圍剪切速率增大，表觀黏度逐漸下降，彈性能量釋放后的應力也會減小，因此顆粒會有一個逐漸加速的過程，剪切稀化效應越強，顆粒沉降終速越大。

圖9 Giesekus流體遷移因子α對顆粒沉降瞬時速度的影響Fig.9 The effects of mobility factor α on transient settling velocity of a particle in Giesekus fluids

對FENE-P 流體，隨著伸長率L的增加，剪切稀化效應減弱。由圖10 可以看出，拉伸率越小，顆粒沉降速度波動幅度越小，當L=2時，顆粒速度不出現過沖現象，因此彈性應力松弛效應增強顆粒沉降速度的波動，而剪切稀化削弱顆粒沉降速度的波動。

圖10 FENE-P流體拉伸率L對顆粒沉降瞬時速度的影響Fig.10 The effects of extensibility L on transient settling velocity of a particle in FENE-P fluids

在實驗中觀察到顆粒速度的持續擬周期性波動，可以推斷顆粒速度波動是由流體的彈性應力松弛效應觸發，彈性越強，顆粒沉降波動效應越強，而剪切稀化效應會削弱顆粒的波動，但是采用Giesekus 和FENE-P 本構方程模擬的結果發現顆粒在開始沉降后，速度波動很快消失，達到穩定的終端速度。黏彈性流動的高維森貝格數會引起求解不穩定，容易發散，本文在模擬中發現De超過10 后計算變得難以收斂，使得模擬更強彈性效應的沉降變得困難，可以推斷，隨著彈性進一步增加，顆粒沉降速度波動周期數和幅度均增加。另一方面，表面活性劑分子形成的蠕蟲狀膠束溶液由于其流變性質比共價骨架連接的聚合物分子更為復雜，其是一種活性聚合物，膠束之間不斷地斷裂和重構，其斷裂-重構效應會對顆粒運動產生不可忽略的影響，需要從流體中表面活性劑溶液形成的微觀膠束結構的動態演變規律去揭示彈性和膠束斷裂-重構的相互作用機理。

5.2.2 顆粒沉降負尾跡特性 表面活性劑溶液中顆粒沉降速度的波動可能與顆粒負尾跡有關，因此對顆粒負尾跡特性進行了研究，為了獲得顆粒負尾跡的流場特征，將模擬結果轉換到慣性參考系，即靜止流體中小球的自然下落，結果如圖11 所示，可以看出，顆粒在沉降過程中其后緣滯止點（相對于顆粒速度為零的點）附近的流體由于黏性作用隨著顆粒向下運動，前緣滯止點附近的流體被顆粒推動向下運動，但是在顆粒兩側表面附近以及向下游延伸的區域里，流體速度反而向上，這是由于此區域的流體受剪切作用較強，由于流體的剪切稀化效應，顆粒周圍向下的流體對周圍流體的黏滯作用力減小，因此速度最大的區域對外圍的影響區域也減小，強剪切帶把此區域分成兩個部分，就像顆粒帶動尾跡區將流體劈開一樣。此外，在尾跡區向下游延伸處，出現了一個速度方向與顆粒沉降方向相反的區域，如圖11(a)中“負尾跡”所示，這部分流體不斷地卷入尾跡區中的流體向上運動，經過一段距離后又被主流攜帶進入負尾跡旁邊的主流區，轉而速度向下，在顆粒下游處形成一個循環區域，形成負尾跡區，此處對附近流體的剪切和拉伸會造成大量膠束的斷裂和重構，使得流體內部的應力變得不穩定，是造成流場不穩定及顆粒沉降速度波動的根本原因。圖11(b)、(c)為沿對稱軸的顆粒尾跡區中心線上流體速度分布，縱坐標為速度大小，橫坐標為從顆粒后緣滯止點開始的距離，在沉降過程中顆粒表面附近流體由于流體黏滯作用隨著顆粒一起向下運動，速度為負，在某些區域出現速度大于零，流體速度與顆粒沉降速度相反，即為負尾跡。由圖11(b)可以看出，隨著α增大，剪切稀化和拉伸硬化效應增強，負尾跡區增大，并且距離顆粒更遠。對Giesekus 流體，當α為0.1 時，負尾跡區反向速度值最大，但是負尾跡向下游延伸較短，隨著α逐漸增加至0.7，反向速度值減小，負尾跡向下游延伸加長，隨著α繼續增加，負尾跡反向速度又逐漸增加。對FENE-P 流體，隨著拉伸率減小，剪切稀化增強，反向速度峰值和尾跡區向下延伸長度也是先增大后減小[圖11(c)]。

圖11 以速度矢量表示的顆粒下游負尾跡(a)；Giesekus模型遷移因子α對負尾跡的影響(b)；FENE-P模型拉伸率L對負尾跡的影響(c)Fig.11 The negative wake downstream of the particle presented using velocity vector(a);The effects of mobility factor α on negative wake for Giesekus fluids(b);The effects of extensibility L on negative wake for FENE-P fluids(c)

因此，流體剪切稀化和拉伸硬化共同導致顆粒周圍流場出現負尾跡，其對顆粒沉降速度不穩定性有著重要作用。

此外，研究了在同樣剪切稀化和拉伸硬化時彈性對負尾跡的影響，如圖12 所示，可以看出，對于同樣的剪切稀化和拉伸硬化（對Giesekus 流體遷移因子α取0.1，對FENE-P 流體伸長率L取5），當松弛時間增加，即彈性效應增強時，顆粒負尾跡反向速度峰值和負尾跡向下游延伸長度都增加，說明了彈性效應越強，負尾跡越強，因此流體的高彈性會強化顆粒負尾跡，其與剪切稀化和拉伸硬化共同導致了顆粒周圍流體流動以及顆粒沉降的不穩定，這與剪切稀化和拉伸硬化導致的微觀蠕蟲狀膠束的動態斷裂-重構有關，為完全準確描述蠕蟲狀膠束溶液中顆粒沉降的波動現象，需要進一步探討適用于微觀膠束動態斷裂-重構機制的本構模型。

圖12 黏彈性對顆粒沉降過程中負尾跡的影響Fig.12 The effects of viscoelasticity on negative wake during the particle settling

6 結 論

本文對黏彈性表面活性劑溶液中球形顆粒沉降特性進行了研究，同時考慮了溶液黏彈性效應和非線性應力應變響應特性，發現顆粒在沉降初期的不穩定性主要是由流體自身的彈性效應引起，彈性效應越強，顆粒沉降速度隨時間不穩定性越強，而剪切稀化效應會削弱顆粒沉降速度的波動。

此外發現顆粒在表面活性劑溶液中沉降時在流體剪切稀化和拉伸硬化作用下顆粒下游形成的“負尾跡”，而在純黏彈性流體無這一現象發生，驗證了實驗中對蠕蟲狀表面活性劑溶液中顆粒沉降時負尾跡存在的推測，并且發現剪切稀化和拉伸硬化效應增強會強化負尾跡效應。彈性效應增加時，負尾跡區反向最大速度和向下游延伸長度均增加，負尾跡區與顆粒下沉方向相反的速度分布可能是引起顆粒沉降速度持續波動的原因。

采用Giesekus 和FENE-P 模型雖然可以描述顆粒在經歷初始階段的沉降速度的波動性，但是在解釋顆粒連續周期性波動時仍有一定的局限性，這是因為表面活性劑溶液中形成的蠕蟲狀膠束是一種“活性聚合物”，不同于高分子聚合物溶液，這種流體在剪切拉伸時微觀尺度的膠束處于不斷地斷裂-重構的動態過程，而傳統的黏彈性本構模型不能描述在拉伸過程的斷裂-重構特性，因此下一步研究應考慮流體中微觀膠束的斷裂-重構動態過程對于顆粒沉降過程的影響，研究能更好地描述表面活性劑蠕蟲狀膠束溶液流變性質的本構方程。此外，對于懸浮顆粒體系，顆粒群效應和非球形不規則顆粒均會對顆粒沉降和運移產生重要影響，因此后續應對顆粒群效應和顆粒形狀的影響開展研究，獲得黏彈性表面活性劑溶液中顆粒沉降的普適規律。

符 號 說 明

a——小球半徑，m

D——圓柱通道直徑，m

d——小球直徑，m

ez——z方向單位矢量

F——慣性力，N

Fg——小球重力，N

Fz——慣性力在方向z的分量，N

fp——黏度系數

fr——松弛函數

G0——彈性模量，Pa

I——單位張量

L——伸長率

n——球表面法向矢量

p——壓力，Pa

R——容器半徑，m

Re——雷諾數

t——時間，s

Us——顆粒沉降速率，m·s-1

Ut——顆粒斯托克斯沉降速率，m·s-1

u——流體速度矢量，m·s-1

Vs——球形顆粒體積，m3

Wi——維森貝格數

α——遷移因子

γ˙——流體應變率，s-1

λ——松弛時間，s

μP——溶質黏度，Pa·s

μS——溶劑黏度，Pa·s

μ0——零剪切黏度，Pa·s

ρ——流體密度，kg·m-3

ρs——顆粒密度，kg·m-3

τP——黏彈性應力，Pa