許慶彥，夏鵠翔

（清華大學材料學院先進成形制造教育部重點實驗室，北京100084）

0 引言

高溫合金渦輪葉片是現(xiàn)代航空發(fā)動機極為重要的熱端部件，其服役環(huán)境惡劣，需在高溫環(huán)境下承受高應力載荷，且要求保持組織穩(wěn)定，這對高溫合金的成分、組織及制備工藝均提出了嚴苛要求。常見的高溫合金有Fe基、Co基和Ni基3大體系[1]，其中鎳基高溫合金由于其綜合性能優(yōu)異，在航空發(fā)動機葉片上應用最為廣泛。

葉片的微觀組織對其最終性能有直接影響，為獲得具有明顯各向異性的單晶組織，在生產(chǎn)過程中一般采用定向凝固工藝。高溫合金葉片生產(chǎn)中常見的定向凝固方法有功率降低法[2]、高速凝固法[3]和液態(tài)金屬冷卻法[4]等，這些方法能建立垂直的溫度梯度，使晶粒沿熱流方向定向生長。

鎳基單晶高溫合金的鑄態(tài)組織主要包括枝晶組織和枝晶間組織[5]。枝晶組織為單一γ相，面心立方晶體；枝晶間組織主要由γ相和γ'相組成，γ'相為具有L12結(jié)構(gòu)的有序相。由于合金中C、B、Hf等晶界強化元素的存在，在實際成形中，往往還會形成MC型碳化物[6]。此外在定向凝固成形過程中，容易產(chǎn)生各類缺陷，如雀斑[7-8]、雜晶[9-10]、條帶晶[11]、小角度晶界[12-13]等。這些鑄造缺陷的形成與單晶高溫合金的成分、組織形態(tài)、定向凝固工藝密切相關(guān)。研究各類缺陷的形成機理，有助于優(yōu)化定向凝固工藝，減少乃至消除缺陷，提高葉片合格率。

目前針對各類缺陷機理的相關(guān)研究主要通過試驗和數(shù)值模擬2種方法進行。但由于鎳基單晶高溫合金價格昂貴，通過試驗研究耗費大，且試驗周期長，工序繁瑣，難以對工藝進行針對性優(yōu)化。而隨著計算機技術(shù)的發(fā)展，數(shù)值模擬技術(shù)在航空發(fā)動機葉片的成形研究中發(fā)展十分迅速，扮演著越來越重要的角色。早期數(shù)值模擬研究聚焦定向凝固過程中各類宏觀物理過程，如傳熱、流體流動、受力變形等，并開發(fā)了有限差分法、有限體積法和有限元法等計算方法。隨著各類表征技術(shù)的發(fā)展，研究者們對發(fā)動機葉片的宏微觀組織有了更進一步的認知。伴隨著計算能力的提高，數(shù)值模擬技術(shù)也不斷朝著模擬晶體宏微觀組織演變的方向發(fā)展。近30年來，分子動力學方法[14]、元胞自動機方法[15]以及相場法[16]在鎳基高溫合金組織演變數(shù)值模擬研究中發(fā)展迅猛。

本文回顧了航空發(fā)動機鎳基高溫合金葉片定向凝固過程數(shù)值模擬方面的研究工作。從定向凝固的物理過程出發(fā)，描述宏觀物理場如溫度場、流場和應力應變場的數(shù)學模型，并介紹元胞自動機法和相場法2種組織計算方法；重點介紹微觀組織數(shù)值模擬技術(shù)在鎳基高溫合金成形過程的具體應用；對比分析了不同組織模擬方法的優(yōu)缺點，闡述不同組織模擬方法的適用范圍；展望了航空發(fā)動機鎳基高溫合金葉片定向凝固過程數(shù)值模擬技術(shù)的未來發(fā)展方向。

1 宏微觀數(shù)理建模

1.1 定向凝固工藝的物理建模

