阮春蕾，劉春太

（1河南科技大學數(shù)學與統(tǒng)計學院，河南 洛陽 471023；2橡塑模具國家工程研究中心，河南 鄭州 450001）

?

剪切流場中聚乙烯結晶過程的建模與模擬

阮春蕾1，2，劉春太2

（1河南科技大學數(shù)學與統(tǒng)計學院，河南 洛陽 471023；2橡塑模具國家工程研究中心，河南 鄭州 450001）

基于Eder模型推導了剪切流場中球晶、串晶形態(tài)演化的數(shù)學模型，將第一法向應力差作為串晶成核的驅動，并引入兩相懸浮模型描述體系，認為其由無定形相和半結晶相組成，分別用FENE-P模型和剛性啞鈴模型描述。基于上述數(shù)學模型，分別構造了捕捉球晶、串晶生長的Monte Carlo法與體系控制方程求解的有限差分法，成功模擬了二維剪切流場中聚乙烯的結晶過程，給出了球晶、串晶的形態(tài)演化，分析了剪切時間、剪切速率對串晶成核密度、結晶速率、流體黏度等的影響。數(shù)值結果表明：所構造的Monte Carlo法合理有效，不僅成功捕捉了晶體的生長與碰撞，而且較為準確地預測了結晶速率。此外，提高剪切時間或剪切速率，將增加串晶成核密度、提高結晶速率、使流體黏度突增的時間點提前。

流動誘導結晶；形態(tài)學；動力學；Monte Carlo模擬；結晶

DOI：10.11949/j.issn.0438-1157.20151472

引　言

聚合物的結晶過程及其形態(tài)演化是決定制品力學性能的關鍵因素［1］。廣義而言，聚合物必須經過成型加工才能成為有用的制品。成型過程中，聚合物受到流場的作用，高分子鏈發(fā)生運動，在低溫條件下形成凝聚態(tài)結構，并最終形成球晶、串晶等復雜結晶形態(tài)［1］。由于成型中以剪切流場居多，因此，本文將研究限定在簡單剪切流場下，并以聚乙烯為代表，研究結晶形態(tài)與動力學，流體性質等的改變，揭示外場因素與結晶過程間的相關規(guī)律，并為復雜成型條件下的結晶過程提供指導。

流場中的聚合物結晶又被稱為流動誘導結晶（flow induced crystallization，F(xiàn)IC）［1-2］。流動誘導結晶的實驗研究表明［1］：與靜態(tài)結晶相比，流場中的結晶發(fā)生了實質性變化，不僅提高了結晶速率，而且改變了結晶形態(tài)（不同于靜態(tài)條件下的球晶，流場作用下可獲得高度異相的串晶結構）?；谏鲜鰧嶒灲Y果，許多學者提出了不少關于流動誘導結晶的動力學模型，其建立大多基于Nakamura方程［3］和Kolmogorov方程［4］。如，Doufas等［5］、Tanner等［6-7］通過修正Nakamura方程得到了流動誘導結晶動力學模型；Koscher等［8］、Zheng等［9］認為流動提高了晶體的成核數(shù)，并將其反映在Kolmogorov方程中，得到了流動誘導結晶動力學模型。該類方法的主要缺點在于不能提供較為詳細的微觀結構信息。Eder等［10］首次提出了基于結晶形態(tài)的數(shù)學模型。他們把球晶看成是不斷長大的球，串晶看成是不斷長大的圓柱體，得到一系列微分方程組。Zuidema等［11］對其進行了改進，用可恢復形變替代Eder模型中流場應變率作為串晶成核的驅動。他們的模型較為成功地反映了串晶的微觀信息，但相關研究并沒有顯式地采用模擬方法給出其成核、生長與碰撞的演化細節(jié)。因此，仍然擺脫不了使用結晶動力學模型。Guo等［12］引入了大分子鏈形變因子，他們認為當形變因子大于某臨界值時，結晶按串晶形式生長；反之，按球晶形式生長。他們利用該方法成功預測了iPP注塑制件的皮層厚度，但未得到結晶形態(tài)細節(jié)。

國內關于流動誘導結晶數(shù)值模擬方面的研究也不多見。王錦燕等［13］、榮彥等［14］對剪切流場中的聚合物結晶進行了模擬，定性分析了剪切對結晶過程及流體性質等的影響。他們的研究并未涉及到結晶形態(tài)。Zhou等［15］采用Guo等［12］的方法對成型中的結晶形態(tài)進行了顯示，區(qū)分了球晶和串晶，但相關研究并沒有將其與結晶動力學建立聯(lián)系，而只是給出了結晶形態(tài)的可能表示。

