王家鋆 敖勇 鄧志成

摘要：為研究不同網(wǎng)格劃分方案對(duì)汽輪機(jī)內(nèi)缸應(yīng)力場計(jì)算結(jié)果的影響，對(duì)某汽輪機(jī)中壓內(nèi)缸進(jìn)行建模，選擇熱固耦合的有限元法采用4種不同的網(wǎng)格劃分方案進(jìn)行計(jì)算，得到額定負(fù)荷工況下典型應(yīng)力集中部位的等效應(yīng)力。根據(jù)計(jì)算結(jié)果，分析網(wǎng)格整體控制和局部細(xì)化過程的不同設(shè)置參數(shù)對(duì)整個(gè)求解過程的影響。綜合考慮求解成本和計(jì)算精度，最終確定一種理想的四面體網(wǎng)格劃分方案，該方案可在結(jié)果準(zhǔn)確合理的前提下提高計(jì)算效率。

關(guān)鍵詞：汽輪機(jī)內(nèi)缸;應(yīng)力集中;網(wǎng)格劃分;優(yōu)化;有限元

中圖分類號(hào)：TK263.1;TB115.1文獻(xiàn)標(biāo)志碼：B

0 引言

有限元法是求解復(fù)雜力學(xué)問題的有效數(shù)值分析方法之一。隨著計(jì)算機(jī)技術(shù)的不斷發(fā)展，在針對(duì)汽輪機(jī)零部件結(jié)構(gòu)強(qiáng)度的計(jì)算方面，有限元法得到廣泛應(yīng)用。在有限元強(qiáng)度分析中，只要原問題的數(shù)學(xué)模型正確、算法穩(wěn)定可靠，隨著單元數(shù)目增加，有限元解與精確解的近似程度就被不斷改善，并最終逼近原數(shù)學(xué)模型的精確解。因此，有針對(duì)性地選擇網(wǎng)格劃分策略尤為重要。

對(duì)幾何結(jié)構(gòu)復(fù)雜、工作環(huán)境惡劣的汽輪機(jī)中壓內(nèi)缸進(jìn)行網(wǎng)格劃分時(shí)，若選用六面體單元，則需要在保證計(jì)算精度的前提下對(duì)模型進(jìn)行大量的簡化和分割;四面體單元具有劃分方便、對(duì)模型邊界擬合能力較強(qiáng)的特點(diǎn)。文獻(xiàn)[3]和[4]均對(duì)汽輪機(jī)中壓內(nèi)缸進(jìn)行完整的結(jié)構(gòu)強(qiáng)度研究，并采用四面體單元進(jìn)行網(wǎng)格劃分，但網(wǎng)格劃分分別采用自動(dòng)劃分和設(shè)置全局尺寸的方法，對(duì)強(qiáng)度薄弱部位未作區(qū)分且得到的網(wǎng)格數(shù)量較多。針對(duì)汽輪機(jī)中壓內(nèi)缸的具體結(jié)構(gòu)，研究網(wǎng)格劃分對(duì)應(yīng)力計(jì)算結(jié)果影響的文章較少。

本文以某電廠運(yùn)行中的汽輪機(jī)中壓內(nèi)缸為研究對(duì)象，應(yīng)用有限元計(jì)算分析軟件worLbench對(duì)額定負(fù)荷工況下的部件強(qiáng)度進(jìn)行有限元分析，并采取不同的四面體網(wǎng)格劃分方案進(jìn)行應(yīng)力計(jì)算和分析對(duì)比，從而總結(jié)高效準(zhǔn)確的汽輪機(jī)內(nèi)缸網(wǎng)格劃分方法，以期為汽輪機(jī)啟動(dòng)和停機(jī)等歷時(shí)較長的瞬態(tài)過程分析提供參考。

1 基本理論

1.1 溫度場與應(yīng)力場數(shù)學(xué)模型

