馬氏體相變的相場模型研究進展

于繼東，姚松林，吳 強

（中國工程物理研究院流體物理研究所，四川 綿陽 621999）

相場（Phase field）模型是目前模擬材料微結構演化的重要理論與計算工具[1–5]。其優(yōu)勢主要體現在以下兩點：（1）避免了三維空間中顯示追蹤界面在數學處理上的困難，相場模型通過場變量的連續(xù)變化將尖銳界面擴散為具有一定厚度的界面，界面的演化通過場變量的連續(xù)變化來實現，無需顯示追蹤界面，因此在處理復雜微結構演化問題上具有獨特的優(yōu)勢；（2）相場模型基于能量泛函的表述形式，其物理模型具有良好的擴展性，可統(tǒng)一地在能量函數中引入力、熱、電、磁等多場物理機制。目前，相場模型已廣泛應用于計算材料學、計算固體力學等諸多領域，如凝固[6]、馬氏體相變[7–9]、鐵電鐵磁相變[10]、孿晶變形[11–12]、位錯動力學[13–14]、裂紋萌生與擴展[15–19]等。本文將側重于介紹相場模型應用于馬氏體相變中的相關進展。

馬氏體相變是一種以切變?yōu)橹鞯臒o擴散型、固固、一階位移型相變，相轉變過程極快（新相長大速度接近聲速），并伴隨著豐富的微結構生成，這些微結構演化與其宏觀物理力學性能緊密相關，因此針對馬氏體相變微結構的形成與演化開展了大量的研究工作。Ball 等[20–21]、Kohn 等[22]、Bhattacharya[23]發(fā)展了零應力平衡態(tài)下微結構形成的理論描述，對于理解微結構模式的形成具有重要的指導意義。然而，這些理論無法回答馬氏體相變的動力學問題。相變動力學關注相變的轉變機制和轉變速率，以及其他因素（如應力、溫度等）產生的影響[24]。馬氏體相變過程中微結構演化與彈性應變場的相互作用非常復雜，數值模擬是研究此類問題的主要方法[7]。相關研究主要集中在不同壓力、溫度等外在條件下不同材料的馬氏體相變的微結構演化模式，以揭示物理機制，并指導材料設計[25–31]。

本文首先概述相場模型的基本概念以及馬氏體相變相場模型的理論框架，圍繞相場模型的兩個獨特優(yōu)勢，介紹近年來馬氏體相變相場模型的主要進展：一是基于相變路徑理論建立相場模型，該模型對于加卸載條件下重構型馬氏體相變的復雜微結構演化問題具有特別的適用性；二是基于相場模型在處理多物理過程中的良好擴展性，介紹馬氏體相變與塑性變形耦合的相場模型的系列進展。

1 相場模型

相場模型是為研究非均勻系統(tǒng)的微結構演化問題而發(fā)展起來的唯象理論模型，其基本思想可追溯至19 世紀末van der Waals 提出的擴散邊界的思想[32]，Gibbs[33]和Cahn 等[34]的經典工作奠定了相場模型的基礎物理概念。類似于Landau 在處理連續(xù)相變問題時引入序參量的概念，相場模型引入了一系列序參量來描述非均勻系統(tǒng)中不同的相，系統(tǒng)能量不僅是序參量的函數，也是序參量空間梯度的函數，序參量的空間梯度描述的是不同相之間的界面。序參量可以有明確的物理意義，如金屬合金相變問題中的質量分數、鐵電及鐵磁問題中的電偶極矩和磁矩，其熱力學和動力學參數均有明確的物理量相對應；場變量也可以沒有明確的物理意義，例如研究凝固問題中的固液界面時，場變量的引入是為了避免顯示界面追蹤，在這類問題中定義具有嚴格物理意義且數學上易于處理的序參量仍是一個挑戰(zhàn)。守恒型序參量（如濃度）和非守恒型序參量（如馬氏體變體取向）的控制方程分別為Cahn-Hilliard 方程[34]和Allen-Cahn 方程[35]。隨著計算機能力的提升，從20 世紀90 年代開始，相場模型逐漸應用于凝固問題的模擬，并在涉及微結構演化一類的問題中得到廣泛應用，相關研究進展可見綜述[1–5]。