本文基于Eder模型推導了剪切流場中球晶、串晶的形態(tài)演化模型，并采用Monte Carlo法［16］成功模擬了二維剪切流場中聚乙烯的結晶過程，給出了球晶、串晶的形態(tài)演化，分析了剪切時間、剪切速率對串晶成核密度、結晶速率、流體黏度、體系應力的影響。

1　模型與算法

1.1結晶動力學模型

Avrami方程［17］是聚合物結晶動力學模型中應用最廣的模型，它也可用于描述流場中的聚合物結晶，其具體表達式為

其中，0φ為球晶的“虛幻體積”，0ψ為串晶的“虛幻體積”。Eder等［10］給出了0φ和0ψ的表達式，其中，描述三維球晶生長的方程組為

這里，Ns、Rtot、Stot、Vtot分別為球晶的晶核數(shù)、總體直徑長度、總體表面積、總體體積；而a為成核速率相關的參數(shù)，Gs為球晶的生長速率。三維串晶生長的方程組為

這里，Ns-k、Ltot、S~tot、V~tot分別為串晶的晶核數(shù)、總體長度、總體表面積、總體體積；而τn為與溫度相關的串晶成核率參數(shù)，為串晶成核的驅動參數(shù)，γ.為應變率，為擬合參數(shù)；τl為與溫度相關的串晶軸向生長率參數(shù)，為串晶軸向生長率的驅動參數(shù)，為擬合參數(shù)，Gs-k，r為串晶徑向生長速率。

Eder模型中將球晶視為不斷長大的球，而串晶視為不斷長大的圓柱體。本文研究二維情況下的晶體生長，此時，認為球晶按圓形生長，而串晶按長方形生長。

由式（2）可導出三維球晶的等價微分方程組，將其降為二維，則有

從中可知，球晶的生長由成核數(shù)sN及生長速率sG所決定。采用如下的成核密度公式［2］

其中，sN是過冷度TΔ的函數(shù)，過冷度TΔ定義為為平衡熔點，0N和?是經驗參數(shù)。球晶的生長速率用Hoffman-Lauritzen方程表示［18］

其中，G0、Kg為參考因子，U?為聚合物的分子活化能，Rg為氣體常數(shù)，T∞=Tg-30，Tg為玻璃化轉變溫度，

由式（3）可導出三維串晶的等價微分方程組，同樣，將其降為二維，有

這里，認為參數(shù)nτ=∞，lτ=∞［11］，則有

由式（8）可知，串晶的生長由串晶成核密度Ns-k、串晶軸向生長速率Gs-k，l以及串晶徑向生長速率Gs-k，r決定。串晶軸向生長速率Gs-k，l，可記為［10］

一般情況下，認為串晶徑向生長速率Gs-k，r與球晶生長速率Gs一致［11］，即

問題的關鍵在于求解串晶成核密度。采用如下公式［8］

其中，C為參數(shù)，N1為體系的第一法向應力差。

1.2無定形相與半結晶相的模型

式（11）中出現(xiàn)了第一法向應力差，因此，必須對晶體所在體系做一下闡述。這里，采用Zheng等［9］的思想，認為半結晶相懸浮于無定形相之中，無定形相采用FENE-P彈性啞鈴模型描述，而半結晶相采用棒狀啞鈴模型來描述。

無定形相的FENE-P彈性啞鈴模型為［19］

其中，C為大分子鏈構型張量，λa（T）為分子松弛

其中，，0aλ為溫度為0T時的分子松弛時間，ag/ER為擬合常數(shù)。

無定形相對應力的貢獻為［9，19］

其中，aτ為無定形相應力，n為分子密度，k為Botzmann常數(shù)。

半結晶相中的分子不能拉伸，只能取向，故采用硬棒模型來描述［9］，其動力學方程為

其中，為分子構型張量，γ.為形變率張量，λsc（x，T）為棒狀啞鈴的松弛時間，采用如下公式［9］α為相對結晶度，A為經驗值。式（15）中的四階張量<>RRRR，采用如下的二次封閉近似［9］半結晶相對體系應力的貢獻為［9］

這里，ηsc（x，T）為半結晶相的黏度，與無定形相的黏度間滿足如下關系［9］

至此，整個體系就完整描述出來。由式（14）、式（18）可得到體系的整體應力表達式而式（11）中的第一法向應力差N1=τ11-τ22由式（20）計算得到。

