基于有限元分析的大型水輪機(jī)軸聯(lián)銷(xiāo)釘拔銷(xiāo)強(qiáng)度及壽命預(yù)測(cè)方法

錢(qián) 祥，何亞文，張寶慶，李友才，吳秩嵩，周海林，黃俊瑋，劉桂英

（1.湖南五凌電力工程有限公司，湖南 長(zhǎng)沙 410004；2.長(zhǎng)沙理工大學(xué)電氣與信息工程學(xué)院，湖南 長(zhǎng)沙 410114）

0 引言

銷(xiāo)釘拔銷(xiāo)系統(tǒng)是大型水輪機(jī)軸聯(lián)檢修中拆解聯(lián)軸的必要設(shè)備，它的狀態(tài)對(duì)大型水輪機(jī)的安全檢修起著舉足輕重的作用[1，2]。疲勞失效是機(jī)械元件失效中的一種主要形式，銷(xiāo)釘拔銷(xiāo)系統(tǒng)斷裂主要發(fā)生在高應(yīng)力區(qū)域[3，4]。為了提高大型水輪機(jī)檢修拔出銷(xiāo)釘?shù)陌踩?，降低機(jī)械元件的運(yùn)維成本，減少檢修事故的發(fā)生次數(shù)，對(duì)于銷(xiāo)釘拔銷(xiāo)系統(tǒng)強(qiáng)度疲勞程度及壽命的準(zhǔn)確預(yù)測(cè)就具備了至關(guān)重要的工程意義[5]。

已有大量學(xué)者對(duì)于銷(xiāo)釘損傷程度及壽命預(yù)測(cè)進(jìn)行了研究，姜豐 等[6]通過(guò)計(jì)算不同負(fù)荷條件下剩余壽命的分布情況，建立了載荷數(shù)-疲勞壽命的干涉模型。劉儉輝 等[7]采用剩余強(qiáng)度作為控制參數(shù)的傳統(tǒng)應(yīng)力強(qiáng)度干涉模型來(lái)對(duì)機(jī)械元件的疲勞可靠性進(jìn)行分析。余天堂 等[8]提出了一種含有位置參數(shù)的對(duì)數(shù)正態(tài)分布來(lái)進(jìn)行銷(xiāo)釘剩余壽命預(yù)測(cè)的方法，同時(shí)分析了在恒幅值載荷和變幅值載荷作用下銷(xiāo)釘?shù)氖Ｓ鄩勖植夹问?。楊嘉順[9]提出了一種在復(fù)雜應(yīng)力狀態(tài)下仍能適用的剩余壽命預(yù)測(cè)模型，該模型是基于臨界平面法搭建的?，F(xiàn)有的研究中，對(duì)于疲勞壽命的描述主要采用正態(tài)分布以及對(duì)數(shù)正態(tài)分布的方式[10，11]。在工程應(yīng)用可靠性分析中，該分布在產(chǎn)品壽命以及可靠性分析預(yù)測(cè)中應(yīng)用廣泛[12]。但是，專(zhuān)門(mén)針對(duì)銷(xiāo)釘拔銷(xiāo)系統(tǒng)的研究文獻(xiàn)極少。

有鑒于此，本文針對(duì)銷(xiāo)釘拔銷(xiāo)系統(tǒng)的強(qiáng)度和壽命預(yù)測(cè)問(wèn)題進(jìn)行研究，基本思想是：根據(jù)有限元理論構(gòu)建起銷(xiāo)釘拔銷(xiāo)系統(tǒng)的有限元模型，在有限元模型基礎(chǔ)上，根據(jù)各有限單元的應(yīng)力和應(yīng)變分布特征，運(yùn)用Smith-Watson-Topper（SWT）模型理論，建立起包含銷(xiāo)釘拔銷(xiāo)系統(tǒng)的材料屬性及其疲勞參數(shù)、損傷參數(shù)的銷(xiāo)釘疲勞強(qiáng)度與剩余壽命的SWT預(yù)測(cè)模型。在銷(xiāo)釘拔銷(xiāo)系統(tǒng)有限元模型中，簡(jiǎn)化了非線性預(yù)測(cè)問(wèn)題，將其轉(zhuǎn)化為線性預(yù)測(cè)問(wèn)題，以提高基于有限元模型的強(qiáng)度與壽命預(yù)測(cè)模型的求解速度。

1 銷(xiāo)釘拔銷(xiāo)的有限元模型