式中：f為局域自由能密度， αi和 βij為與界面厚度相關的梯度項系數。界面處的序參量梯度不為零，式(1)積分項中的梯度項與界面能相關。不同相場模型的主要區(qū)別體現在如何構建式(1)所示的系統(tǒng)自由能。例如：相的失穩(wěn)分解和凝固問題中構建雙勢阱形式的能量函數[36–37]，位錯動力學問題中構建周期性勢阱能量函數[38]，馬氏體相變問題中構建滿足晶格對稱性操作的能量函數[39]。

守恒型場變量的演化方程為Cahn-Hilliard 方程

非守恒型場變量的演化方程為Allen-Cahn 方程

式中：Mi j和Lpq為動力學系數，與原子或界面的遷移率相關，Mij和Lpq滿足Onsager 關系，即Mij=Mji，Lpq=Lqp。Cahn-Hilliard 方程通過ci的守恒方程推導得到，體現了守恒型場變量的守恒性質。Allen-Cahn 方程是個唯象方程，其物理意義是場變量的演化朝著系統(tǒng)自由能減小的方向進行。

物質界面具有一定的寬度，如相界面厚度的典型尺度為納米量級，當相場模型的擴散界面寬度與界面的物理尺寸相當時，相場模型的界面具有明確的物理意義，相關參數也可通過第一性原理計算獲得[40–41]（如層錯能等），對缺陷等微結構的演化具有預測性[5,8,42]，此類相場模型被稱為微觀相場模型（Microscopic phase field）[5]。在相場模型的數值模擬中，界面至少需要3～4 個計算網格，因此微觀相場模型適用的空間尺度為微觀尺度。目前，相場模型也廣泛應用于介觀尺度上的相變、位錯結構、枝晶生長等問題中，此時擴散界面的寬度遠大于物理界面的寬度（一般為1 個數量級），相場模型并不具備預測缺陷基本特性的能力，但在描述微結構演化過程中的復雜幾何結構變化及微結構模式形成等問題時，仍能給出與實驗結果一致的有價值的認識，此類相場模型稱為粗粒化相場模型（Coarse grained phase field）[5]。理論上，粗?；鄨瞿Ｐ蛢H在尖銳界面極限（即通過尖銳界面的漸進分析對模型參數及控制方程進行修正）下才能給出定量化結果。Karma 等[43–44]將尖銳界面漸進分析方法應用于枝晶生長問題，將計算網格由納米量級提升至0.1 μm量級仍可保證模擬結果的定量程度，極大地拓展了枝晶生長問題的相場模擬能力。最近，Finel 等[45]提出了一種新的方法，以提高粗粒化相場模型的界面動力學模擬精度，通過在離散形式的能量方程中引入晶格的平移及旋轉不變性，從而在能量方程中引入對稱性不變的限制條件，使得在一個計算網格內高精度地捕捉界面動力學特性，相關方法對于提升固固相變相場模擬的研究尺度具有重要意義。另一方面，也可通過自適應模擬技術[46–47]提高相場模擬的計算效率。有效提升粗?；鄨瞿Ｐ偷念A測能力及計算效率將是相場模型研究中持續(xù)關注的問題之一[48–51]。

2 馬氏體相變的相場模型