葉片生產(chǎn)中主流的定向凝固方法如圖1所示[17]。從圖中可見，3類裝置均有1個加熱區(qū)，定向凝固開始之前鑄型被放置在布置了感應線圈的加熱區(qū)中，并向內(nèi)注入熔煉后的高溫合金金屬液。放置鑄型的平臺為激冷盤，目的是為了讓鑄型底部金屬液率先過冷形核。3類定向凝固方法最主要的區(qū)別在于其冷卻方式不同。功率降低法僅通過調(diào)節(jié)感應線圈的閉合來控制凝固過程，得到的溫度梯度較低；液態(tài)金屬冷卻法采用的裝置中冷區(qū)有大量液態(tài)金屬冷卻液，在定向凝固過程中，隨著激冷盤平臺緩慢向下抽拉，型殼與冷區(qū)金屬冷卻液接觸，通過對流換熱獲得較大的溫度梯度；高速凝固法的裝置與液態(tài)金屬冷卻法的相似，但冷區(qū)沒有液態(tài)金屬冷卻液，冷卻方式為向外輻射換熱，溫度梯度處于上述2種方法之間。采用不同定向凝固方法可實現(xiàn)的溫度梯度和生長速度范圍如圖2所示[18]。圖中箭頭表示改變工藝參數(shù)時，溫度梯度和生長速度的改變方向。

圖1 主要的定向凝固方法

圖2 采用不同定向凝固方法得到的溫度梯度和生長速度范圍

1.2 數(shù)學模型（溫度場、流場、應力場）

在定向凝固過程中，葉片各部分的溫度梯度是極為重要的變量，直接影響葉片微觀組織和殘余應力。為了準確描述并模擬定向凝固過程的宏觀溫度變化，可以通過傳熱方程[19]確保體系能量守恒

式中：ρ為材料密度；cp為比熱容；T為熱力學溫度；t為時間；λ為導熱系數(shù)；ΔH為結(jié)晶潛熱；f為凝固時網(wǎng)格固相率；Qr為與環(huán)境之間的熱流密度，在不同的定向凝固方法中，Qr的形式和計算方法不同。

在求解流體流動過程中，需對計算單元進行連續(xù)性方程和N-S方程求解，來保持系統(tǒng)質(zhì)量守恒和動量守恒

傳統(tǒng)的彈塑性本構(gòu)模型和基于滑移理論的晶體塑性有限元方法在高溫合金成形過程的應力應變模擬中均有廣泛應用。在晶體塑性有限元模型中將總變形梯度F分解為塑性剪切Fe和彈性變形Fp2部分

依據(jù)2階張量的標準分解法則，將速度梯度L張量分為對稱部分（變形率D）和偏對稱部分（旋轉(zhuǎn)率R）之和，二者分別代表晶格的變形與旋轉(zhuǎn)

1.3 形核與生長模型

目前在高溫合金組織模擬領(lǐng)域常用的研究方法有蒙特卡洛方法、分子動力學方法、元胞自動機方法以及相場法。蒙特卡洛方法隨機性太強，分子動力學方法研究尺度太小，所以二者的研究應用范圍受到了限制；而元胞自動機方法和相場法由于其優(yōu)勢明顯，發(fā)展迅速，是目前高溫合金組織模擬領(lǐng)域的主流方法。

1.3.1 元胞自動機方法

元胞自動機是一種應用比較廣泛的模型理論，由馮·諾依曼建立，經(jīng)多位數(shù)學家和物理學家補充后迅速發(fā)展，目前被廣泛應用于通信、交通、材料學等多學科和領(lǐng)域。Rappaz[20]首先將元胞自動機方法引入凝固過程組織生長模擬中，按照特定的局部規(guī)則，在離散的時間維度上求解晶體組織變化。

定向凝固過程組織模擬主要包括晶粒形核與晶體生長2部分。晶粒形核在算法上有2種實現(xiàn)方式，一種是瞬時形核，直接在特定位置部署晶體核心。在定向凝固模擬中，一般會在模擬的底端設(shè)置1層帶有噪聲的晶體籽晶，用來模擬遇到激冷盤時葉片底部初始過冷形核。另一種是連續(xù)形核，主要發(fā)生在過冷區(qū)域。當局部位置滿足一定的形核條件時，將自發(fā)形核。連續(xù)形核過程常用Gauss模型[21]描述為