對(duì)于中壓內(nèi)缸在額定負(fù)荷工況下溫度場和應(yīng)力場的求解問題，其穩(wěn)態(tài)溫度場的場變量在由直角坐標(biāo)系中應(yīng)滿足拉普拉斯方程，即

利用加權(quán)余量的伽遼金法可以得到確定各節(jié)點(diǎn)溫度φ_i的有限元求解方程，即

Kφ=P （3）

式中：K為熱傳導(dǎo)矩陣;φ為節(jié)點(diǎn)溫度矩陣;P為溫度載荷矩陣。這是一組以節(jié)點(diǎn)溫度為獨(dú)立變量的線性方程組。

求得彈性體的溫度場后，進(jìn)一步通過節(jié)點(diǎn)位移求出彈性體各部分的熱應(yīng)力，其有限元求解方程為

K'a=P（4）

P'=P_T+P_ε0（5）

式中：K'為結(jié)構(gòu)整體剛度矩陣;a為單元節(jié)點(diǎn)位移矩陣;P為結(jié)構(gòu)節(jié)點(diǎn)載荷矩陣;P_ε0和P_T分別為表面載荷和溫度應(yīng)變引起的載荷矩陣。

1.2 局部網(wǎng)格細(xì)化依據(jù)

圣維南原理是彈性力學(xué)中的經(jīng)典理論，該原理表明：如果物體在局部邊界表面所受的表面力是平衡力系，則這個(gè)平衡表面力只在受力點(diǎn)附近產(chǎn)生顯著的應(yīng)力，在受力點(diǎn)遠(yuǎn)處其應(yīng)力可忽略不計(jì)。在網(wǎng)格劃分中也存在著類似的規(guī)律：在距離所關(guān)注部位的遠(yuǎn)處，其網(wǎng)格劃分的疏密程度對(duì)該點(diǎn)的應(yīng)力造成的影響可以忽略不計(jì)。因此，在具體計(jì)算中可以對(duì)模型的不同部位分別進(jìn)行網(wǎng)格細(xì)化設(shè)置而互不產(chǎn)生影響。

2 網(wǎng)格劃分目標(biāo)方案計(jì)算分析

2.1 汽輪機(jī)中壓內(nèi)缸有限元模型

計(jì)算對(duì)象為某電廠亞臨界汽輪機(jī)的中壓內(nèi)缸，該結(jié)構(gòu)主要由圍繞轉(zhuǎn)子旋轉(zhuǎn)的殼體結(jié)構(gòu)、高窄法蘭、缸體與法蘭間的過渡結(jié)構(gòu)和進(jìn)汽管道等組成。內(nèi)缸由外缸的水平中分面支承，頂部和底部由定位銷導(dǎo)向，內(nèi)缸凸臺(tái)與外缸槽配合實(shí)現(xiàn)內(nèi)缸定位。左、右半缸結(jié)構(gòu)對(duì)稱，故取整個(gè)中壓內(nèi)缸的1/2作為研究對(duì)象建立三維有限元模型，計(jì)算傳熱系數(shù)（作為熱載荷）和內(nèi)缸表面不同區(qū)域處的蒸汽溫度。

進(jìn)行穩(wěn)態(tài)額定負(fù)荷工況溫度場計(jì)算，中壓內(nèi)缸額定負(fù)荷工況溫度場云圖見圖1。將該溫度場作為應(yīng)力場計(jì)算的溫度邊界條件導(dǎo)人靜力分析模塊;將連接上、下半缸所需的螺栓預(yù)緊力，由于進(jìn)氣管在上、下半缸非對(duì)稱布置所引起的管道應(yīng)力，以及在內(nèi)缸不同區(qū)域處的蒸汽壓力作為等效應(yīng)力場計(jì)算的應(yīng)力邊界條件。