馬氏體相變有不同的分類[52]。Bhattacharya 等[53]根據逆相變的微結構可恢復性，將馬氏體相變分為弱（Weak）馬氏體相變和重構型（Reconstructive）馬氏體相變。弱馬氏體相變指的是外載卸載時，逆相變使相變微結構消失的同時，在母相中不會或僅殘留少量位錯和孿晶缺陷，如圖1 所示的TiNi 形狀記憶合金的掃描電鏡（Scanning electron microscope，SEM）圖像，降溫時發(fā)生高溫奧氏體到低溫馬氏體相變，形成復雜的V 形及六角形形貌，再升溫后馬氏體相變微結構基本消失[54]；重構型馬氏體相變指的是逆相變不會使馬氏體相變微結構消失，外載卸載過程會形成一系列微結構，導致母相中形成大量的位錯和孿晶缺陷，如金屬鐵在沖擊加載條件下的相變過程，圖2 所示的沖擊加載金屬鐵回收樣品的透射電鏡（Transmission electron microscope，TEM）圖像[55]顯示，回收樣品中存在大量孿晶及二次孿晶組織，而初始樣品中不存在這些孿晶組織，這些特征結構是材料經歷相變-逆相變后形成的不可恢復的微結構。Bhattacharya 等[53]從晶格對稱群的角度分析發(fā)現，當新相的對稱群是母相對稱群的子群時，馬氏體相變微結構是可逆的，并稱此類相變?yōu)槿躐R氏體相變，其他類型的馬氏體相變?yōu)橹貥嬓婉R氏體相變。采用相場模型模擬兩種類型的馬氏體相變時，理論模型上是有明顯區(qū)別的。

圖1 TiNi 合金降溫-升溫過程中的原位SEM 圖像：(a) 290 K，(b)降溫至222 K，(c)再升溫至265 K，(d)再升溫至290 K[54]Fig. 1 Series of in situ SEM images of TiNi alloy at (a) 290 K and (b) 222 K upon cooling, and (c) 265 K and (d) 290 K upon heating[54]

相場模型的能量函數一般基于Landau 相變理論，將自由能在相變點附近展開為序參量的多項式函數。對于馬氏體相變問題，按序參量的選擇方式分為兩類相場模型。第一類是將序參量與總應變直接關聯，能量函數為序參量的多項式函數[56–58]，或者進一步考慮能量函數在晶格對稱性操作下保持不變[59–61]。此類相場模型的前提條件是晶格尺度的變形與宏觀尺度的應變一致，即Cauchy-Born 假設成立，對于存在塑性變形等復雜變形的問題，此類相場模型不適用。另一類是將序參量與相變應變（也稱本征應變或特征應變）關聯[7,26,62–65]，相變應變?yōu)榱銘l件下相變導致的應變，序參量為內變量，在能量函數構建中可考慮不同內變量之間的耦合，使相場模型更易于向多場問題等復雜情況拓展。目前馬氏體相變相場模型研究中主要采用第2 類方法，下面將對其進行簡要介紹。

圖2 沖擊加載下金屬鐵回收樣品的TEM 圖像[55]Fig. 2 TEM image showing the microstructure of shock-compressed iron[55]

馬氏體相變相場模型的能量函數由化學能和彈性能兩部分組成?；瘜W能描述系統(tǒng)的熱力學特性，不依賴于應力，可表示為[8]

例如，對于“hcp→正交相”馬氏體相變，滿足晶格對稱性操作且展開至6 次冪的化學能表示為[66]

彈性能可表示為

Cijkl(η)

式中：彈性常數為母相彈性常數和新相彈性常數的差值函數

對于多序參量問題，兩種限制條件使得能量函數的構建較為復雜：一個是需要考慮不同馬氏體變體之間轉變的重取向問題；另一個是能量函數需在數學形式上保證兩個序參量不能同時到1，即馬氏體變體不能相互重疊。一般通過引入罰函數，使2 個及以上的序參量同時到1 的能量極高，從而使能量函數滿足限定條件，然而隨著序參量的增加，數學形式變得非常復雜，限制了這種做法向更多序參量情況的拓展。具體做法可見Steinbach[4]和Levitas 等[63]的工作。

3 基于相變路徑理論的相場模型

3.1 模型簡介