式中：n為形核密度；nmax為最大形核密度；ΔT為當前單元的過冷度；ΔTσ為過冷度標準差；ΔTN為平均形核過冷度。

利用Kurz-Giovanola-Trivedi（KGT）解析模型[22]，可計算固液界面推進速度，描述晶體生長過程

式中：vn為界面法向推進速度；α和β為速度動力學系數(shù)。

1.3.2 相場法

經(jīng)過30多年的發(fā)展，相場法已成為模擬相變過程的主流研究方法。與其他隨機性模擬方法不同，相場法是以金茲堡朗道2級相變理論為基礎(chǔ)的確定性數(shù)值模擬方法，具備堅實的物理理論基礎(chǔ)。Kobayas?hi[23]首先提出了描述純物質(zhì)凝固的相場模型，發(fā)展出了多元多相體系相場算法。針對計算量大、計算空間小等問題，各國學者通過改進計算方式，實現(xiàn)了目前主流的非均勻網(wǎng)格加速計算和GPU加速計算，極大地提高了計算效率。相場模型中金茲堡朗道類型的自由能泛函F表達式為

式中：φα為α相的相場變量為α相的各元素百分比；Ω為計算域；f為自由能密度函數(shù)fintf（{φα}），由界面自由能密度函數(shù)和化學自由能密度函數(shù)fchem（{φα}，組成

式中：n為相的個數(shù)；σαβ為α相與β相之間界面能；為α相的體積自由能密度函數(shù)，與其當前成分相關(guān)。

由體系自由能泛函對相求變分可以得到相場變量隨時間的變化規(guī)律

為求解化學自由能，還需對體系吉布斯自由能G進行分析

式中：Nα為相的個數(shù)；u為合金元素種數(shù)。

2 定向凝固過程的宏觀物理場模擬

溫度場模擬是高溫合金定向凝固過程中其他物理場模擬的基礎(chǔ)，應力應變模擬、組織演變、缺陷形成等數(shù)值模擬分析手段的合理應用都需要依賴準確的溫度場模擬結(jié)果。

梁作儉[24]、崔鍇[25]、于靖[26]和王潤楠等[27]分別對各種簡單形狀和實際復雜形狀單晶葉片的定向凝固過程進行模擬，研究了凝固過程葉片溫度分布、縱向溫度梯度和冷卻速度的變化規(guī)律。某復雜葉片的定向凝固過程溫度場模擬結(jié)果如圖3所示。從圖中可見，抽拉速度越大，生長界面曲率和等溫線曲率越大，縱向溫度梯度越小，越容易產(chǎn)生結(jié)晶缺陷。其結(jié)晶缺陷預測結(jié)果與試驗結(jié)果大體一致，模擬模型準確性較高。

圖3 某復雜單晶葉片定向凝固過程中的溫度場演變

為了克服單晶葉片緣部雜晶缺陷，確保凝固過程葉片葉身段溫度分布平穩(wěn)、糊狀區(qū)連續(xù)生長，于靖[26]、潘冬等[28]模擬研究了定向凝固過程中變抽拉速度工藝對葉片溫度及糊狀區(qū)分布的影響，結(jié)果如圖4所示。從圖中可見，相比于定拉速工藝，在變拉速工藝條件下，葉片溫度分布良好，縱向溫度梯度均勻，糊狀區(qū)穩(wěn)定且水平，孤立過冷區(qū)大幅減少，有利于抑制雜晶的形成。

圖4 某復雜葉片在變拉速工藝下的模擬結(jié)果

許自霖[29]采用商用有限元軟件ProCAST模擬了定向凝固過程中單晶葉片的溫度場，然后將溫度場結(jié)果導入通用有限元力學軟件Abaqus中，對葉片受力應變過程進行分析，其應力分布模擬結(jié)果如圖5所示。從圖中可見，在圓柱坐標系中，葉片的σmises應力與Z方向σzz應力變化趨勢大致相同。等效塑性應變出現(xiàn)在排氣邊與上部緣板交匯處，最大值超過2%。該位置與工廠實際生產(chǎn)中葉片裂紋的位置相符，驗證了模型的準確性。

