高英俊， 羅志榮， 鄧芊芊， 黃禮琳， 林 葵

（1.廣西大學廣西有色金屬及特色材料加工重點實驗室，南寧 530004；2.廣西大學物理科學與工程技術學院，南寧 530004）

韌性材料的微裂紋擴展與分叉的晶體相場模擬

高英俊1，2， 羅志榮1，2， 鄧芊芊1， 黃禮琳1， 林 葵1

（1.廣西大學廣西有色金屬及特色材料加工重點實驗室，南寧 530004；2.廣西大學物理科學與工程技術學院，南寧 530004）

采用晶體相場方法研究韌性單晶材料在雙軸拉伸條件下微裂紋擴展與分叉的演化過程，分析應變、溫度、初始裂口形狀等因素對裂紋擴展和分叉的影響.結果表明：對于簡單的單向拉伸，應變需要達到一定的臨界值，裂紋擴展才會啟動.對于二組相互垂直的雙軸拉伸作用，當應變達到臨界值后，裂紋擴展過程中會發生分叉現象.溫度越高，裂紋擴展越快且分叉越多.裂紋在擴展過程中，體系能量不斷降低，當裂紋出現分叉時，體系能量降低得更快.在裂紋擴展過程中，有時是會在裂尖處前方出現微小的空洞，類似在裂紋尖端前方出現位錯發射情況，這些微小的空洞逐漸擴大連成裂紋.本文所得結果與相關模擬結果和實驗結果符合.

晶體相場方法；微裂紋擴展；應變；韌性材料

0 引言

材料在介觀和宏觀尺度上的性能在很大程度上由其內部的微觀缺陷，例如，空位、空洞、位錯、晶界和微裂紋所決定.這些缺陷起因于在原子尺度發生的復雜非平衡動力學過程.就目前所采用的多尺度模擬方法來說，對這些微觀特征結構的實時演化過程的建模和模擬是一個重大挑戰.相場方法是當今研究微觀結構演化的強有力的數值計算方法［1－2］.傳統的相場方法［3－7］是建立在平衡態均勻場的基礎上的，忽略了許多由原子的周期排列結構產生的物理特性，難以反映晶體學結構特性以及原子尺度的行為信息，因而無法從根本上闡明微觀組織演化過程中原子尺度上的動力學機理.雖然，目前可用分子動力學方法［8－9］模擬預測這些微結構演化動力學過程，且能準確地描述這些現象，但該方法適用的特征時間尺度主要在原子振動時間尺度（10－12到10－14秒），不易拓展到擴散時間尺度（10－6秒），并且，分子動力學要求加載的應變率很大，通常達到107～1010s－1，與實際工業的應變速度102s－1相差甚遠［10，11］.最近，Elder等［11－13］基于密度泛函理論提出了晶體相場（phase-field-crystal，PFC）模型.該模型給出了新的自由能函數形式，引入的序參量為局域密度場，它將無序相（例如液相等）的密度場定義為一常數，將固相的密度場表示成周期性函數（波）的形式，進而通過周期的原子密度函數表現晶體的晶格結構.這樣的周期結構的密度場就很自然地與彈性效應、晶粒取向和位錯的運動等由周期結構產生的物理特性緊密地關聯起來［13］.PFC模型既可以描述晶體學結構特性以及原子尺度的行為，又可以揭示時間尺度為10－6秒量級的原子、缺陷運動的行為特征.由此可以預見，PFC方法在今后的微結構演化的研究中將展現其強大的優勢，能夠廣泛地應用到材料微結構研究的各個領域，為研究原子尺度組織結構的10－6秒量級的行為特征，例如研究晶體內部的空位、位錯的運動、晶界的遷移和微裂紋擴展及其對材料宏觀性能的影響，提供先進和強有力的計算方法.

目前，PFC方法已經有一些具體的應用，例如，模擬位錯滑移和攀移［10］、相結構轉變［14］、異質外延生長［15］、納米晶粒長大［16］等.PFC方法最新的研究前沿和發展方向［17－21］之一主要集中在：將加工變形的外界作用因素引入PFC模型中，建立形變和位錯動力學的PFC模型，用于超塑性、疲勞、斷裂、蠕變等行為和性能的機理研究.盡管傳統相場法和其它方法對材料的裂紋擴展研究已有較多的報道［22－25］，但是應用先進的PFC方法開展微裂紋擴展和分叉的研究工作還少見有報道.本文應用PFC方法模擬在雙軸拉伸應變作用下納米級微裂紋擴展與分叉的演化過程，詳細分析應變、溫度、初始缺口形狀等因素對微裂紋擴展和分叉的影響.