1.3算法分析

1.3.1Monte Carlo法本文采用Monte Carlo法［16］來捕捉結晶形態(tài)的演化。這里，研究對象設為空間中某塊特定區(qū)域上［0，1］mm×［0，1］mm的聚合物。給定溫度、剪切速率及剪切時間，并按式（5）、式（6）計算球晶成核密度Ns、球晶生長速率Gs；按式（9）～式（11）計算串晶軸向生長速率Gs-k，l、串晶徑向生長速率Gs-k，r以及隨時間變化的串晶成核密度Ns-k。

實施時，先將區(qū)域剖分為Ntot個小單元，這里取Ntot=106，并對每個單元賦底色。在t=0時生成密度為Ns的球晶晶核和密度為Ns-k的串晶晶核，并對晶核賦顏色。在t時刻，計算球晶的半徑R=Gst，對于每一個球晶，生成足夠多的隨機點，若該點的顏色為底色，則判斷其坐標是否落入該球晶范圍內，若是則認為該點被此球晶覆蓋，賦上該球晶的顏色；對于串晶，由于其成核與時間相關，先由式（11）計算串晶的成核密度，隨機產生新生核的坐標并賦值，同時，計算串晶的長度l=2Gs-k，l及半徑R=Gs-k，r，對于每一個串晶，生成足夠多的隨機點，判斷其是否落入某串晶長度及半徑覆蓋的長方形中，若是則賦上該串晶的顏色。如此循環(huán)，直到所有單元都被晶體覆蓋，則結束。需要指出的是，上述t~的計算與串晶的晶核生成時間相關，從晶核生成后開始計時。此外，相對結晶度的計算，這里采用式（21）來表示。

其中，晶體單元數(shù)可根據(jù)單元顏色來統(tǒng)計，而總體單元數(shù)為totN。Monte Carlo法的優(yōu)勢在于可以避開使用結晶動力學模型。常見的結晶動力學模型，如Avrami方程和Kolmogorov方程，在處理球晶、串晶共混體系時遇到困難，很難確定形狀因子mC和n。而Monte Carlo法直接從形態(tài)演化中統(tǒng)計獲得可信的相對結晶度，避免了由結晶動力學模型中參數(shù)的不確定而導致的誤差。

1.3.2有限差分法無定形相與半結晶相方程的計算采用有限差分法。演化方程式（12）和式（15）采用時間一階向前進行離散，即

2　結果與討論

2.1材料參數(shù)

本文以聚乙烯為代表，考察結晶過程及流體性質。模擬所用參數(shù)為：N0=17.4×106，?=0.155， G0=2.83×102m/s，，Kg=5.5×105K2，Tg=269 K，=2.69×10-7， C=105Pa-1·s-1·m-1，λa，0=4.00×10-2s，T0=476.15 K，Ea/Rg=5.602×103K，b=5，n= 1.26×1026/m3，k=1.38×10-23，β=9.2，β1=0.05，A=0.44。其中，結晶形態(tài)參數(shù)見文獻［8，11］，無定形相及半結晶相參數(shù)見文獻［9］。

2.2算法有效性驗證

為了驗證算法的有效性，將模擬結果與Eder模型預測結果進行比較。這里，假定球晶和串晶晶核瞬時產生，且有成核密度Ns=107m-2，Ns-k=107m-2，球晶按生長速率Gs=10-6m·s-1進行生長，串晶按徑向生長速率Gs=10-6m·s-1及軸向生長速率Gs-k，l=10-5m·s-1進行生長。圖1給出了模擬結果與理論模型預測所得相對結晶度的比較。由圖可知，模擬結果與理論預測值吻合得較好。模擬結果主要在結晶后期發(fā)生偏離，這可能是由于Eder模型的建立是基于Avrami方程的，而該方程在建模時并未考慮到晶體間的碰撞。這也是Avrami方程的主要局限性。

圖1　數(shù)值解與Eder模型的比較Fig.1　Comparison of simulation result with Eder model

2.3剪切時間的影響本節(jié)討論剪切時間的影響。為了便于分析，這里固定剪切速率γ.=10 s-1，溫度T=137℃。

2.3.1結晶速率及形態(tài)的變化圖2為剪切時間t=0、1、5、10 s時所生成的串晶成核密度。其中，t=0代表不對流體進行剪切，即靜態(tài)條件。施加剪切后，串晶的成核密度將有所增加；且剪切時間越大，串晶成核密度越大。當剪切停止后，成核密度將不再隨時間的變化而變化。