圖5 渦輪葉片定向凝固過程的σmises、σzz應力分布及實際葉片裂紋

李忠林[30]針對單晶空心渦輪葉片進行了應力應變分析，其等效塑性應變分布的模擬結(jié)果如圖6所示。從圖中可見，在葉片排氣邊上下兩側(cè)、型芯內(nèi)部截面突變處、緣板根部以及緣板的引晶條連接處均出現(xiàn)了較大的塑性變形。模擬得到的較大的塑性變形區(qū)域與實際生產(chǎn)中葉片的再結(jié)晶位置大體吻合。

圖6 單晶空心渦輪葉片定向凝固1h時塑性變形分布模擬結(jié)果（800℃）

3 葉片定向凝固組織模擬

3.1 晶粒結(jié)構(gòu)模擬

在凝固過程組織生長模擬中應用最為廣泛的2種方法——元胞自動機方法和相場法中，相場法計算量大，對計算機內(nèi)存要求高，現(xiàn)階段還無法實現(xiàn)宏觀尺度模擬；元胞自動機方法由于其在計算量和計算尺度上的優(yōu)勢，在凝固過程宏觀組織生長模擬方面的應用無可替代。

于靖[26]使用元胞自動機有限差分數(shù)值模擬對鎳基高溫合金葉片定向凝固過程的晶粒組織進行了研究，模擬結(jié)果如圖7所示。圖中不同的顏色代表不同的晶粒。從圖中可見，在定向凝固過程中，底部晶粒激冷形核后以柱狀晶的方式生長，隨著生長過程的進行，僅余最優(yōu)生長方向與縱軸結(jié)晶的晶粒得以繼續(xù)生長。

圖7 定向凝固葉片組織模擬

張航等[31]在此基礎(chǔ)上對鎳基高溫合金選晶器定向凝固過程的晶粒組織進行了詳細的研究，模擬結(jié)果如圖8所示。從圖中可見，隨著凝固分數(shù)fs增大，固/液界面傾斜程度增大，這種生長傾斜是溫度場的傾斜對晶體競爭生長和選晶過程產(chǎn)生不利影響的結(jié)果。當溫度場發(fā)生傾斜，部分晶粒橫向生長速度加快，造成晶粒大小不均，嚴重影響晶粒密度分布。同時，可能造成部分粗大晶粒直接進入到選晶器的過渡段，使螺旋段作用失效，從而對晶粒取向分布產(chǎn)生不利影響。

圖8 引晶段溫度場與凝固組織模擬結(jié)果

3.2 枝晶生長模擬

相場法是一種基于金茲堡朗道2級相變原理以及吉布斯自由能體系的確定性組織模擬方法。相比于其他隨機性模擬方法（如元胞自動機方法和蒙特卡洛方法），在物理理論基礎(chǔ)上更為扎實，結(jié)果相對更為可信。元胞自動機方法[32]與相場方法[33]的模擬對比結(jié)果如圖9所示。從圖中可見，相場方法采用彌散界面模型代替了傳統(tǒng)的尖銳界面處理方式來描述界面狀態(tài)，這與真實材料一致，在真實材料中，相界面或晶界并不是嚴格的零厚度界面，而是具有一定厚度的邊界層。

圖9 元胞自動機方法和相場方法的模擬結(jié)果比較

相場法在鎳基高溫合金定向凝固過程中組織形態(tài)變化的數(shù)值模擬結(jié)果已經(jīng)與實際試驗結(jié)果十分接近，如圖10所示[34]。楊聰?shù)萚35]采用相場法系統(tǒng)定量地研究合金中每種成分的偏析行為，如圖11所示。其研究還對多元鎳基高溫合金定向凝固過程中初生枝晶間距進行了模擬。結(jié)果表明，當晶粒傾角θ＜20°時，一次枝晶間距與θ=0°時沒有明顯區(qū)別。當θ＞20°時，一次枝晶間距隨θ的增大而顯著增大。

圖10 相場模擬的組織形態(tài)與試驗結(jié)果

圖11 在50 μm/s的抽拉速度和5 K/mm的溫度梯度、穩(wěn)態(tài)下的枝晶形態(tài)和不同成分的溶質(zhì)分布