1 PFC模型與方法

1.1 體系的能量密度函數

與傳統的相場模型不同，PFC模型的自由能密度函數是通過對原子密度場變量取極小值得到的［11］.因此，該模型能夠揭示晶體學結構特性以及原子尺度的行為.對于固態金屬材料，其原子的位置不依賴于時間，呈規則排列.因此，要求相場變量必須能夠反映原子周期性排列的特征.引入周期性相場變量，其局域位置的最大值對應于原子的位置；另一方面，液相中的原子位置隨時間隨機變化，取時間平均可看成為常量.符合這兩方面要求的相場變量定義，可用保守的原子密度場函數作為相場變量，其表達式［11］可寫成

式中，等號右邊第一項反映晶格原子的周期排列結構特征，第二項反映液相的原子無規則運動，平均值為一個常量.此時體系無量綱的自由能函數［11］可以寫成

式中，γ是反映體系溫度的參數，?2為Laplace算子.該能量模型自洽地包含了晶體結構的物理特征，例如，多晶取向、彈性效應、塑性變形等.

1.2 二維體系的極小自由能密度函數

其中，ρ0和γ為體系自由能函數的兩個重要參數.由體系的極小自由能密度函數，可以計算并畫出體系不同相結構區的相圖.對于二維體系，相區有均勻無序相和固相，且固相又有三角格子相和條狀相兩種，文獻［11］已給出三角格子相和條狀相的自由能函數形式.利用這些相的自由能函數，再按照平衡相圖的計算方法得到的二維相圖如圖1（a）所示.其中，L、T和S分別代表均勻無序相、三角格子固相和條狀固相.

1.3 演化動力學方程與數值化處理

保守的原子密度場變量的演化可用與時間相關的Cahn-Hilliard動力學方程描述

式中，ζ為高斯隨機噪聲項，具有零平均值．本文不需考慮ζ的作用.

圖1 （a）PFC模型的二維相圖［11］；（b）（a）中方框的放大圖Fig.1 （a）2D phase diagram of PFCmodel（b）Enlarged map of the rectangle block in（a）

為求解復雜的動力學方程式（7），還必須將動力學方程在時間和空間上進行離散化處理，即采用數值求解的辦法.在本文的數值求解中，采用顯型Euler迭代公式［26］

式中，Δx為空間步長，j和n分別代表i的最近鄰格點與次近鄰格點.

1.4 施加應變的方法

在受到拉伸作用后，樣品中所有格點的坐標位置都會發生相應的變化.因此，可通過格點的坐標位置變換來反映應變的作用所引起的樣品變形.設變形前的樣品中格點位置的坐標用（x，y）表示，受拉伸作用引起變形后，格點位置的坐標用（x′，y′）表示．樣品在x方向和y方向分別受到應變εx和εy的拉伸作用下，變形前后坐標變換關系為x′＝x／（1＋εx）和y′＝y／（1＋εy）．由于應變作用引起樣品變形，其內部原子密度函數變為ρ（x′，y′） ＝ρ（x／（1＋εx），y／（1＋εy）），因此，在考慮應變作用引起變形的情況下，原子密度函數表達式（4）將變換為ρ（x′，y′），體系的能量也隨之發生相應的變化，但控制演化的動力學方程的形式保持不變．將樣品變形后的原子密度函數ρ（x′，y′）代入動力學方程（7）中，按照數值化計算步驟進行數值求解，然后進行可視化處理［5］，即可得到裂紋擴展的演化細節圖.本文模擬不涉及具體材料的物性參數，所用參數均已無量綱化處理［11］，并將連續空間離散為四方格子，采用周期性邊界條件.計算區域網格為1 024×512 gp（gp表示格子點數），空間步長設為Δx′＝Δy′＝π／3，時間步長Δt＝0.05，其它參數在下節中介紹.

1.5 初始裂口的設置