有限元分析（Finite Element Analysis，F(xiàn)EA）是一種利用數(shù)學(xué)近似對(duì)真實(shí)物理系統(tǒng)進(jìn)行模擬的方法?；舅枷胧怯幂^簡(jiǎn)單的問(wèn)題代替復(fù)雜問(wèn)題后再求解，利用簡(jiǎn)單而又相互作用的有限微小單元及其有限數(shù)量的未知量去逼近無(wú)限未知量的真實(shí)系統(tǒng)，將整體系統(tǒng)的求解域看成是由許多稱(chēng)為有限元的小的互連子域組成，對(duì)每一單元假定一個(gè)合適的近似解，然后推導(dǎo)求解這個(gè)域總的滿(mǎn)足條件（如結(jié)構(gòu)的平衡條件），從而得到整體系統(tǒng)的解[13]。

銷(xiāo)釘拔銷(xiāo)系統(tǒng)的有限元建模原理如下：

（1）建立單元位移模式模型

將銷(xiāo)釘拔銷(xiāo)系統(tǒng)中的銷(xiāo)釘劃分為有限數(shù)量的N個(gè)微小單元u，微小單元u為剛性六面體，每個(gè)單元有8個(gè)節(jié)點(diǎn)，每個(gè)節(jié)點(diǎn)有6個(gè)自由度的位移，其中3個(gè)距離位移，3個(gè)相位位移。因此，第n（n=1,2,…,N。下同）單元un的第i節(jié)點(diǎn)的未知量包括6個(gè)應(yīng)力變量、6個(gè)應(yīng)變變量、3個(gè)位移變量，分別構(gòu)成應(yīng)力向量、應(yīng)變向量、位移向量。位移未知量也稱(chēng)基本未知量。第n單元un的第i節(jié)點(diǎn)的應(yīng)力向量可以表示為

第n單元un的第i節(jié)點(diǎn)的應(yīng)變向量可以表示為

第n單元un的第i節(jié)點(diǎn)的位移向量可以表示為

第n單元的全部節(jié)點(diǎn)的位移向量構(gòu)成的節(jié)點(diǎn)位移擴(kuò)展向量可以表示為

將微小單元u看成一個(gè)質(zhì)點(diǎn)，則第n單元un有6個(gè)應(yīng)力變量、6個(gè)應(yīng)變變量、3個(gè)位移變量，因此，第n單元un的應(yīng)力向量可以表示為

第n單元un的應(yīng)變向量可以表示為

第n單元un的位移模式θun表示為

設(shè)第n單元的位移模式使用線性差值表示，則θnj可以表示為

寫(xiě)成矩陣形式，有

對(duì)于任意單元n的任意節(jié)點(diǎn)i的位移變量λnij應(yīng)滿(mǎn)足：

寫(xiě)成矩陣形式，有

其中，xni、yni、zni為位移向量λni在xyz坐標(biāo)系中的x、y、z軸上的坐標(biāo)。

根據(jù)式（7）~（11），可以得到

其中，

稱(chēng)為形態(tài)函數(shù)（形函數(shù)）矩陣。

（2）建立單元應(yīng)力向量與應(yīng)變向量關(guān)系模型

第n單元un的單元應(yīng)力向量δun與應(yīng)變向量關(guān)系εun可以寫(xiě)成

式（13）中的矩陣D為

其中，E為楊氏彈性模量，μ為柏松比。

（3）建立單元應(yīng)變向量與位移向量關(guān)系模型

第n單元un的單元應(yīng)變向量與節(jié)點(diǎn)位移擴(kuò)展向量關(guān)系可以寫(xiě)成

式（15）中，B=[Bn1Bn2… Bn8]，其中，

（4）建立單元?jiǎng)偠染仃嚹Ｐ团c載荷向量模型

根據(jù)式（13）、（15），可以得到單元應(yīng)力向量δun與節(jié)點(diǎn)位移擴(kuò)展向量的關(guān)系為

其中，S為應(yīng)力矩陣，S=DB。

基于虛功原理，可得到單元的節(jié)點(diǎn)力擴(kuò)展向量（即為單元載荷向量）與節(jié)點(diǎn)位移擴(kuò)展向量的關(guān)系為

其中，F(xiàn)n也稱(chēng)為節(jié)點(diǎn)力擴(kuò)展向量，

式（17）也稱(chēng)單元特性方程；Kn稱(chēng)為單元?jiǎng)偠染仃嚕琄n=∫VnBTDBdV=BTDBVn。Vn表示單元的體積。

（5）建立整體的總剛度矩陣模型與總載荷向量模型

整體結(jié)構(gòu)的總剛度矩陣的總載荷向量分別是單元?jiǎng)偠染仃嚭蛦卧?jié)點(diǎn)力擴(kuò)展向量的疊加，因此有

