功能梯度熱障涂層熱應力的有限體積法研究

龔京風，張文平，宣領寬，明平劍，韓春旭

（1.武漢理工大學能源與動力工程學院，湖北武漢430063；2.哈爾濱工程大學動力與能源工程學院，黑龍江哈爾濱150001）

功能梯度熱障涂層熱應力的有限體積法研究

龔京風1，2，張文平2，宣領寬2，明平劍2，韓春旭2

（1.武漢理工大學能源與動力工程學院，湖北武漢430063；2.哈爾濱工程大學動力與能源工程學院，黑龍江哈爾濱150001）

為了研究功能梯度熱障涂層（TBC）熱應力問題，改進了交錯格點型有限體積方法（SCV-FVM）。運用交錯網格技術，待解變量定義在單元節點上，物性參數定義在單元中心上。計算功能梯度圓盤的應力場，結果與理論解吻合良好，驗證方法的正確性。通過與有限元結果及高階功能梯度理論（HOTFGM）結果對比發現，SCV-FVM可以有效避免物性參數變化引起的人工數值不連續。研究功能梯度TBC熱應力問題，分析涂層厚度變化對熱應力場的影響，通過與ANSYS計算結果對比驗證數值模擬的正確性。

熱障涂層；功能梯度；有限體積法；交錯網格；應力場

在內燃機部件表面噴涂熱障涂層（TBC）是發展低散熱發動機的重要手段。TBC通常由粘結底層和陶瓷頂層組成，由于材料屬性的突變，在分層界面上產生應力集中。采用功能梯度材料（FGM），使材料屬性逐漸變化，可以有效緩解應力集中問題，同時增強材料粘結強度。準確預測功能梯度TBC熱力性能是設計優化涂層結構和材料分布的前提。

目前，關于功能梯度TBC熱力性能的研究主要基于試驗［1-5］和數值模擬［6-11］。其中，數值模擬主要利用有限元商用軟件，如ADINA［8］、ANSYS［7，10-11］。由于有限元方法（FEM）需要計算并存儲大型剛度矩陣，計算量和內存消耗較大［12］。另外，TBC屬于非均勻材料，對物性參數的不合理處理可能產生人工數值不連續問題［13-15］。為了克服這一問題，Santare和Lambros［16］將梯度單元引入傳統FEM中，而Aboudi等人［17］提出高階功能梯度理論（HOTFGM）。2種方法及其發展形式均能夠有效減輕不連續問題，但不能完全避免［10，13-20］。

本文旨在發展一種能有效預測功能梯度TBC熱力性能的數值方法。格點型有限體積法（CVFVM）兼有FVM的守恒性和FEM的靈活性，是FEM的一種有效補充方法。文獻［12］將交錯網格技術引入CV-FVM提出SCV-FVM，采用時間推進法解離散的平衡方程，因此該方法僅適用于動力學問題。本文進一步發展SCV-FVM用于研究穩態熱應力問題，利用隱格式耦合法求解離散的熱彈性方程。通過數值算例驗證發展的SCV-FVM預測熱應力場的準確性和對人工不連續問題的有效性，進而研究功能梯度TBC熱應力問題。

1 熱應力控制方程

考慮二維各向同性線彈性材料穩態熱應力問題，基于計算域內任意控制體微元S給出熱應力基本方程及邊界條件。

1.1 熱傳導控制方程

不考慮內熱源，穩態熱傳導方程為

式中：k和T分別為熱傳導系數和溫度。考慮3種邊界條件

式中：TB為邊界LD上給定的溫度（Dirichlet邊界）；qB為邊界LN上給定的法向熱流密度（Neumann邊界）；hB和T∞分別為邊界LR上給定的對流換熱系數和環境溫度（Robin邊界）。nα（α＝x，y）為邊界單位外法線矢量n的分量。

1.2 平衡方程

忽略體積力，根據牛頓第二定律得到平衡方程：

式中：σ為應力張量。本構方程和幾何方程為

式中：σαβ，εαβ和uα分別為應力張量σ、應變張量ε、位移矢量u的分量。a為熱膨脹系數，Tr為參考溫度。G和λ為拉姆系數，Γ為熱系數，根據不同條件由楊氏模量E和泊松比μ計算得到：

平面應變

平面應力

考慮2種邊界條件

式中：uαB為邊界LD上給定的位移（Dirichlet邊界）；tαB為邊界LN上給定的法向應力（Neumann邊界）。簡支邊界可由上述2種邊界條件組合得到。

2 數值方法