根據圖1給出的相圖，選取體系自由能密度函數的溫度參數γ＝－1.0，原子密度參數ρ0＝0.49（圖1（b）中的A點），位于液相與三角格子固相共存區.本文模擬單晶樣品，其原子的排列方向與x軸夾角θ為15°.在樣品的中心位置，挖出一個缺口作為初始裂紋的位置.為了便于比較，設置了兩種初始裂口，分別為20×10 gp的長方形缺口和直徑為16 gp的圓形缺口，缺口面積近似相等，其局部放大圖如圖2所示.缺口處周圍的原子環境如圖2（a）和2（c）所示，圖2（b）和2（d）為缺口周圍原子密度分布圖.缺口處的參數設為γ＝－1.0，ρ0＝0.79（圖1（b）中的B點），其正好處在液相與三角格子固相共存區附近的液相區.選取該點的參數的好處在于有利于三角晶相結構向裂紋結構轉變.

圖2 初始裂口的局部放大圖（εx和εy分別為x方向和y方向的應變.）Fig.2 Local enlarged zone of initial crack rim（εxandεyare strains in x and y direction.）

2 結果與分析

2.1 應變對裂紋的影響

設置初始裂口為長方形，如圖2（a）所示.圖3為樣品在y方向單軸拉伸作用下經過5.0×105ts（ts表示時間步數）演化的裂紋擴展情況，其中y方向的應變分別為εy＝7%、7.5%和8%.由圖3（a）可見，對于單軸拉伸情況，當εy＜7.5%時，裂紋不擴展；當εy≥7.5%時，裂紋才開始擴展．由此可見，裂紋沿x方向擴展的臨界應變約為7.5%.

圖3 5.0×105ts時裂紋擴展情況（εx＝0%）Fig.3 Crack propagation at t＝5.0×105ts andεx＝0%

設置初始裂口為長方形，如圖2（a）所示.樣品同時受到x和y方向的雙軸拉伸作用，設沿y方向的應變保持εy＝10%，沿x方向的應變分別為εx＝0%、3%、5%、7%、7.5%和8%，經過5.0×105ts演化的裂紋擴展分布如圖4所示.由圖4（a）可見，當εx＝0%時，裂紋擴展方向沿著原子間距最寬的方向進行，與x方向的夾角為15°，個別地方出現轉折（圖4（a1）），然后又恢復到原來的前進方向；圖4（a1）給出了圖4（a）的裂紋波折區域的放大圖.由圖4（a1）可見，裂紋的邊緣與原子的排列方向平行，在轉折處轉折方向與原方向夾角約為60°.當εx＜5%時，裂紋不出現分叉，存在鋸齒邊緣，裂紋擴展的形狀像一把“劍”.當εx＝5%時（圖4（c）），在裂紋前端轉折處，開始出現微小空洞（見圖4（c）中圓形區域）.當εx＞7%時，裂紋在擴展中開始出現分叉，εx越大，分叉越多．例如，當εx＝8%時（圖4（f）），由主裂紋分叉為子裂紋，子裂紋又分叉為孫子裂紋，形成了多枝狀結構，且裂紋越來越寬.由圖4中還可以看到在裂紋周圍出現灰白波紋分布，這是在裂紋周圍形成的原子位移的細微變化，反映了裂紋周圍彈性應變場分布發生的畸變.圖5是體系能量隨時間的變化曲線.由此可見，應變εx越大，能量降低就越快，表明裂紋擴展也越快，分叉也越多．當εx較小時，能量變化曲線近似為直線，表明這一過程中，裂紋擴展過程不分叉．裂紋要出現分叉，施加的εx存在一個臨界值，在本文設定的條件下，由圖5可見，該εx的臨界值接近7%．當εx＞7% 時，曲線下降速率明顯加快，表明開始出現裂紋顯著分叉．由圖4的裂紋變化情況，可以看到裂紋演化的形狀變化特征：開始時，施加的εx較小，裂紋形狀像一把“劍”；隨著應變εx增加，裂紋出現分叉擴展，逐漸長成“蝦”狀；最后，當εx≥8%時（圖4（f）和（g）），出現多個分支結構，裂紋形狀變成了“蟹”狀.文獻［28］也給出了類似的現象.

圖4 5.0×105ts時裂紋擴展情況（εy＝10%）Fig.4 Crack propagation asεy＝10%at t＝5.0×105ts

圖6給出了圖4（b）中虛線方形區域內裂紋尖端B的擴展速率曲線.由該圖可見，裂紋沿直線擴展過程中，速率較快，在轉折處，裂紋擴展速率減小；從轉折到恢復直線擴展，速率有一個加速突變的過程，然后又恢復到近似勻速運動過程.在轉折處，出現速率振蕩現象［29］，見圖6中速率的峰和谷.由圖5的變化趨勢和圖6速率的變化特點可以預測，當應變εx大于臨界值后繼續增加，則裂紋擴展速度加快，導致分叉更快、更多，將使得裂紋擴展速率出現更復雜的振蕩情況.類似文獻［29］指出的裂紋擴展速率超過臨界值后，發生裂紋分叉，裂紋擴展速率出現振蕩現象.