對(duì)中壓內(nèi)缸進(jìn)行較為粗糙的自動(dòng)網(wǎng)格劃分。根據(jù)初步計(jì)算結(jié)果，在該內(nèi)缸表面選取4處典型的應(yīng)力集中部位（見圖2），分別為高壓側(cè)內(nèi)缸表面圓角（A）、進(jìn)汽管內(nèi)表面上部圓角（B）、進(jìn)汽管外表面圓角（C）和法蘭凸臺(tái)上部圓角（D）。

2.2 目標(biāo)方案的網(wǎng)格無關(guān)性驗(yàn)證

有限元網(wǎng)格劃分直接影響后續(xù)數(shù)值計(jì)算分析結(jié)果的精確性，是進(jìn)行數(shù)值模擬分析至關(guān)重要的一步。在本算例中，四面體網(wǎng)格的生成采用Workbench軟件自帶的Patch conforming算法，即協(xié)調(diào)分片算法。對(duì)單元進(jìn)行不斷細(xì)化，直至前后兩次的等效應(yīng)力計(jì)算結(jié)果的相對(duì)誤差小于0.2%，且應(yīng)力集中部位的最大單元平均應(yīng)力與節(jié)點(diǎn)平均應(yīng)力差別不大，相對(duì)誤差不超過5%，則認(rèn)為有限元計(jì)算結(jié)果穩(wěn)定。以此時(shí)的網(wǎng)格劃分方案為目標(biāo)方案，將其應(yīng)力場計(jì)算結(jié)果作為其他方案的參考值，考察相對(duì)誤差大小。目標(biāo)方案的中壓內(nèi)缸應(yīng)力集中部位的有限元網(wǎng)格見圖3。

在結(jié)構(gòu)場后處理分析中，VON Mises應(yīng)力是一項(xiàng)重要的反映指標(biāo)，其在Workbench中有多種表現(xiàn)形式，最常見的是節(jié)點(diǎn)平均應(yīng)力，即在對(duì)所有單元進(jìn)行計(jì)算得到其節(jié)點(diǎn)應(yīng)力后，再對(duì)共享節(jié)點(diǎn)的應(yīng)力進(jìn)行平均，即為該點(diǎn)的應(yīng)力。單元平均應(yīng)力是指在對(duì)節(jié)點(diǎn)應(yīng)力進(jìn)行平均計(jì)算后，對(duì)單元內(nèi)所有的節(jié)點(diǎn)應(yīng)力再一次進(jìn)行平均，得到單元內(nèi)部的平均應(yīng)力。隨著求解精度的不斷提高，單元平均應(yīng)力會(huì)逐步趨近于節(jié)點(diǎn)平均應(yīng)力，因此兩者的近似水平可作為判斷計(jì)算結(jié)果是否穩(wěn)定的標(biāo)準(zhǔn)。以部位A為例，其節(jié)點(diǎn)應(yīng)力與單元應(yīng)力的對(duì)比見圖4。部位A的最大節(jié)點(diǎn)平均應(yīng)力和最大單元平均應(yīng)力分別為62.29和61.54MPa，且均出現(xiàn)在部位A的相近位置處。

3 網(wǎng)格劃分改進(jìn)方案的確定

在Workbench軟件中，難以直接對(duì)具有復(fù)雜結(jié)構(gòu)的幾何體進(jìn)行掃掠法網(wǎng)格劃分。在工程設(shè)計(jì)和計(jì)算中，一般選用Automatic方法自動(dòng)生成四面體網(wǎng)格。采用Automatic方法劃分網(wǎng)格時(shí)，較密的網(wǎng)格參數(shù)設(shè)置會(huì)極大地增加計(jì)算工作量和求解時(shí)間，而粗糙的網(wǎng)格設(shè)置又無法準(zhǔn)確捕捉零件在應(yīng)力集中部位的應(yīng)力分布，因此一般采用網(wǎng)格重新劃分的思路，即首先進(jìn)行整體網(wǎng)格控制，隨后在應(yīng)力集中部位進(jìn)行局部網(wǎng)格細(xì)化，從而達(dá)到兼顧計(jì)算成本和求解精度要求的目的。