圖2　不同剪切時間下的串晶成核密度Fig.2　Nucleation density of shish-kebabs with different shear time

圖3　不同剪切時間下的相對結晶度比較Fig.3　Relative crystallinity comparison with different shear time

圖3為不同剪切時間下的相對結晶度演化。與靜態(tài)條件（t=0）相比，剪切時間越大，結晶速率越快。結晶的加速作用主要是由剪切誘導的串晶所貢獻的。由于流場的剪切作用，促使串晶晶核的產生，并提供軸向的生長速率，使串晶成長，進而加速結晶過程。在剪切時間較小的情況下，所提供的串晶晶核相對較少，其加速作用不甚明顯。文獻［9］也報道了大剪切時間對結晶速率的加速作用。一些學者［13-14］認為這種加速不是無限制的，當剪切時間達到某臨界值時，加速作用會變得緩慢。

圖4給出了剪切時間為0、5、10 s時的結晶形態(tài)比較。由圖可知，隨著剪切時間的增加，結晶形態(tài)發(fā)生較大變化：不僅提高了串晶成核密度，而且提高了串晶的各向異性。這種變化趨勢與Zhou等［15］的實驗結果也是一致的。

2.3.2流體黏度及應力的變化圖5給出了剪切時間t=0、1、5、10 s時流體黏度的變化。黏度隨著時間的增加在某臨界值時會發(fā)生突增現(xiàn)象。這是由相對結晶度引起的。由式（19）知，ηsc（x，T）=，當α接近常數(shù)A時，ηsc將趨向于無窮；故體系黏度急劇變化。此外，剪切時間越大，流體黏度越早發(fā)生突增。這與Zheng等［9］的結果一致。

圖6為剪切時間t=10 s時的剪切應力和第一法向應力差。應力值在前期變化較為緩慢，而到臨界值附近變化劇烈，該臨界值對應于相對結晶度。當結晶完成時，這些應力將鎖在材料內部，形成殘余應力［20］。

2.4剪切速率的影響

本節(jié)討論剪切速率的影響。為了便于分析，這里固定剪切時間t=10 s，溫度T=137℃。

2.4.1結晶速率及形態(tài)的變化圖7給出了剪切速率0γ=.、1、5、10 s-1時串晶成核密度的比較。其中，γ=.0代表不對流體進行剪切，即靜態(tài)條件。由圖可知，剪切速率越大，串晶成核密度越大。當剪切停止后，成核密度將不再隨時間的變化而變化。

圖8為剪切速率γ=.0、1、5、10 s-1時的相對結晶度發(fā)展。由圖可知，由于剪切作用的存在，結晶速率將比靜態(tài)條件下（0γ=.）的結晶速率有所增加，但這種增加在小剪切速率下（γ=.1 s-1）較微弱。同樣，結晶的這種加速作用也是由剪切誘導的串晶所貢獻的。

結晶形態(tài)的演化受剪切速率的影響與受剪切時間的影響類似。隨著剪切速率的增加，串晶成核密度將有所增加，并且串晶的各向異性也更為明顯。這里將不再作圖顯示，圖形與圖4類似。

2.4.2流體黏度的變化圖9為不同剪切速率下流體黏度的變化。流體隨時間的變化首先較為緩慢，但到某臨界值時發(fā)生突增，這也是由結晶度所引起的。此外，剪切速率越大，流體黏度突增的時間點越早。這與Zheng等［9］的結果也是一致的。

圖4　不同剪切時間下結晶形態(tài)的比較Fig.4　Morphology comparison with different shear time

3　結　論

本文給出了一個二維簡單剪切流場下聚乙烯結晶過程的建模與模擬?；贓der模型推導了球晶、串晶的形態(tài)演化模型；并用Zheng的兩相懸浮模型描述聚合物流體，將體系的第一法向應力差作為串晶成核的驅動；與此同時，構建了模擬球晶、串晶形態(tài)演化的Monte Carlo法，分析了剪切時間、剪切速率對結晶速率、結晶形態(tài)、體系黏度等的影響。得出如下結論。

（1）本文所建立的Monte Carlo法合理有效，成功區(qū)分了晶體生長模式，捕捉了晶體生長前沿，較為直觀地給出了晶體的演化發(fā)展；同時，統(tǒng)計獲得了可信的相對結晶度。該方法的優(yōu)勢在于成功避開了使用結晶動力學模型，避免由結晶動力學模型中參數(shù)的不確定而導致預測的失真。