圖5 不同εx對應的能量曲線（εy＝10%）Fig.5 Energy curves atεy＝10%and differentεx

圖6 在圖4（b）中虛線方形區域內裂紋尖端B的裂紋擴展速率振蕩曲線Fig.6 Oscillations of crack propagation in Fig.4（b）rectangle block at crack tip B

2.2 溫度參數對裂紋擴展的影響

設置初始裂口為長方形，如圖2（a）所示，在裂口處取ρ0＝0.79，溫度參數分別取γ＝－1.0、－0.95和－0.90（分別對應圖1（b）中的B點、B1點和B2點）得到3類樣品.對樣品施加的應變為εx＝8% 和εy＝10%.圖7給出了演化4.0×10 ts時不同溫度參數取值情況下的裂紋擴展形貌圖.由圖7（c）可見，γ＝－0.9對應溫度較高的情況，裂紋擴展較快，分叉也較多；γ＝－1.0對應溫度較低的情況，裂紋擴展較慢，分叉也較少．圖8給出了三種溫度參數值情況下的體系能量隨時間的變化曲線．由圖可見，體系能量隨著演化時間增加而減小．γ＝－0.9時，能量曲線減小最快；γ＝－1.0時，能量減少最慢.在最初的演化時間里，能量曲線隨時間線性地減小，經過一段時間后，曲線開始向下彎曲，能量降低有加速趨勢，表明從此時開始，裂紋開始分叉，隨后的裂紋擴展過程中，出現更多的裂紋分叉；裂紋擴展分叉越多，體系能量降低得越快.能量曲線降低的本質原因是由于彈性應變能的釋放要比裂紋界面能（表面能）引起的能量增加要大得多，使得體系總能量隨著裂紋長大不斷降低，裂紋得以擴展.

圖7 4.0×105ts時取不同溫度參數值γ對應的裂紋擴展形貌圖Fig.7 Crack propagation with different temperature parameterγat t＝4.0×105ts

2.3 初始裂口形狀及周圍原子密度的影響

設置的初始裂口為圓形，如圖2（c）所示.同時施加的應變為εx＝8%和εy＝10%，溫度參數γ＝－1.0．為研究初始裂口處的原子密度的影響，裂口處的初始原子密度ρ0設置了兩種值，分別為0.79（圖1（b）中的B點）和0.90（圖1（b）中的C點）.圖10（a）和（b）分別給出了這兩種ρ0參數下演化5.0×105ts時的裂紋分叉形貌圖.從圖10可見，ρ0＝0.90時，裂紋分叉較多，分枝較豐富，沿x方向和沿y方向都發展較快，且發展較均衡；而ρ0＝0.79時，相對而言裂紋沿x方向的發展較快，而沿y方向擴展較慢．由圖9能量曲線變化可見，體系能量減小較快的是ρ0＝0.90情況．因為初始裂口處原子密度ρ0＝0.90比ρ0＝0.79偏離相圖的晶相與均勻無序相共存區更遠，不利于三角格子固相的穩定，而有利于裂紋的擴展.由圖10還可見，在裂紋的擴展過程中，多個裂紋的子分枝出現鋸齒折線結構.在分叉的地方，有的是先在分叉處前方出現微小的空洞，類似在裂紋尖端前出現位錯發射，然后這些微小的空洞逐漸連成裂紋.這是由于裂紋尖端處存在非線性彈性與塑性變形作用而引起的現象，文獻［30］也觀察到類似的結果.

圖8 取三種不同溫度參數值（γ＝－1.0，－0.95和－0.90）對應的能量曲線Fig.8 Energy curves with different parameter（γ＝ －1.0， －0.95 and －0.90）

圖9 初始圓形裂口處ρ0分別為0.79和0.90時的能量曲線Fig.9 Energy curves asρ0＝0.79 and 0.90 at circular crack rim

圖10 圓形初始裂口在應變εx＝8%和εy＝10%情況下的裂紋形貌圖（5.0×105ts）Fig.10 Crack propagation asεx＝8%andεy＝10%at circular crack rim