整體網(wǎng)格控制主要包括默認(rèn)（Defaults）設(shè)置和網(wǎng)格尺寸函數(shù)設(shè)置。Defaults設(shè)置中確定物理場選項(xiàng)為結(jié)構(gòu)場。網(wǎng)格相關(guān)度（Relevance）選項(xiàng)數(shù)值從-100到+100，代表網(wǎng)格由疏到密，不同的值對(duì)應(yīng)不同的網(wǎng)格數(shù)和節(jié)點(diǎn)數(shù)，在劃分過程中可根據(jù)需要進(jìn)行調(diào)整。在目標(biāo)網(wǎng)格劃分方案的確定過程中，可通過不斷增大關(guān)聯(lián)度使得計(jì)算結(jié)果穩(wěn)定。在網(wǎng)格尺寸函數(shù)設(shè)置中，主要考慮高級(jí)尺寸函數(shù)的選取、關(guān)聯(lián)中心的設(shè)定和全局單元尺寸的調(diào)整等。這里用到的高級(jí)尺寸函數(shù)為自適應(yīng)網(wǎng)格函數(shù)和曲率網(wǎng)格函數(shù);關(guān)聯(lián)中心包括描述網(wǎng)格粗糙、中等和細(xì)化等3種模式，和網(wǎng)格相關(guān)度一起對(duì)網(wǎng)格整體尺寸產(chǎn)生影響;全局單元尺寸取各劃分方法下的默認(rèn)值。對(duì)于該汽輪機(jī)中壓內(nèi)缸而言，用到的局部網(wǎng)格細(xì)化方法主要包括插入局部尺寸控制（Sizing）選項(xiàng)和插入直接細(xì)化分割（Refinement）選項(xiàng)。Refinement直接基于細(xì)化倍數(shù)對(duì)初始單元邊界進(jìn)行不同程度的分割，網(wǎng)格間的過渡性較好、劃分方便。調(diào)用Sizing選項(xiàng)時(shí)可以根據(jù)所關(guān)注部位的圓角大小分別設(shè)定單元尺寸，靈活性較強(qiáng)。

通過以上討論并進(jìn)行反復(fù)試驗(yàn)，最終確定的網(wǎng)格劃分目標(biāo)方案和改進(jìn)方案見表1，其中R3表示采用細(xì)化等級(jí)為3的Reflnement選項(xiàng)，S8表示在該處插入Sizing選項(xiàng)控制局部單元大小且平均單元特征長度為8mm，以此類推。目標(biāo)方案是網(wǎng)格細(xì)化程度最高的計(jì)算方案，以此作為其他網(wǎng)格劃分方案得到的應(yīng)力計(jì)算結(jié)果的參考。

4 結(jié)果與分析

4.1 計(jì)算結(jié)果

網(wǎng)格劃分的目標(biāo)方案和其他4種改進(jìn)方案的計(jì)算結(jié)果對(duì)比見表2。根據(jù)表2數(shù)據(jù)計(jì)算得到的各方案最大等效應(yīng)力與目標(biāo)值的相對(duì)誤差見表3。由表3得到不同方案在各部位的最大等效應(yīng)力與目標(biāo)值相對(duì)誤差絕對(duì)值對(duì)比的柱形圖，見圖5。

圖5可以直觀地反應(yīng)各方案計(jì)算結(jié)果相對(duì)誤差絕對(duì)值在不同部位的分布情況。通過試算可知，以上各方案的最大等效應(yīng)力在同一部位的出現(xiàn)位置相近，可保證不同方案在同一部位的計(jì)算結(jié)果對(duì)比有意義。同時(shí)，各計(jì)算結(jié)果也均未超過材料在該溫度場下的屈服極限。