本文基于CV-FVM發展適用于功能梯度TBC熱應力問題的數值方法。圍繞單元節點依次連接相鄰單元中心及邊長中點建立控制體，如圖1所示。單元由實線圍成，其節點用空心圓點表示。控制體由虛線圍成，其中心及邊長中點用實心圓點表示。

圖1 網格單元及控制體示意圖Fig.1 Sketch map of grids and control volumes

采用交錯網格技術，從而將物性參數的空間變化引入離散過程。待解變量定義在單元節點上，并假設變量在控制體內均勻分布。物性參數定義在單元中心并假設均勻分布，則物性參數在控制體內是變化的，如圖1所示。

2.1 熱傳導方程的離散

采用高斯公式將方程寫為式中：L為控制體的邊。采用FEM的思想，利用型函數離散方程

式中：nc為當前節點周圍的單元總數，ncni為第i個單元內的節點個數。三角形單元內，ncni＝3；四邊形單元內，ncni＝4。N為型函數。下標i代表第i個單元中心的變量值，下標ij代表第i個單元內第j個節點上的變量值。

對于邊界上的節點，相應的控制方程離散需要考慮邊界條件的影響。對于Dirichlet邊界上的節點，其溫度直接根據式給定。對于Neumann邊界和Robin邊界上的節點，將式（3）和（4）代入方程（13）得

2.2 平衡方程的離散

將式（6）和（7）代入方程（5）得利用高斯公式及型函數，方程（15）可離散寫為

式中：Liα為積分線Li在α方向上的投影長度。

對于Dirichlet邊界上的節點，其位移直接根據式（10）給定。對于Neumann邊界上的節點，平衡方程的積分利用式（11）得

離散方程（14）和（17）中的型函數積分及積分線長度的計算可參考文獻［12］。

2.3 數值求解

將僅與網格幾何形狀有關的型函數積分及積分線長度作為幾何常數處理，計算一次并存儲，減少計算量。離散的平衡方程（17）中涉及2個方向的位移，采用耦合解法同時求解，避免迭代，從而一次計算獲得整個熱應力場。采用Fortran語言編程實現數值模擬。

3 FGM應力問題

適用于功能梯度TBC熱應力問題的數值方法，必須能夠準確計算FGM應力場。為了評價SCVFVM對FGM應力問題的適用性，計算文獻［15］中的平面應變FGM圓盤應力場，幾何模型見圖2。圓盤內徑ri＝0.025 4 m，受到法向壓力pi＝68.95 MPa。圓盤外徑ro＝0.050 8 m，受到法向壓力po＝6.89 MPa。圓盤泊松比為常數μ＝0.22，楊氏模量沿徑向變化E＝E0r2，其中E0＝300 GPa，r為半徑。

圖2 FGM圓盤示意圖Fig.2 Sketch map of the FGM disk

選取圓盤的1／4作為計算域，網格劃分與文獻［15］相同。圖3為不同方法得到的計算結果，其中Q4代表4節點四邊形單元，Q8代表8節點四邊形單元。通過比較可以看出SCV-FVM結果與理論解吻合良好，且能夠有效避免應力場的不連續，而HOTFGM理論及傳統FEM不能。另外，SCV-FVM結果與理論解誤差在邊界附近相對較大，這是由插值誤差引起的。SCV-FVM以位移為求解變量，單元中心的應力根據式（6）和（7）計算，然后插值得到節點應力，從而引入誤差。

圖3 FGM圓盤計算結果Fig.3 Results of the FGM disk

4 功能梯度TBC熱應力問題

將本文發展的SCV-FVM用于研究含FGM層的平面應變TBC熱力性能，計算模型見圖4，材料屬性見表1。其中第3層（L3）為FGM層，其物性參數變化規律為

式中：P代表任意物性參數，Pa、Pb分別為陶瓷頂層L4和粘結層L2的物性參數。

表1 功能梯度TBC材料屬性Table 1 Material properties of the functionally graded TBC

圖4 功能梯度TBC計算模型示意圖Fig.4 Sketch map of the functionally graded TBC

TBC下表面溫度Td＝25℃，上表面溫度沿x方向變化Tu＝1 325cos20（|x-1|）+25℃，左右表面絕熱。參考溫度Tr＝0℃。TBC下表面簡支，其余自由。采用均勻四邊形單元劃分計算域，網格總數為5 000。