圖10（a）與圖4（f）所對應的設置條件，除了初始裂口形狀不同，其它條件相同.對比兩圖可以發現，圓形裂口的裂紋沿x方向的擴展占有較明顯的優勢，而長方形缺口的裂紋較均勻的向四周擴展.主要原因是PFC方程的高度非線性和初始裂口形狀不同將會導致裂紋周圍的原子作用環境不同，可見，初始裂口形狀，對裂紋擴展也有明顯的影響［24］.

3 結論

采用PFC方法能夠在原子空間尺度和擴散時間尺度較好地研究分析模擬韌性材料的裂紋擴展和分叉的詳細過程和細節.主要結論如下：

1）裂紋的擴展方向沿著原子間距最寬的方向進行，施加的應變存在一個臨界值，超過該臨界值，裂紋擴展才得以啟動；對于裂紋分叉，需要x和y方向的應變達到臨界應變值才能出現裂紋分叉.

2）體系的溫度參數取值較小時，則裂紋發展較快，分叉也較多.裂紋在擴展過程中，體系能量在不斷下降；當裂紋出現分叉時，體系能量減少得更快.裂紋的子分枝出現鋸齒折線結構，可看成裂紋擴展過程的振蕩現象.在裂紋擴展過程中分叉的地方，有時是先在裂尖處前方出現微小的空洞，類似在裂紋尖端前方出現位錯發射，然后逐漸連成裂紋.裂口處原子密度對裂紋擴展有明顯影響.

3）裂紋擴展的形狀變化特征：當施加單軸拉伸，應變較小時，裂紋形狀像一把“劍”；同時施加x和y方向的雙軸拉伸作用，隨著應變增加，裂紋形狀轉變成“蝦”狀；最后，裂紋出現多個分支結構，裂紋形狀變成了“蟹”狀.

［1］Moelans N，Blanpain B，Wollants P.An introduction to phase-fieldmodeling ofmicrostructure evolution［J］.Calphad，2008，32：268－294.

［2］Wang Y，Li J.Phase field modeling of defects and deformation［J］.Acta Mater，2010，58：1212－1235.

［3］Chen LQ.Phase-field models formicrostructure evolution［J］.Annu Rev Mater Res，2002，32：113－140.

［4］Steinbach I.Phase-field models in materials science［J］.Modell Simul Mater Sci Eng，2009，17：73001.

［5］高英俊，羅志榮，胡項英，等.相場方法模擬AZ31鎂合金的靜態再結晶過程［J］.金屬學報，2010，46：1161－1172.

［6］羅炳池，王海鵬，魏炳波.液態三元Ni-Cu-Pb偏晶合金相分離的相場模擬 ［J］.科學通報，2009，54：7－11.

［7］王錦程，李俊杰，楊玉娟，等.定向凝固界面形態演化及其穩定性的相場法研究 ［J］.中國科學：E輯，2009，38：16－23.

［8］Millett P C，Selvam R P，Saxena A.Molecular dynamics simulations of grain size stabilization in nanocrystallinematerials by addition of dopants［J］.Acta Meterialia，2006，54：297－303.

［9］Holm E A，Foiles SM.How grain growth stops：A mechanism for grain-growth stagnation in pure materials［J］.Science，2010，328：1138.

［10］楊濤，陳錚，董衛平.應力誘發雙位錯組亞晶界湮沒的晶體相場模擬［J］.金屬學報，2011，47：1301－1306.

［11］Elder K R，Grant M.Modeling elastic and plastic deformations in nonequilibrium processing using phase field crystals［J］. Phys Rev E，2004，70：051605－20.

［12］Elder K R，KatakowskiM，Haataja M，GrantM.Modeling elasticity in crystal growth［J］.Phys Rev Lett，2002，88：245701－5.

［13］Berry J，GrantM，Elder K R.Diffusive atomistic dynamics of edge dislocations in two dimensions［J］.Phys Rev E，2006，73：031609－16.

［14］高英俊，羅志榮，黃創高，等.晶體相場法研究二維六角結構向正方結構轉變 ［J］.物理學報，2013，62（5）：050507－010.

［15］Yu Y M，Backofen R，Voigt A.Morphological instability of heteroepitaxial growth on vicinal substrates：A phase-field crystal study［J］.JCryst Growth，2011，318：18－22.

［16］Gao Y J，Wang JF，Luo ZR，Lu Q H，Liu Y.Nano-twin structure simulation with phase field crystalmethod［J］.Chinese JComputational Physics，2013，30（4）：577－581.