第2 節(jié)介紹的相場方法已成功應用于形狀記憶合金等功能材料的弱馬氏體相變[5,9]，但應用于重構型馬氏體相變時卻遇到較大困難。圖3 顯示了金屬鐵bcc（ α相）→hcp（ ε相）重構型馬氏體相變加卸載過程[70]，加載過程的 α→ε相變過程存在6 個對稱性相關的ε相變體（若考慮晶格內的shuffle 變形，則存在12 個等效的ε相變體），卸載時的逆相變過程不會沿著正相變的路徑返回，而是形成新的bcc 相變體，存在12 個等效的 α′相變體，這種“相變-逆相變”路徑的不唯一性導致大量孿晶及二次孿晶等組織結構的出現（見圖2）。采用相場模型描述金屬鐵在加載過程中的馬氏體相變問題時，需要引入19 個序參量，并且考慮多次加卸載過程時，變體數目呈幾何級數增長，若采用對每個變體引入一個序參量的做法，將導致相場模型能量函數的數學形式非常復雜，數值求解變得非常困難。

圖3 鐵的α→ε正相變過程（紅線）及ε→α′逆相變過程（藍線）中對稱性相關的變體示意圖[70]Fig. 3 Schematic illustration of the multiple symmetry-related variants for the forward α→ε (red) and the reverse ε→α′ (blue) phase transitions in iron[70]

式(12)的分解方式并不唯一，不同分解順序的差別在文獻[74]中已有分析。在小應變近似下，式(12)可退化為類似于式(8)的分解形式。將式(12)對時間求導，得到

系統(tǒng)的Helmholtz 自由能為彈性變形梯度、相變變形梯度和相變變形梯度的空間梯度的函數，可表述為[70]

系統(tǒng)的耗散能大于或等于零，即Clausius-Duhem 不等式

式中：第1 個積分項為彈性能（不包括黏彈性），不發(fā)生耗散，由超彈性本構描述；第3 個積分項為塑性耗散能，可由塑性屈服準則或晶體塑性模型描述；第2 個積分項和第4 個積分項為相變耗散能，借助Gauss 積分定理，可得與對偶的相變驅動力

在界面各向同性的假設下，梯度項自由能表示為

式中： β類似于式(1) 中的梯度項系數。相變變形梯度的演化依賴于相變驅動力，表述為Ginzburg-Landau 方程形式[70]

式中：L類似于式(3)中的動力學系數。與經典相場模型中每個變體對應于一個序參量的做法不同，式(20)中相變變形梯度Ft（Ft為非對稱矩陣，有9 個分量）為序參量。相變變形梯度演化的關鍵是相變驅動力的求解，需結合超彈性本構、塑性本構及相變本構耦合求解[70]。小變形下，應力-應變關系可近似為線性關系，大變形條件下該線性關系不再滿足，需采用非線性超彈性本構模型[75]。非線性超彈性本構會涉及如何確定高階彈性模量中的相關參數，若再考慮熱彈性的貢獻以及彈性模量對溫度和壓力的依賴關系，則參數量進一步增加，不便于應用。對于高壓加載情況，位錯滑移、孿晶變形等塑性變形機制將剪應力水平松弛至材料的屈服強度附近，使得偏應力與靜水壓相比為小量，此時在本構模型中引入標量形式的狀態(tài)方程，高壓下仍具有較好的精度，并減少了參數量，是一種良好的近似處理。但是當加載壓力與材料的屈服強度可比擬時，需采用非線性超彈性本構。

Vattré等[70,72]采用相變路徑理論構建了式(20)中的相變能ψt(Ft)，其基本思路是在應變空間中基于相變本征應變及晶格對稱群構建相變路徑，相變遵循馬氏體相變的Bain 應變機制（原子在馬氏體相變過程中移動的距離最?。?，每個相變變體對應于相變路徑中能量穩(wěn)定的節(jié)點，以相變路徑為基礎進一步構建相變能函數 ψt，使相變沿著最接近本征應變的方向進行時能量最低。反復加卸載過程中變體數目增加，表現為應變空間中相變路徑的擴展，并不增加序參量數目，避免了經典相場模型中由于序參量數目的增加導致自由能函數構建復雜性增加的問題。

對于金屬材料在爆炸與沖擊等強動載條件下的馬氏體相變問題，反復加卸載下的馬氏體相變與塑性變形互相耦合形成大量新的變體，采用基于相變路徑理論的相場模型研究此類加載條件的微結構演化問題將具有良好的適用性。