（2）剪切時間對結晶速率、結晶形態(tài)影響重大；剪切時間越長結晶速率越快，串晶對結晶形態(tài)的影響越明顯，各向異性更為顯著。剪切時間對流體黏度也有顯著影響，這種影響主要體現(xiàn)在對結晶度的依賴上，當結晶度接近于某臨界值時，體系黏度將發(fā)生突增現(xiàn)象。

（3）剪切速率對結晶速率、結晶形態(tài)、體系黏度等的影響與剪切時間的影響類似；增加剪切速率能提高結晶速率、提高串晶的各向異性以及使體系黏度發(fā)生突增的時間點提前。

圖5　不同剪切時間下流體黏度的變化Fig.5　Evolution of melt viscosity with different shear time

圖6　剪切應力和第一法向應力差（剪切時間t=10 s）Fig.6　Shear stress and first normal stress （shear time t=10 s）

圖7　不同剪切速率下的串晶成核密度Fig.7　Nucleation density of shish-kebabs with different shear rate

圖8　不同剪切速率下的相對結晶度比較Fig.8　Relative crystallinity comparison with different shear rate

圖9　不同剪切速率下流體黏度的變化Fig.9　Evolution of melt viscosity with different shear rate

符號說明

b——啞鈴的有限拉伸參數(shù)

C——彈性啞鈴構型張量，m

Gs，Gs-k，l，Gs-k，r——分別為球晶、串晶軸向、串晶徑

向生長速率，m·s-1

Ns，Ns-k——分別為球晶、串晶成核密度，m-2

N1——第一法向應力差，N·m-2

R——剛性啞鈴取向向量，m

T——熱力學溫度，K

t——時間，s

Δt——時間步長，s

α——相對結晶度

ηa，ηsc——分別為無定形相、半結晶相黏度，Pa·sa

λ，scλ——分別為無定形相、半結晶相分子

松弛因子時間，s

τ，aτ，scτ——分別為總體、無定形相、半結晶

相應力，N·m-2

?——梯度算子

〈〉——系綜平均

下角標

a ——無定形相

s ——球晶

s-k ——串晶

sc ——半結晶相

References

［1］KENNEDY P， ZHENG R. Crystallization//Flow Analysis of Injection Molds ［M］. Munich： MunichCarl Hanser Verlag， 2013： 141-153.

［2］ZHE M， LUIGI B， GIUSEPPE P， et al. Flow induced crystallization in isotactic polypropylene during and after flow ［J］. Polym.， 2014， 55：6140-6151.

［3］NAKAMURA K， WATANABE T， KATAYAMA K， et al. Some aspects of nonisothermal crystallization of polymers （Ⅰ）：Relationship between crystallization temperature， crystallinity， and cooling conditions ［J］. J. Appl. Polym. Sci.， 1972， 16： 1077-1091.

［4］KOLMOGOROV A N. On the statistic of crystallization development in metals ［J］. Bull. Akad. Sci. USSR， Class Sci.， Math. Nat.， 1937， 1：355-359.

［5］DOUFAS A K， MCHUGH A J， MILLER C. Simulation of melt spinning including flow-induced crystallization （1）： Model development and predictions ［J］. J. Non-Newtonian Fluid Mech.， 2000， 92： 27-66.

［6］TANNER R I， QI F Z. A comparison of some models for describing polymer crystallization at low deformation rates ［J］. J. Non-Newtonian Fluid Mech.， 2005， 127： 131-141.

［7］WO D L， TANNER R I， FLETCHER D F. A numerical treatment of crystallization in tube flow ［J］. Prog. Polym. Sci.， 2012， 56： 1356-1366.

［8］KOSCHER E， FULCHIRON R. Influence of shear on polypropylene crystallization： morphology development and kinetics ［J］. Polym.，2002， 43： 6931-6942.

［9］ZHENG R， KENNEDY P. A model for post-flow induced crystallization： general equations and predictions ［J］. J. Rheol.， 2004， 48： 823-842.

［10］EDER G， JANESCHITZ-KRIEGL H， LIEDAUER S. Crystallization processer in quiescent and moving polymer melts under heat transfer conditions ［J］. Prog. Polym. Sci.， 1990， 15： 629-714.

［11］ZUIDEMA H， PETERS G W M， MEIJER H E H. Development and validation of a recoverable strain-based model for flow induced crystallization of polymers ［J］. Macromol. Theory Simul.， 2001， 10：447-460.

