似線性黏彈性模型 (Quasi-Linear Viscoelastic Model)：Part 1

Revised: 2021/12/23

原始的 Jeffreys model (Eq. 5.2-9) 只能定性描述線性黏彈行為

(5.2-9)

我們於是基於上式，將之推廣成似線性黏彈模型 (quasi-linear viscoelastic model)。做法是將微分型式的 Jeffreys model (Eq. 5.2-9) 中的時間偏導數 (partial time derivatives) 取代成對流時間導數 (convected time derivatives)，以得到 convected Jeffreys model，或稱為 Oldroyd's fluid B

(7.2-1)

[註：Equation 7.2-1 亦被稱為 contravariant convected Jeffreys model，其與 covariant convected Jeffreys model (或 Oldroyd's fluid A) 不同]

Convected Jeffreys model 有三個參數，零剪切率黏度 η₀、鬆馳時間 λ₁、延遲時間 λ₂。其運動學張量 (kinematic tensors) 已定義於 §6.1。以下三種模型為 convected Jeffreys model 的特殊例子:
(a) 如果 λ₂ = 0，變成 convected Maxwell model (也稱 contravariant convected Maxwell model 或 upper convected Maxwell model)，因其簡單的原故，常被用於黏彈流場之計算。
(b) 如果 λ₁ = 0，變成第二正向力係數為零的二階流體 (second-order fluid with a vanishing second normal stress coefficient)。
(c) 如果 λ₁ = λ₂ = 0，變成黏度 η₀ 的牛頓流體。

雖然我們知道 convected Jeffreys model 是可接受的，我們無法事先知道它是否能描述典型高分子的物質函數。因此，必需將之儘可能測試於不同的流變實驗 (如同對任何一種模型)。

以下針對 convected Jeffreys model，推導其於時間相依的剪切流場之預測 (time-dependent shearing flows)，包括： (a) 簡單剪切流場下之時間相依本質方程式；(b) 穩態剪切流場 (steady shear flow)；(c) 小振幅振盪剪切流場 (small-amplitude oscillatory shearing flow)；(d) 穩態剪切流場之前的起始過程 (start-up of steady shear flow)；(e) 穩態剪切流場之後的應力鬆馳 (stress relaxation following steady shear flow)。

A. 剪切流下之時間相依本質函數

根據 Appendix C，τ₍₁₎、γ₍₁₎、γ₍₂₎ 於非穩態剪切流場 (unsteady shearing flow) 的張量表示式分別如 Eqs. 1 到 3 所示

(1) 

(2) 

(3)

將 Eqs. 1 至 3 代入 convected Jeffreys model，可得到下方的矩陣方程式 (matrix equation)

(7.2-2)

從矩陣方程式中，我們可以得到一組應力張量分量的耦合微分方程式 (a set of coupled differential equations for the stress tensor components)

(7.2-3)

(7.2-4)

(7.2-5)

(7.2-6)

從這些方程式可得知，對於簡單的時間相依剪切流場，正向應力 τ_yy 和 τ_zz 均為零，故 Eq. 7.2-6 虛底線的項為零，我們接下來使用 Eqs. 7.2-3 至 6 計算特定的物質函數。

B. 穩態剪切流之物質函數

對於穩態剪切流，Eqs. 7.2-3 至 6 可剪化成代數方程式，解 τ_xx 和 τ_yx 的兩個方程式可得到以下函數

(7.2-7)

Convected Jeffreys model 預測黏度、第一正向力係數為定值，第二正向力係數為零。

C. 小振幅振盪剪切流之物質函數

為了得到小振幅振盪剪切的性質，我們假設應變 γ_yx 隨時間 t 的變化為

(4)

其中，γ₀ = γ̇ ₀/ω 是應變的振幅。接著找出在小應變極限下，Eq. 7.2-3 到 6 的解。在此流場下，剪切應力的微分方程式為

(7.2-8)

因為我們在求一個穩態週期解，一階線性常微分方程式的非齊性部分 (nonhomogeneous part of the first-order, linear, ordinary differential equation) 建議我們嘗試下方型式 τ_yx 的解 (或者說，在小應變的極限下，τ_yx 可表示成與應變異相、同相的應力疊加)

(7.2-9)

將 Eq. 7.2-9 代入 Eq. 7.2-8 可得

(7.2-10)

(7.2-11)

由上兩式可得 η' 和 η"  (或 G' = η"ω 和 G" = η'ω) 的表示式，兩者之值與應變大小無關 (η' 和 η" 的結果相同於 linear viscoelastic Jeffreys model)。

D. 穩態剪切流之前的起始過程

對於穩態剪切流場之起始，剪切率可表示為

(5)

其中，H(t) 是 Heaviside 單位階梯函數 (Heaviside unit step function)。Equations 7.2-3 和 7.2-6 (τ_xx 和 τ_yx)可因此表示成微分方程式如下 

(7.2-12)

(7.2-13)

其中，Dirac delta 函數 δ(t) 的引入是因為其為階梯函數的導數 (dH/dt = δ(t))。我們將 Eq. 7.2-12 乘上積分因子  exp(t/λ₁) (integrating factor)，並自 t = 0^- (τ_yx = 0) 積分至任一時間 t > 0，可得

(7.2-14)

接著將最後的結果結合 Eq. 7.2-12 可得

(7.2-15)

由 Eqs. 7.2-14 和 7.2-15 可知，convected Jeffreys model 預測應力和第一正向力的成長函數均不隨剪切率而變，然而這個預測僅與高分子溶液線性黏彈區間之響應吻合。再者，模型預測應力是隨時間呈單調遞增至穩態值，並無應力過衝 (stress overshoot)，此與實驗觀察不符。最後，應力在 t = 0 時有一個跳躍，這個特徵與模型的延遲時間 λ₂ 有關，但實驗上並未有此現象；不過，也因為有這個跳躍特徵，預測的應力比正向力成長更加快速，此趨勢與實驗結果定性相符。

E. 穩態剪切流之後的應力鬆馳

對於穩態剪切流之後的應力鬆馳，我們讓剪切率的切換以下式表之

(6)

接著將 Eq. 6 代入 Eq. 7.2-6 可得應力的微分方程式如下

(7.2-16)

若將 Eq. 7.2-16 自 t = 0^- (τ_yx = - η₀γ̇ ₀) 積分至任一時間 t > 0，可得應力鬆馳函數

(7.2-17)

 為了得到正向力的鬆馳函數，我們注意 Eq. 7.2-3 中，當時間大於等於零時且剪切率為零時，正向力 τ_xx 以 exp(-t/λ₁) 的型式自穩態值鬆馳，故得

(7.2-18)

如同前面已發現應力成長的物質函數不受剪切率而影響，應力鬆馳也是呈現同樣的結果。此外，τ_yx 一開始的鬆馳 (t = 0) 也是有一個跳躍特徵，這導致應力鬆馳到零的速度較正向力來得快，這點與實驗定性吻合。


Reference: RB Bird, RC Armstrong, O Hassager, Dynamics of Polymeric Liquids, Vol. 1, Fluid Mechanics, 2nd ed (Wiley-Interscience 1987).

靜態與流動態之高分子結晶 (Polymer Crystallization in Quiescent and Flow Conditions)

高分子加工 (polymer processing) 過程總是涉及流動、熱傳等現象，最終將高分子熔體製造成不同形狀的產品。而產品在不同位置的特性 (例如，機械性質)，將取決於材料之微結構 (microstructures)。因此，當我們模擬加工的過程，必需對結構的發展 (structure development) 加以考慮。不管是恆溫 (isothermal) 還是非恆溫 (non-isothermal) 狀態，高分子的靜態結晶 (quiescent state crystallization) 已被廣泛探討。反觀，由於理論和實驗的困難，鮮少學者對流動誘發造成的結晶 (flow-induced crystallization) 進行完整的研究。

本文根據下方文獻 [Guo 2000]，在此介紹一個流動誘發的結晶模型 (flow-induced crystallization model)，此模型基於一個假說，也就是當高分子處於受應力狀態下，平衡熔化溫度 T_m⁰ (equilibrium melting temperature) 將受影響而偏移至較高的溫度。值得注意的是，流動誘發的結晶模型仍處於發展階段，這裡僅提供其一做為參考。以下分別介紹 (I) 靜態結晶 (crystallization in quiescent state)；(II) 靜態結晶之參數決定 (determination of parameters for quiescent state crystallization)；(III) 剪切誘發結晶 (shear-induced crystallization)。

I. 靜態結晶 (Crystallization in Quiescent State)

Avrami 方程式 [1941] 經常被用於描述恆溫狀態的結晶動力學

(1)

其中，θ(t) 是時間 t 的相對結晶度 (relative crystallinity)，Χ(t)、Χ∞ 分別是時間 t 的絕對結晶度和最終絕對結晶度，n、k 分別是 Avrami 指數和 Avrami 恆溫結晶速率常數 (isothermal crystallization rate constant)。在另一方面，對於非恆溫的情況，Nakamura 等人 [1972] 根據 Avrami 理論得到下方的關係式

(2)

這是基於等動力 (isokinetic) 狀態且假設活化的晶核數目 (number of activated nuclei) 為常數。K(T) 是非恆溫結晶速率常數 (non-isothermal crystallization rate constant)，其與 Avrami 恆溫結晶速率常數的關係為

(3)

t_1/2 是結晶半數時間 (crystallization half-time)，1/t_1/2 是整體恆溫結晶速率 (overall rate of isothermal crystallization)。若假設成核點的數目 (number of nucleation sites) 與溫度無關且所有點均在同一時間被活化，Hoffman 等人 [1976] 得到溫度相依之整體恆溫結晶速率如下式

(4)

其中，T 是結晶溫度，R 是通用氣體常數，△T = T_m⁰ - T 是過冷度，f = 2T/(T + T_m⁰) 是熔化潛熱 (latent heat of fusion) 因溫度下降而減少的修正因子，U* 是結晶單元通過相邊界之活化能，T_∞ 是結晶完全停止的溫度，K_g 是球晶成長速率。根據 Hoffman 等人，U* 和 T_∞ 可分別被給定通用的值 6,284 J/mol 和 T_g - 30 K。(1/t_1/2)₀ 是與溫度無相依性的前因子 (prefactor)。

若對 Eq. 2 進行微分可得到微分型式的 Nakamura 方程式

(5)

對非恆溫的模流計算而言，Equation 5 較 Eq. 2 更加方便使用。我們可以注意到，一個非零的初始值需給定相對結晶度，值必需夠小且不會對最終結果造成影響，一般介於 10^-10 到 10^-15 之間。

II. 靜態結晶之參數決定 (Determination of Parameters for Quiescent State Crystallization)

於 Eqs. 3 和 4 中，如果我們給定 U* = 6,284 J/mol、T_∞ = T_g - 30 K，則僅需決定兩個未知參數 (1/t_1/2)₀、K_g。Figure 4.1 為恆溫結晶速率常數對溫度的作圖，圓點為 PET 的 DSC 實驗數據，實線為 Eqs. 3 和 4 的結合，透過 Fig. 4.1 可決定 (1/t_1/2)₀ = 4.25×10⁴ 1/s、K_g = 3.83×10⁵ K²。

Figure 4.1 恆溫結晶速率常數對溫度之作圖，實線為 Eqs. 3 和 4 的結合

至於 Avrami 指數 n 的決定，由 Eq. 1 可知，如果我們將 ln(- ln(1-θ)) 對 ln t 作圖，其斜率即為 n，如 Fig. 4.2 所示，指數 n 介於 1.5 到 2.5 之間。

Figure 4.2 相對結晶度對時間之作圖，實線為 Eq. 1

由前面的 Figs. 4.1 和 4.2，我們已自恆溫結晶的實驗決定了主要參數，即 (1/t_1/2)₀、K_g、n。接著我們可以將微分型式的 Nakamura 方程式 (Eq. 5)，對不同降溫速率的結晶過程進行預測，可於 Fig. 4.3 發現 Eq. 5 與實驗結果具高度吻合 。

Figure 4.2 相對結晶度對溫度之作圖，實線為不同降溫速率下 Eq. 5 之預測

III. 剪切誘發結晶 (Shear-Induced Crystallization)

於塑料加工過程中，由於高分子經歷很高的剪切率，外加上高分子長鏈之特性，將造成分子鏈具特殊方向性 (chain orientation)，規則的分子排列將加速高分子的結晶速度，並引響最終成品的形態 (morphology)。至於如何以最簡單的方式將剪切誘發結晶考慮在內？Guo [2000] 假設剪切力 (τ) 只會改變靜態下的平衡熔化溫度 (T_m⁰)，剪切造成的偏移量 T_shift 可近似為

(6)

所以，流動狀態之平衡熔化溫度 (T_m) 等於靜態之平衡熔化溫度 (T_m⁰) 加上 Eq. 6 的偏移量 (T_shift) ，即

(7)

Equation 6 的材料參數值 c₁ 和 c₂ 可以透過適當的流變實驗設計取得。若將 Eq. 4 中的 T_m⁰ 以 T_m 取代後，再結合微分型式的 Nakamura 方程式 (Eq. 5)，我們可對高分子流動狀態下的結晶動力進行預測。

Reference: J Guo, "Numerical simulation of stress-induced crystallization of injection molded semicrystalline thermoplastics," PhD dissertation (NJIT 2000).

非線性黏彈性模型 (Nonlinear Viscoelastic Model): Part 4

我們已於先前的兩篇文章中 (詳見 "非線性黏彈性模型: Parts 2 & 3")，介紹如何於本質方程式中納入不變量 (inclusion of invariants) (Part 2) 以及納入速度梯度的平方項 (inclusion of quadratic terms in velocity gradient) (Part 3)，藉以提升模型對真實流體數據之預測能力。這裡我們將介紹第三種方法，也就非線性應力項之納入 (inclusion of nonlinear terms in stress)。

Oldroyd 曾提出，似乎沒有太多流變的考量需要把本質方程式中的應力項侷限於線性。因此，接下來，我們將介紹一個具有非線性應力項的本質方程式 (a constitutive equation with nonlinear stress terms)，即 Giesekus 模型如下 (四個參數)

(1,2,3)

在這個模型中，應力張量 (τ) 是溶劑 (τ_s) 和高分子 (τ_p) 貢獻的疊加，此疊加型式相同於動力學中 (kinetic theory) 推導高分子溶液本質方程式的作法。

Giesekus 模型有四個參數，分別為鬆弛時間 λ₁、溶液和高分子對零剪切率黏度的貢獻 η_s 和 η_p、無因次的遷移率因子 α (mobility factor) 。涉及 α 的項與作用於組成高分子之非均向布朗運動 (anisotropic motion) 或非均向流體動力阻力 (anisotropic hydrodynamic drag) 有關。如果我們將 Eq. 3 的 τ_p 以 τ - τ_s = τ + η_s γ̇  取代之，Equations 1 至 3 可以合併寫成單一本質方程式如下

(4)

其中，η₀ 為零切率黏度、λ₂ 為延遲時間、a 為修正的遷移率，它們可由 η_s、η_p、α 等參數表之

(5)

由 Eq. 4 可看出，Giesekus 模型是 Oldroyd 8-constant 模型的一個特殊例子，差別僅在於額外加入一個 τ ‧ τ 的項。如果 a = 0，則 Eq. 4 還原成 convected Jeffreys 模型。因此 convected Jeffreys 模型可以表示成 Eqs. 1 至 3 的疊加且 α = 0。Table 7.3-4 列出幾個文獻常出現的本質方程式，這些方程式均來自於原始 Giesekus 模型之簡化。由此可知，原始完整的 Giesekus 模型對流變性質的預測較具多樣性。

Table 7.3-4 Giesekus 模型 (Eqs. 1 至 3 或 Eq. 4) 的簡化模型

由於 τ ‧ τ 的項之納入，使得 Giesekus 模型相較於 Oldroyds 8-constant 模型更加合理預測物質方程式。例如，(i) 對於 α ≠ 0 或 1，當 λ₂ = 0，黏度的冪次律斜率是 - ，這樣的結果有點不合理的陡峭 (即黏度對剪切率的相依程度過大)；但是，如果我們增加一個值很小的延遲項 (例如，λ₂/λ₁ = 10^-3)，則 dlogη/dlog(γ̇ ) 的值可以變成大於 -1，導致應力的大小總是隨著剪切率增加，較符合實驗趨勢。(ii) 第二正向力係數不是零且其值可相對於第一正向力係數而變，例如，Ψ_2,0 = - (α/2)Ψ_1,0。(iii) 如果 α ≠ 0，拉伸黏度是有界的 (bounded)，並於大的拉伸率時達一個常數。Figures 7.3-5 到 7.3-8 是 Giesekus 模型對各種流場下的物質方程式預測，一般而言，為了得到合理的性質，α 需落於 0 < α < 1/2。

非線性應力項的納入，確實使物質方程式的解析解 (analytical solutions) 難以得到，Table 7.3-5 列出幾個具解析解的物質方程式。然而，對於暫態 (transient) 的物質方程式，必需透過數值方法取得，下方參考文獻的附錄 C (Table C.2) 提供統衘方程式，它們可用於計算暫態物質方程式。

Figure 7.3-5

Figure 7.3-6

Figure 7.3-7

Figure 7.3-8

Table 7.3-5

Reference: RB Bird, RC Armstrong, O Hassager, Dynamics of Polymeric Liquids, Vol. 1, Fluid Mechanics, 2nd ed (Wiley-Interscience 1987).

非線性黏彈性模型 (Nonlinear Viscoelastic Model)： Part 3

於先前的文章 "非線性黏彈性模性 (nonlinear viscoelastic model): Part 2" ，我們已介紹如何於本質方程式中納入不變量 (inclusion of invariants)，以提升模型對真實流體數據之預測能力。這裡我們將介紹第二種方法，也就是速度梯度的平方項之納入 (inclusion of quadratic terms in velocity gradient)。

Oldroyd 認知到很多模型 (例如，convected Jeffreys 模型)，雖已納入速度梯度的平方項，但是做法並不是那麼的有系統。因此他將原始的 convected Jeffreys 模型進行廣義化 (generalization)，於該模型加入所有可能的平方項，像是 τ 與 γ₍₁₎ 的乘積、γ₍₁₎ 與 γ₍₁₎ 的乘積。Oldroyd 於是提出 Oldroyd 8-constant 模型如下 (八個參數為 η₀、λ₁、λ₂、...、λ₇)

(1)

Equation 1 中虛底線的項為原始的 convected Jeffreys 模型。很多源自於 Oldroyd 8-constant 模型 (Eq. 1) 的簡化版本已被大量用於文獻，特殊的值也已被指定於 Eq. 1 中的常數。Table 7.3-2 匯整了 Oldroyd 8-constant 模型之主要簡化模型。Table 7.3-3 則整理 Oldroyd 8-constant 模型的常見物質函數 (material functions)。為了讓模型預測與實驗物質函數至少有定性上的吻合，常數的值必需加以額外限制 (請參閱下方文獻更深入的探討)。整體而言，雖然 Oldroyd 8-constant 模型不偌 White-Metzner 模型那般成功地描述黏度和第一正向力係數，但是透過適當的參數選定後，Oldroyd 8-constant 模型卻能補捉穩態剪切流的起始過程發生的應力過衝行為 (stress overshoot  in the start-up of steady shear flow)。Figures 7.3-1 到 7.3-4 是 Oldroyd 4-constant 模型所預測的物質函數，此模型至少定性地描述大部分的流變性質。

Table 7.3-2 Oldroyd 8-constant 模型

Table 7.3-3 Oldroyd 8-constant 模型的物質函數

Figure 7.3-1

Figure 7.3-2

Figure 7.3-3

Figure 7.3-4

Reference: RB Bird, RC Armstrong, O Hassager, Dynamics of Polymeric Liquids, Vol. 1, Fluid Mechanics, 2nd ed (Wiley-Interscience 1987).