由表2中各方案網(wǎng)格劃分節(jié)點(diǎn)總數(shù)與對(duì)應(yīng)的求解時(shí)間，可以得到不同方案節(jié)點(diǎn)總數(shù)與求解時(shí)間的關(guān)系，見圖6。網(wǎng)格劃分的目標(biāo)方案、方案3和方案4在部位D的局部網(wǎng)格放大圖對(duì)比見圖7。

4.2 結(jié)果分析

從求解時(shí)間來看，目標(biāo)方案的求解耗時(shí)最長，達(dá)147.58s，計(jì)算成本較高;方案2和方案3的求解耗時(shí)較短，均未超過30.00s。

各方案在部位A、B的計(jì)算結(jié)果與目標(biāo)值差異較小，在部位C、D差異較大，最大達(dá)38.13%。原因在于部位C、D的圓角較小，而隨著圓角的減小，各方案的網(wǎng)格劃分策略和精密程度對(duì)該部位應(yīng)力計(jì)算結(jié)果的影響也越明顯。因此，對(duì)于結(jié)構(gòu)有圓角的部位，采用自動(dòng)網(wǎng)格劃分后若不再對(duì)其進(jìn)行局部細(xì)化，計(jì)算結(jié)果會(huì)有較大誤差，無法準(zhǔn)確反映該部位的真實(shí)應(yīng)力狀況?？傮w而言，方案3和4的相對(duì)誤差要小于方案l和2。

方案1在部位D的計(jì)算結(jié)果與目標(biāo)值的相對(duì)誤差（絕對(duì)值）遠(yuǎn)小于方案2和4。原因是該方案采用CurVaturc的高級(jí)尺寸函數(shù)控制方法。該方法由曲率法確定細(xì)化邊緣和曲面處的網(wǎng)格大小，在有曲率變化的地方網(wǎng)格會(huì)自動(dòng)加密。因此，對(duì)于具有小圓角結(jié)構(gòu)（如部位D）的應(yīng)力計(jì)算結(jié)果誤差較小。但是，該函數(shù)同時(shí)也在其他并不重要的曲率變化處自動(dòng)進(jìn)行網(wǎng)格加密，所以計(jì)算工作量和求解時(shí)間增加。

對(duì)于該模型的整體網(wǎng)格劃分來說，各方案所產(chǎn)生的節(jié)點(diǎn)總數(shù)與最終求解時(shí)間基本均呈線性關(guān)系，因此兩者均可以作為反映求解成本的參數(shù)。

從應(yīng)力結(jié)果和圖7可以看出，對(duì)于部位D而言，采用根據(jù)圓角大小設(shè)定單元平均尺寸的Sizing控制方法比采用Refinement直接進(jìn)行局部細(xì)化形成的網(wǎng)格更細(xì)密規(guī)整，因此方案3對(duì)細(xì)小圓角處的真實(shí)應(yīng)力狀況反映得更好。通過以上分析，綜合考慮有限元計(jì)算的時(shí)間成本和各應(yīng)力集中部位的求解精度，最終認(rèn)為方案3為最佳方案。

5 結(jié)束語

針對(duì)汽輪機(jī)內(nèi)缸的三維結(jié)構(gòu)模型，分析2種局部細(xì)化方法對(duì)應(yīng)力集中部位應(yīng)力計(jì)算的影響，總結(jié)汽輪機(jī)內(nèi)缸在穩(wěn)態(tài)額定負(fù)荷工況下應(yīng)力計(jì)算的網(wǎng)格劃分思路和方法。該方法同樣適用于瞬態(tài)分析，可為同類產(chǎn)品在其他工況下的結(jié)構(gòu)強(qiáng)度分析提供參考。對(duì)于汽輪機(jī)啟動(dòng)和停機(jī)等歷時(shí)較長的瞬態(tài)過程，在保證精度的前提下，選取合理高效的網(wǎng)格劃分方案，對(duì)減少求解時(shí)間、提高求解效率更重要。