［12］GUO X， ISAYEV A I， GUO L. Crystallinity and microstructure in injection modelings of isotactic polypropylenes （Ⅰ）： a new approach to modeling and model parameters ［J］. Polym. Eng. Sci.， 1999， 39（10）： 2096-2114.

［13］王錦燕， 陳靜波， 劉春太， 等. 聚合物流動誘導結晶數(shù)值模擬 ［J］.化工學報， 2011， 62 （4）： 1150-1156. WANG J Y， CHEN J B， LIU C T， et al. Numerical simulation of shear-induced crystallization of polymer ［J］. CIESC Journal， 2011， 62（4）： 1150-1156.

［14］榮彥， 賀惠萍， 曹偉， 等. 基于兩相模型的聚合物流動誘導結晶數(shù)值模擬 ［J］. 化工學報， 2012， 63 （7）： 2252-2257. RONG Y， HE H P， CAO W， et al. Numercial simulation of flow-induced crystallization of polymer based on two-phase model ［J］. CIESC Journal， 2012， 63 （7）： 2252-2257.

［15］ZHOU Y G， TURNG L， SHEN C Y. Modeling and prediction of morphology and crystallinity for cylindrical-shaped crystals during polymer processsing ［J］. Polym. Eng. Sci.， 2010： 1-10.

［16］楊玉良， 張紅東. 高分子科學中的Monte Carlo方法 ［M］. 上海：復旦大學出版社， 1993. YANG Y L， ZHANG H D. Monte Carlo Method in Polymer Science［M］. Shanghai： Fudan University Press， 1993.

［17］AVRAMI M. Kinetics of phase change （Ⅰ）： General theory ［J］. J. Chem. Phys.， 1939， 7： 1103-1112.

［18］HOFFMAN J D， MILLER R L. Kinetics of crystallization from the melt and chain folding in polyethylene fractions revisited： theory and experiment ［J］. Polym.， 1997， 38： 3151-3212.

［19］NASIR A， MUHAMMAD A J， MUHAMMAD S. Theoretical analysis of the exiting thickness of sheets in the calendaring of FENE-P fluid ［J］. J. Non-Newtonian Fluid Mech.， 2015， 225： 28-36.

［20］ABOUHAMZEH M， SINKE J， JANSEN K， et al. Closed form expression for residual stresses and warpage during cure of composite laminates ［J］. Compos. Struct.， 2015， 133： 902-910.

Numerical simulation of morphology and kinetics of polyethylene in shear flow

RUAN Chunlei1，2， LIU Chuntai2
（1School of Mathematics and Statistics， Henan University of Science and Technology， Luoyang 471023， Henan， China；2National Engineering Research Center for Advanced Polymer Processing Technology， Zhengzhou 450001， Henan， China）

The mathematical model of morphology evolution of spherulites and shish-kebabs in the shear flow is deduced based on the Eder model. The model considers that the nucleation and growth of spherulites are determined by the static temperature， while the nucleation and growth of shish-kebabs are determined by the flow，which is depended on the first normal stress difference of the system and the shear rate respectively. In order to calculate the nucleation density of the shish-kebabs， the two-phase suspension model of Zheng is introduced. The model treats the stress as the combination of the amorphous phase and semi-crystalline phase. The amorphous phase is described by FENE-P model while the semi-crystalline phase is depicted by a rigid dumbbell model. Based on the mathematical model， the Monte Carlo method and the finite difference method are constructed，respectively. The former is to capture the crystal growth while the latter is to calculate the equation of the system. By using these methods， the 2D crystallization in the shear flow is simulated. The evolution of spherulites and shish-kebabs is given. Also， the effects of shear time and shear rate on the nucleation density of shish-kebabs，crystallization rate， the viscosity of the fluid and the system stress are discussed. Numerical results show that the Monte Carlo method is valid which not only captures the morphology evolution of crystals successfully， but alsopredicts the crystallization rate well. In addition， the increase of the shear time or the shear rate will increase the nucleation density of shish-kebabs and the crystallization rate.

date： 2015-09-21.

RUAN Chunlei， ruanchunlei622@mail.nwpu. edu.cn

supported by the National Natural Science Foundation in China （11402078， 51375148）.

flow induced crystallization； morphology； kinetics； Monte Carlo simulation； crystallization

O631.+3

A

0438—1157（2016）05—2144—08

2015-09-21收到初稿，2016-01-11收到修改稿。

聯(lián)系人及第一作者：阮春蕾（1983—），女，博士，副教授。

國家自然科學基金項目（11402078，51375148）。