3.3 共晶生長模擬

了解鎳基單晶高溫合金凝固過程中的凝固行為對于確定后續(xù)熱處理工藝至關(guān)重要，而試驗研究卻發(fā)現(xiàn)了相互矛盾的結(jié)論。一方面，枝晶間細小的γ/γ'結(jié)構(gòu)表明了共晶反應的發(fā)生；另一方面，粗大的γ'沉淀表明了包晶反應的發(fā)生。楊聰?shù)萚36]通過耦合了熱力學數(shù)據(jù)庫的多元多相場模型對其生長過程進行模擬，并在γ'相附近觀察到γ相的重熔，表明γ/γ'結(jié)構(gòu)是共晶反應的結(jié)果，如圖12所示。圖中紅色代表液態(tài)，藍色代表γ相，黃色代表γ'相。

圖12 γ′相附近觀察到γ相的重熔

3.4 宏-微觀多尺度耦合模擬

微觀組織模擬計算量大，計算時間長，現(xiàn)階段還無法對整個宏觀尺度部件進行全局微觀組織模擬計算。在此局限下，宏-微觀耦合模擬應運而生，逐漸成為一種模擬趨勢。在定向凝固宏-微觀耦合模擬中，常使用宏觀模擬模型計算整個計算域流場、溫度場等宏觀物理場，截取局部區(qū)域的溫度場等物理場作為物理邊界條件，導入至微觀組織模擬模型進行計算，最后通過組合宏觀和微觀計算結(jié)果，得出整個零件的微觀組織發(fā)展規(guī)律。

張航[37]耦合了宏觀溫度計算結(jié)果及介微觀晶體生長模型，模擬結(jié)果如圖13所示。溫度場計算結(jié)果為介觀和微觀晶粒組織模擬提供基本溫度分布信息。介觀尺度晶粒競爭生長模擬結(jié)果如圖13（b）所示，該結(jié)果為微觀枝晶生長模擬提供晶粒取向及結(jié)構(gòu)分布等信息。耦合考慮宏觀溫度分布及介觀晶粒生長模擬，最終實現(xiàn)微觀枝晶3維生長過程的模擬。

圖13 引晶段多尺度耦合模擬計算

閆學偉等[3]通過模擬和試驗研究了液態(tài)金屬冷卻定向凝固法制備單晶渦輪葉片過程。建立了耦合溫度場、晶粒生長和溶質(zhì)擴散的多尺度模型，模擬了大尺度2維枝晶在整個截面上的分布，模擬和試驗結(jié)果具有良好的一致性，其結(jié)果如圖14所示。

圖14 枝晶形態(tài)的模擬和試驗結(jié)果

4 總結(jié)與展望

本文回顧了鎳基高溫合金航空發(fā)動機葉片成形過程數(shù)值模擬研究進展情況。經(jīng)過幾十年的發(fā)展，目前已經(jīng)成功地實現(xiàn)了宏觀傳熱、流體流動、應力應變以及宏微觀組織生長模擬，數(shù)值模擬技術(shù)已經(jīng)成為工藝優(yōu)化必不可少的關(guān)鍵手段。近些年來，高溫合金定向凝固過程宏微觀組織跨尺度模擬成為新的研究熱點，多種不同尺度的數(shù)值模擬技術(shù)均得到迅速發(fā)展和廣泛應用。在未來，高溫合金航空發(fā)動機葉片定向凝固數(shù)值模擬研究還有望從如下幾方面進一步加強：

（1）耦合熱力學數(shù)據(jù)庫，建立具有各向異性的鎳基單晶高溫合金的力學模型；

（2）開發(fā)多物理場耦合、跨尺度組織計算方法，進一步研究宏微觀組織演變過程；

（3）利用數(shù)值模擬技術(shù)輔助研究小角度晶界等缺陷機理，建立小角度晶界等缺陷數(shù)值模擬計算模型，并構(gòu)建相應的工藝優(yōu)化方案；

（4）開發(fā)構(gòu)建高效的并行計算方法以及內(nèi)存分配方案，實現(xiàn)3維大尺度相場組織模擬。