［17］Stefanovic P，Haataja M，Provatas N.Phase field crystal study of deformation and plasticity in nanocrystallinematerials［J］. Phys Rev E，2009，80：046107.

［18］Chan PY，Tsekenis G，Dantzig J，et al.Plasticity and dislocation dynamics in a phase field crystalmodel［J］.Phys Rev Lett，2010，105：015502.

［19］Huang Z F，Elder K R，Provatas N.Phase-field-crystal dynamics for binary systems：Derivation from dynamical density functional theory，amplitude equation formalism，and applications to alloy heterostructures［J］.Phys Rev E，2010，82：21605.

［20］Jaatinen A，Ala-Nissila T.Eighth-order phase-field-crystalmodel for two-dimensional crystallization［J］.Phys Rev E，2010，82：061602.

［21］Hakim V，Karma A.Laws of crackmotion and phase-fieldmodels of fracture［J］.JMech Phys Solid，2009，57：342－368.

［22］Spatschek R，Brener E，Karma A.Phase field modeling of crack propagation［J］.Philos Mag，2011，91：75－95.

［23］高英俊，羅志榮，黃禮琳，林葵.韌性材料的微裂紋擴展和連通的相場方法模擬 ［J］.中國有色金屬學報，2013，23（7）：1892－1899.

［24］Huang S，Zhang S，Belytschko T，et al.Mechanics of nanocrack：Fracture，dislocation emission，and amorphization［J］.J Mech Phys Solids，2009，57：840－850.

［25］Livne A，Bouchbinder E，Svetlizky I，et al.The near-tip fields of fast cracks［J］.Science，2010，327：1359－1363.

［26］Chen L Q，YangW.Computer simulation of the domain dynamics of a quenched system with a large number of nonconserved order parameters：The grain-growth kinetics［J］.Phys Rev B，1994，50：15752－15756.

［27］Oono Y，PuriS.Computationally efficientmodeling ofordering ofquenched phases［J］.Phys Rev Lett，1987，58：836－839.

［28］Imbeni V，Kruzic J J，Marshall G W，et al.The dentin—enamel junction and the fracture of human teeth［J］.Nature Materials，2005，4：229－232.

［29］Sharon E，Fineberg J.Confirming the continuum theory of dynamic brittle fracture for fast cracks［J］.Nature，1999，397：333－335.

［30］Buehler MJ，Abraham F F，Gao H J.Hyperelasticity governs dynamic fracture at a critical length scale［J］.Nature，2003，426：141－146.

date：2013－08－21；Revised date：2013－11－25

Phase-Field-Crystal Modeling of Microcrack Propagation and Branching in Ductile Materials

GAO Yingjun1，2，LUO Zhirong1，2，DENG Qianqian1，HUANG Lilin1，LIN Kui1
（1．Guangxi Key Laboratory for Non-ferrousMetal and Featured Materials，Guangxi University，Nanning 530004，China；2．College of Physics Science and Engineering，Guangxi University，Nanning 530004，China）

Morphology ofmicrocrack propagation and crack branching of single crystal ductile materials under tensile are simulated with phase-field-crystalmethod.It shows that shape of notch has a strong effect on crack propagation.With uniaxial tension crack propagation appears if only strain reaches a critical value.With biaxial tension a crack branches if only strain reaches a critical value. Temperature parameter takes significant effects on crack propagation.The smaller the temperature parameter is（the temperature is higher），the faster the cracks propagate and themore the crack branches are.It is observed that free energy of system decreaseswith time and decreases faster after crack starts to branch.During crack propagation，disconnectedly isolated cavities around themain cracks can be found and these defects become new cracks.They grow up along a line and become a branching crack.Simulated results are in agreementwith other simulations and experimental results.

phase-field-crystalmethod；microcrack propagation；strain；ductilematerials

TG111.5

A

1001-246X（2014）04-0471-08

2013－08－21；

2013－11－25

國家自然科學基金（51161003）；廣西自然科學基金重點項目（2012GXNSFDA053001）；廣西大學廣西有色金屬及特色材料加工重點實驗室開放基金（GXKFJ12-01）及廣西大學科研基金（XJZ110611）資助項目

高英俊（1962－），男，博士，教授，博士生導師，主要從事材料微結構設計與計算模擬實驗研究，E-mail：gaoyj＠gxu.edu.cn