當金屬底層L1厚度為Δy1＝0.2 m，L2厚度為Δy2＝0.1 m，L3厚度為Δy3＝0.5 m，L4厚度為Δy4＝0.2 m時，熱應力場云圖如圖5所示。由于L1和L2交界面處物性參數突變，存在應力集中現象（圖5（c））。FGM層L3起到過渡作用，L2、L3及L3、L4交界面處沒有應力集中。SCV-FVM計算結果沒有出現應力不連續現象。比較SCV-FVM與ANSYS在y＝0.5 m處的計算結果（見圖6），兩者吻合良好，驗證本文方法的正確性。

圖5 TBC熱應力場（Δy1＝0.2 m，Δy2＝0.1 m，Δy3＝0.5 m，Δy4＝0.2 m）Fig.5 The thermoelastic field of the TBC（Δy1＝0.2 m，Δy2＝0.1 m，Δy3＝0.5 m，Δy4＝0.2 m）

保持Δy2＝0.1 m，Δy4＝0.2 m，Δy1+Δy3＝0.7 m，改變Δy3研究FGM層厚度對TBC熱應力場的影響，見圖6。隨著Δy3的增加，位移ux及溫度T增大。合理選擇FGM層厚度，能夠有效減小應力幅值，使應力分布更光滑，如Δy3＝0.5 m和Δy3＝0.4 m。不合理的FGM層厚度不但不能改善應力分布，甚至使應力幅值比無FGM層時大得多，如Δy3＝0.3 m和Δy3＝0.2 m。

圖6 TBC中y＝0.5 m處的熱應力曲線Fig.6 The curves of variables at y＝0.5 m in the TBC

4 結束語

本文發展SCV-FVM，采用交錯網格技術，使得SCV-FVM適用于不均勻材料穩態熱應力分析。計算FGM圓盤應力場，結果與理論解吻合良好，通過與文獻對比發現SCV-FGM能夠有效避免不合理的應力不連續現象。基于SCV-FVM進一步研究功能梯度TBC熱應力問題，分析涂層厚度對TBC熱力性能的影響。計算得到的熱應力場不存人工數值不連續現象，結果與ANSYS吻合良好，表明SCV-FVM可以用作TBC設計階段的數值預測工具，作為FEM的有效補充方法。

［1］胡春枝，肖金生.陶瓷／金屬梯度熱障涂層的制備及性能評價［J］.海軍工程大學學報，2006，18（6）：75-78.HU Chunzhi，XIAO Jinsheng.Fabrication and characterization of ceramic／metal gradient thermal barrier coating［J］.Journal of Naval University of Engineering，2006，18（6）：75-78.

［2］原紅星.等離子噴涂法制備鋁基體梯度陶瓷涂層的組織及性能研究［D］.西安：西安理工大學，2009：21-45.YUAN Hongxing.Structure and performance of gradient ceramic coatings on aluminum substrate prepared by plasma spray［D］.Xi’an：Xi’an University of Technology，2009：21-45.

［3］HEJWOWSKI T，WERONSKI A.The effect of thermal barrier coatings on diesel engine performance［J］.Vacuum，2002，65：427-432.

［4］PARLAK A，SAHIN B，YASAR H.Performance optimisation of an irreversible dual cycle with respect to pressure ratio and temperature ratio-experimental results of a ceramic coated IDI diesel engine［J］.Energy Conversion and Management，2004，45：1219-1232.

［5］PARLAK A，YASAR H，ELDOGAN O.The effect of thermal barrier coating on a turbo-charged diesel engine performance and exergy potential of the exhaust gas［J］.Energy Conversion and Management，2005，46：489-499.

［6］劉杰，肖金生，覃峰，等.陶瓷／金屬梯度熱障涂層圓筒的傳熱與熱應力有限差分分析［J］.武漢理工大學學報：交通科學與工程版，2002，26（3）：379-382.LIU Jie，XIAO Jinsheng，QIN Feng，et al.Finite difference analysis for heat transfer and thermal stress of ceramic／metal gradient thermal barrier coating cylinder［J］.Journal of Wuhan University of Technology：Transportation Science and Engineering，2002，26（3）：379-382.

［7］劉春曉.催化／熱障功能梯度涂層的建模、仿真與實驗研究［D］.武漢：武漢理工大學，2003：30-47.LIU Chunxiao.Modeling simulation and experimental study of catalysis／thermal barrier coating functionally gradient coatings［D］.Wuhan：Wuhan University of Technology，2003：30-47.

［8］王素，倪春陽，朱心雄.功能梯度材料活塞三維溫度場分析［J］.北京航空航天大學學報，2005，31（12）：1299-1302.WANG Su，NI Chunyang，ZHU Xinxiong.Temperature field analysis for heterogeneous material piston［J］.Journal of Beijing University of Aeronautics and Astronautics，2005，31（12）：1299-1302.