如果單元載荷比較復(fù)雜，按如下表達(dá)式求取等效載荷作為單元載荷向量。

其中，F(xiàn)c、Fs、FV分別表示作用于單元上的集中力、面力和體力。

（6）建立節(jié)點(diǎn)載荷平衡方程

在整體結(jié)構(gòu)中，一個(gè)節(jié)點(diǎn)為幾個(gè)單元所共有，因此，對(duì)每個(gè)節(jié)點(diǎn)i，均有

其中，ui表示共有節(jié)點(diǎn)i（與節(jié)點(diǎn)i相連）的單元的集合，Tni、Ri分別表示第n單元的第i節(jié)點(diǎn)（假設(shè)共有的節(jié)點(diǎn)在單元中的編號(hào)相同）的節(jié)點(diǎn)力（也稱(chēng)內(nèi)載荷）和外力（也稱(chēng)外載荷），因此，前述的F是T（內(nèi)力）和R（外力）之和。

（7）建立邊界條件

1）力條件

邊界上的單位面積的外載荷的面力Rs等于內(nèi)載荷的面力Ts。

2）位移條件

（8）計(jì)算節(jié)點(diǎn)應(yīng)力

根據(jù)節(jié)點(diǎn)載荷平衡方程和邊界條件，首先，求解式（17）得到位移λn，然后，根據(jù)式（18）計(jì)算各單元的節(jié)點(diǎn)力擴(kuò)展向量Fn，根據(jù)式（16），計(jì)算各單元的單元應(yīng)力向量δun，從而預(yù)測(cè)出銷(xiāo)釘拔銷(xiāo)的受力分布。

2 基于有限分析模型的壽命預(yù)測(cè)SWT模型

銷(xiāo)釘拔銷(xiāo)系統(tǒng)的疲勞主要是以銷(xiāo)釘?shù)奈?dòng)疲勞為主。微動(dòng)疲勞是一種常見(jiàn)的接觸疲勞現(xiàn)象，它發(fā)生在兩個(gè)表面緊密配合并有微小的相對(duì)運(yùn)動(dòng)時(shí)。這種運(yùn)動(dòng)會(huì)在接觸區(qū)域產(chǎn)生高應(yīng)力和高應(yīng)變，從而引起表面的損傷和裂紋。銷(xiāo)釘拔銷(xiāo)系統(tǒng)作為銷(xiāo)釘拔出過(guò)程使用的一種重要器件，反復(fù)使用受到微動(dòng)疲勞的影響，其壽命會(huì)因此降低。本文采用臨界平面法中的SWT模型來(lái)建立銷(xiāo)釘拔銷(xiāo)疲勞受損的剩余壽命預(yù)測(cè)模型[10]。

臨界平面法是一種分析裂紋萌生的方法[14]，基本思想是：將裂紋萌生問(wèn)題轉(zhuǎn)化為尋找臨界平面和綜合參數(shù)的問(wèn)題，在復(fù)雜應(yīng)力狀態(tài)下，存在一個(gè)最有利于裂紋萌生的平面，稱(chēng)為臨界平面，在有限分析模型基礎(chǔ)上，通過(guò)在有限單元積分點(diǎn)計(jì)算一個(gè)綜合參數(shù)，可以比較不同截面的裂紋萌生傾向，從而預(yù)測(cè)工程構(gòu)件的疲勞壽命。

2.1 基于有限分析模型的微動(dòng)損傷預(yù)測(cè)模型

裂紋的萌生通常發(fā)生在一個(gè)特定的危險(xiǎn)平面上，這個(gè)平面可以用一個(gè)綜合參數(shù)來(lái)描述，稱(chēng)為臨界平面參數(shù)。臨界平面參數(shù)可以反映材料在復(fù)雜應(yīng)力狀態(tài)下的疲勞性能，有多種形式，如SWT模型、FP模型、SSR模型、MSSR模型等[15]。本文采用SWT模型建立銷(xiāo)釘拔銷(xiāo)系統(tǒng)微動(dòng)損傷預(yù)測(cè)模型，如式（21）。根據(jù)銷(xiāo)釘拔銷(xiāo)系統(tǒng)的銷(xiāo)釘微動(dòng)損傷度，從而預(yù)測(cè)出銷(xiāo)釘拔銷(xiāo)系統(tǒng)的耐力強(qiáng)度。

式中:D為微動(dòng)損傷參量，δmax為最大應(yīng)力,εa為應(yīng)變幅值。

2.2 基于微動(dòng)損傷模型的疲勞壽命預(yù)測(cè)模型