3.2 相變路徑理論研究

相變路徑描述了馬氏體相變對稱性相關的變體形成基本信息，相場模型的能量函數需體現相變路徑機制，相變路徑和相場模型能量函數的構建是馬氏體相變相場模型研究中持續(xù)關注的基礎問題。相變路徑和能量函數的準確建立，可合理描述相變勢壘、界面能及相變的彈性松弛過程，這些因素對于馬氏體相變的成核至關重要[76]，而成核過程決定了后續(xù)微結構演化的基本模式。

馬氏體相變不僅有馬氏體晶胞的結構變化，還存在晶胞內部的挪動（Shuffling），相場模型需構建包含這些變形自由度的能量曲面。例如：近似情況下考慮6 維的應變空間，加上3 個Shuffle 自由度，則能量函數描述的是這9 維空間中的能量曲面，相變路徑可近似為連接該能量曲面中代表不同亞穩(wěn)相極值點的一維曲線，曲線路徑一般不是9 維空間中的直線，而是曲線甚至糾纏[77]。第一性原理計算可為相變路徑的確定及能量函數的構建提供有價值的信息[78–80]。

Gao 等[81]采用群論描述晶胞的結構變化及晶胞內的原子挪動導致的晶格對稱性的變化，進而定義了變形變體的概念，用于預測結構相變的相變路徑。該理論已應用于不同結構相變的相變路徑描述[82–83]及孿晶變形[84–85]，在相場模型的理論建模中具有較好的應用前景。

4 馬氏體相變與塑性變形耦合的相場模型研究

相變和塑性變形的相互作用是材料變形研究中的基礎問題[86]。馬氏體相變早期共格結構的相界面存在較大的應變能，會通過塑性變形過程進行應力松弛[87]；另一方面，相變與塑性變形都是能量耗散機制，兩個物理過程相互影響、相互競爭，塑性變形會改變相變序列[88]、誘導新的相變變體[70]、降低相變壓力[89]等，使微結構形成模式發(fā)生改變。相場模型基于能量的表述形式具有良好的擴展性，建模中易于考慮塑性變形機制，包括孿晶變形[12]、位錯滑移[13]、剪切帶[90]等塑性變形機制。相變與塑性的耦合是近年來馬氏體相變相場模擬研究中的熱點[50,87,90–96]。

塑性變形機制通過在應變中考慮塑性應變引入。在大變形情況下，采用式(12)的變形梯度分解考慮塑性變形。在小變形近似下，式(12)可近似為

相變與塑性的耦合主要關注塑性變形對相變成核、變體結構演化等的影響。例如：Levitas 等[50]建立了一種尺度無關的、耦合位錯帶與馬氏體相變的相場模型[50]。圖4 為雙晶結構在壓剪條件下馬氏體相變序參量分布的模擬結果，其中：圖4(a)為未考慮塑性變形時，壓剪作用下晶界處的應力集中導致高壓相的成核及變體生長；圖4(b)則考慮了兩個位錯滑移帶，由于位錯運動速度高于相變變體生長速度，因此滑移面的位錯運動松弛了應力狀態(tài)，改變了成核分布，另外滑移帶的存在阻礙了變體的生長。

圖4 塑性變形對雙晶結構在壓剪作用下高壓相成核生長的影響[50]Fig. 4 Effect of dislocation band on nucleation and evolution of the high pressure phase of bicrystal under compression and shear[50]

5 結 束 語

相場模型是研究馬氏體相變的重要理論計算方法，可揭示馬氏體相變不同微結構演化模式的形成機理及影響因素，在材料科學及工程應用中發(fā)揮著積極作用。相變機制和相變路徑是相場模型建立的基礎，第一性原理計算可提供有價值的信息，這對提升相場模型的預測能力具有重要意義。另外，通過在能量函數中考慮塑性變形、斷裂等物理過程，采用相場模型可研究馬氏體相變與其他物理過程的耦合機制，相關研究是馬氏體相變相場模型的重要發(fā)展方向。