［9］杜雙松.梯度熱障涂層活塞的穩態熱分析研究［D］.合肥：合肥工業大學，2008：29-54.DU Shuangsong.Steady-state thermal analysis of functional gradient materials as thermal barrier coating of engine piston［D］.Hefei：Hefei University of Technology，2008：29-54.

［10］BANSAL Y，PINDERA M J.Efficient reformulation of the thermoelastic higher-order theory for functionally graded materials［J］.Journal of Thermal Stresses，2003，26：1055-1092.

［11］BUYUKKAYA E.Thermal analysis of functionally graded coating AlSi alloy and steel pistons［J］.Surface and Coatings Technology，2008，202：3856-3865.

［12］龔京風，明平劍，宣領寬，等.基于有限體積法求解FGM動態響應及固有特性［J］.華中科技大型學報：自然科學版，2013，41（5）：25-30.GONG Jingfeng，MING Pingjian，XUAN Lingkuan，et al.Solving dynamic performance and natural characteristics of FGM by finite volume method［J］.Journal of Huazhong U-niversity of Science＆Technology：Natural Science Edition，2013，41（5）：25-30.

［13］KIM J H，PAULINO G H.Isoparametric graded finite elements for nonhomogeneous isotropic and orthotropic materials［J］.ASME Journal of Applied Mechanics，2002，69：502-514.

［14］CHAROENSUK J，VESSAKOSOL P.Numerical solutions for functionally graded solids under thermal and mechanical loads using a high-order control volume finite element method［J］.Applied Thermal Engineering，2011，31：213-227.

［15］CAVALCANTE M A A，MARQUES S P C，PINDERA M J.Computational aspects of the parametric finite-volume theory for functionally graded materials［J］.Computational Materials Science，2008，44：422-438.

［16］SANTARE M H，LAMBROS J.Use of graded finite elements to model the behavior of nonhomogeneous materials［J］.ASME Journal of Applied Mechanics，2000，67：819-822.

［17］ABOUDI J，PINDERA M J，ARNOLD S M.Higher-order theory for functionally graded materials［J］.Composites Part B，1999，30：777-832.

［18］ZHONG Y，BANSAL Y，PINDERA M J.Efficient reformulation of the thermal higher-order theory for FGMs with locally variable thermal conductivity［J］.International Journal of Computational Engineering Science，2004，5：795-831.

［19］CAVALCANTE M A A，MARQUES S P C，PINDERA M J.Parametric formulation of the finite-volume theory for functionally graded materials—Part I：analysis［J］.ASME Journal of Applied Mechanics，2007，74：935-945.

［20］CAVALCANTE M A A，MARQUES S P C，PINDERA M J.Parametric formulation of the finite-volume theory for functionally graded materials—Part II：numerical results［J］.ASME Journal of Applied Mechanics，2007，74：946-957.

Research on the finite volume method for the thermoelastic stress analysis of thermal barrier coatings

GONG Jingfeng1，2，ZHANG Wenping2，XUAN Lingkuan2，MING Pingjian2，HAN Chunxun2
（1.School of Energy and Power Engineering，Wuhan University of Technology，Wuhan 430063，China；2.College of Power and Energy Engineering，Harbin Engineering University，Harbin 150001，China）

In order to analyze the thermoelastic performance of a thermal barrier coating with a functionally graded material（FGM）layer，the staggered grid finite volume method（SCV-FVM）has been developed.The staggered grid technique was adopted by defining variables at the vertex and material properties at the cell center.SCV-FVM was applied to predict the stress of an FGM disk.The results agree well with the theoretical solution which validates its accuracy.The comparisons between the results from SCV-FVM，the conventional finite element method and the higher order theory functionally graded method（HOTFGM）show that SCV-FVM can avoid the stress discontinuity caused by the change in the matter parameters.The thermoelastic performance of functionally graded TBC was investigated.The obtained results are consistent with the ANSYS results which validates the simulations.The effects of the thickness have been discussed.

thermal barrier coating；functionally graded；finite volume method；staggered grid；stress field

10.3969／j.issn.1006-7043.201306012

http：／／www.cnki.net／kcms／doi／10.3969／j.issn.1006-7043.201306012.html

O343.1

A

1006-7043（2014）06-0713-06

2013-06-03.網絡出版時間：2014-05-14 15：49：09.

中央高校基本科研業務費專項資金資助項目（HEUCF130302）.

龔京風（1986-），女，博士研究生；張文平（1956-），男，教授，博士生導師.

張文平，E-mail：zhangwenping＠hrbeu.edu.cn.