微動(dòng)疲勞壽命預(yù)測(cè)模型的核心在于解決接觸問(wèn)題[13]。接觸表面的應(yīng)力和位移是微動(dòng)疲勞的兩個(gè)主要控制因素。為了提高應(yīng)力和位移計(jì)算的精度，本文采用有限元分析模型進(jìn)行計(jì)算，然后根據(jù)臨界平面法的SWT模型，建立了微動(dòng)疲勞壽命預(yù)測(cè)模型。微動(dòng)損傷疲勞是一種由于微小的位移或變形而導(dǎo)致材料失效的現(xiàn)象。微動(dòng)損傷疲勞的影響因素眾多,包括材料的性質(zhì)、微動(dòng)損傷幅度、作用力幅值及其頻率、溫度、環(huán)境等?；谖?dòng)損傷的疲勞壽命預(yù)測(cè)模型可以表示為:

式中:A和m為材料特性相關(guān)的系數(shù)，δ為應(yīng)力，f為力的頻率。

2.3 強(qiáng)度與壽命預(yù)測(cè)流程

銷(xiāo)釘拔銷(xiāo)系統(tǒng)的強(qiáng)度與壽命預(yù)測(cè)流程如圖1所示，通過(guò)建立銷(xiāo)釘拔銷(xiāo)系統(tǒng)的有限元分析模型，計(jì)算出銷(xiāo)釘拔銷(xiāo)系統(tǒng)的節(jié)點(diǎn)力和單元應(yīng)力分布，在此基礎(chǔ)上，計(jì)算銷(xiāo)釘拔銷(xiāo)系統(tǒng)微動(dòng)損傷度和疲勞壽命。

圖1 強(qiáng)度與壽命預(yù)測(cè)流程圖

3 仿真與試驗(yàn)驗(yàn)證

3.1 拔銷(xiāo)過(guò)程的銷(xiāo)釘應(yīng)力分布仿真

現(xiàn)在常用的有限元軟件有Ansys、Abaqus、Cosmos和Marc等。Abaqus是一款專(zhuān)業(yè)的有限元分析軟件，在國(guó)際上的知名度都是首屈一指的，完美地融合了同類(lèi)軟件中的優(yōu)點(diǎn)，擁有各種類(lèi)型的材料模型庫(kù)，可以模擬典型工程材料的性能，適用于相對(duì)簡(jiǎn)單的線性分析和復(fù)雜的非線性分析。本文使用Abaqus 2022對(duì)銷(xiāo)釘拔銷(xiāo)系統(tǒng)的銷(xiāo)釘拔銷(xiāo)過(guò)程進(jìn)行有限元分析仿真計(jì)算。

仿真的銷(xiāo)釘拔銷(xiāo)系統(tǒng)的結(jié)構(gòu)示意圖如圖2所示，銷(xiāo)釘是構(gòu)成銷(xiāo)釘拔銷(xiāo)系統(tǒng)的核心元件，銷(xiāo)釘為圓柱體金屬構(gòu)件，仿真的大型水輪機(jī)軸聯(lián)的銷(xiāo)釘拔銷(xiāo)系統(tǒng)的銷(xiāo)釘為6個(gè)，呈等間隔圓環(huán)設(shè)置，使用6個(gè)銷(xiāo)釘聯(lián)接的大型水輪機(jī)軸聯(lián)系統(tǒng)的Abaqus 2022仿真模型結(jié)構(gòu)圖如圖3所示，仿真的銷(xiāo)釘拔銷(xiāo)系統(tǒng)的銷(xiāo)釘長(zhǎng)度400 mm，直徑120 mm，彈性模量為206 MPa，泊松比為0.3，銷(xiāo)釘拔銷(xiāo)系統(tǒng)在拔銷(xiāo)過(guò)程中的外載荷（外力）900 MPa，屬于面力。銷(xiāo)釘拔銷(xiāo)系統(tǒng)拔銷(xiāo)過(guò)程的銷(xiāo)釘云圖如圖4所示。由圖4可以看出，在銷(xiāo)釘?shù)?處出現(xiàn)較高應(yīng)力，此應(yīng)力集中現(xiàn)象與環(huán)帶的微動(dòng)疲勞問(wèn)題有關(guān)，符合實(shí)際構(gòu)件破壞情況。

圖2 銷(xiāo)釘拔銷(xiāo)系統(tǒng)的結(jié)構(gòu)示意圖

圖3 銷(xiāo)釘聯(lián)接的大型水輪機(jī)軸聯(lián)系統(tǒng)的Abaqus模型

圖4 銷(xiāo)釘拔銷(xiāo)過(guò)程的銷(xiāo)釘應(yīng)力分布仿真模型

3.2 銷(xiāo)釘拔銷(xiāo)壽命的仿真和實(shí)驗(yàn)

銷(xiāo)釘拔銷(xiāo)系統(tǒng)仿真與實(shí)驗(yàn)的銷(xiāo)釘采用鍛鋼35CrMo材質(zhì)，銷(xiāo)釘?shù)幕瘜W(xué)成分及機(jī)械性能參數(shù)見(jiàn)表1。銷(xiāo)釘拔銷(xiāo)系統(tǒng)仿真，采用Abaqus 2022+Maltlab仿真軟件。銷(xiāo)釘拔銷(xiāo)系統(tǒng)的實(shí)驗(yàn)過(guò)程，采用電子式萬(wàn)能試驗(yàn)機(jī)模擬銷(xiāo)釘拔銷(xiāo)，采用直讀光譜儀進(jìn)行監(jiān)測(cè)，試驗(yàn)在室溫空氣介質(zhì)中進(jìn)行，仿真和實(shí)驗(yàn)結(jié)果如圖5和圖6所示。圖5為銷(xiāo)釘拔銷(xiāo)系統(tǒng)的銷(xiāo)釘SWT參數(shù)與疲勞壽命的關(guān)系曲線，圖6為銷(xiāo)釘疲勞壽命的仿真和實(shí)驗(yàn)預(yù)測(cè)結(jié)果。

表1 銷(xiāo)釘?shù)幕瘜W(xué)成分與機(jī)械性能參數(shù)

圖5 銷(xiāo)釘SWT參數(shù)與疲勞壽命的關(guān)系

圖6 銷(xiāo)釘疲勞壽命的仿真和實(shí)驗(yàn)預(yù)測(cè)結(jié)果

由圖5可以看出，SWT參數(shù)的銷(xiāo)釘疲勞數(shù)據(jù)擬合的效果較好。由圖6可以看出，實(shí)驗(yàn)預(yù)測(cè)結(jié)果基本上均布在試驗(yàn)結(jié)果的兩側(cè)，效果較為理想，銷(xiāo)釘疲勞壽命實(shí)驗(yàn)與仿真預(yù)測(cè)的誤差因子小于1.5，仿真預(yù)測(cè)結(jié)果小于實(shí)驗(yàn)結(jié)果，仿真預(yù)測(cè)結(jié)果偏于保守，但是與實(shí)驗(yàn)結(jié)果的偏差較小，由此可見(jiàn)，本文所研究的壽命預(yù)測(cè)模型具有較高的疲勞壽命預(yù)測(cè)精度。

4 結(jié)論

本文旨在利用有限元法和SWT對(duì)大型水輪機(jī)軸聯(lián)銷(xiāo)釘拔銷(xiāo)系統(tǒng)的強(qiáng)度和壽命進(jìn)行預(yù)測(cè)。通過(guò)理論分析和仿真與實(shí)驗(yàn)，可以得到以下結(jié)論：

（1）基于有限元分析，可以確定銷(xiāo)釘拔銷(xiāo)系統(tǒng)的銷(xiāo)釘是否在承受負(fù)載時(shí)會(huì)發(fā)生塑性變形或破壞。如果超過(guò)了材料的屈服強(qiáng)度，銷(xiāo)釘拔銷(xiāo)系統(tǒng)可能會(huì)發(fā)生塑性變形。如果應(yīng)力達(dá)到了銷(xiāo)釘拔銷(xiāo)系統(tǒng)的破壞強(qiáng)度，銷(xiāo)釘可能會(huì)發(fā)生破壞。

（2）建立的銷(xiāo)釘拔銷(xiāo)系統(tǒng)疲勞壽命預(yù)測(cè)SWT模型，在估算疲勞壽命時(shí)，以臨界平面上的損傷為主，同時(shí)了考慮材料的總體損傷程度，適用于變幅或隨機(jī)多軸載荷情況下的多損傷平面情況的壽命預(yù)測(cè)。以臨界平面上最大剪應(yīng)變幅、正應(yīng)變幅為參數(shù)建立微動(dòng)損傷和疲勞壽命預(yù)測(cè)SWT模型，能夠得到較高的預(yù)測(cè)精度，誤差因子小于1.5。

（3）基于有限元法和SWT的預(yù)測(cè)方法，可以較高精度地評(píng)估銷(xiāo)釘拔銷(xiāo)系統(tǒng)的強(qiáng)度和壽命。該方法為大型水輪機(jī)軸聯(lián)銷(xiāo)釘拔銷(xiāo)設(shè)計(jì)和工程實(shí)踐提供了重要的技術(shù)支持，能夠提高水輪機(jī)的可靠性，降低運(yùn)維成本，并減少檢修事故的發(fā)生。