劉光磊,趙 寧,張瑞庭
(西北工業(yè)大學(xué)機(jī)電學(xué)院,陜西 西安 710072)
在進(jìn)行輪齒加載接觸分析時(shí),如何獲得初始嚙合點(diǎn)處接觸橢圓長(zhǎng)軸上離散點(diǎn)的柔度矩陣是關(guān)鍵技術(shù)之一。許多學(xué)者對(duì)此進(jìn)行了研究[1-7],所提出的方法大致可以分為兩種。一種是直接法[1-5],即直接采用有限元法進(jìn)行分析;另一種是間接法,即先采用有限元法形成齒面網(wǎng)格柔度矩陣,經(jīng)轉(zhuǎn)換,得到輪齒接觸橢圓長(zhǎng)軸上離散點(diǎn)的柔度矩陣[6-7]。直接法計(jì)算結(jié)果精確,但計(jì)算過程復(fù)雜,花費(fèi)時(shí)間長(zhǎng)。間接法計(jì)算精度稍低,但計(jì)算速度快、適應(yīng)面廣。直接法已經(jīng)成熟,間接法尚有待改進(jìn)。
間接法在操作中,文獻(xiàn)[6-7]提出了從齒面網(wǎng)格柔度矩陣向齒面接觸點(diǎn)柔度矩陣轉(zhuǎn)換的二元插值法,但并未給出具體插值的細(xì)節(jié)。為了實(shí)現(xiàn)從齒面網(wǎng)格柔度矩陣到齒面接觸點(diǎn)柔度矩陣的轉(zhuǎn)換,本文以有限元法中的數(shù)值計(jì)算原理為基礎(chǔ),提出從齒面網(wǎng)格柔度矩陣到齒面接觸點(diǎn)柔度矩陣的混合插值法。在求解中發(fā)現(xiàn),在齒面網(wǎng)格結(jié)點(diǎn)施加單位載荷時(shí),相鄰齒面結(jié)點(diǎn)出現(xiàn)齒面畸變。這與彈性力學(xué)中的布辛涅斯克的計(jì)算結(jié)果相左。在網(wǎng)格單元包含兩個(gè)以上的接觸點(diǎn)時(shí),從網(wǎng)格結(jié)點(diǎn)向著它所包圍的接觸點(diǎn)所進(jìn)行的位移插值,有可能出現(xiàn)接觸點(diǎn)柔度矩陣主對(duì)角線不占優(yōu)的情況。為了解決這一問題,提出了齒面網(wǎng)格結(jié)點(diǎn)間距與接觸點(diǎn)間距必須協(xié)調(diào),須保證一個(gè)網(wǎng)格結(jié)點(diǎn)單元只包含一個(gè)接觸點(diǎn)。所提出的方法,物理概念清楚,可操作性強(qiáng)。
本文采用弧齒錐齒輪單齒模型,齒輪副的工作面為小輪凹面和大輪凸面。邊界條件為固定輪轂內(nèi)表面和側(cè)面,計(jì)算齒面網(wǎng)格結(jié)點(diǎn)的柔度矩陣時(shí)采用三維八節(jié)點(diǎn)非協(xié)調(diào)等參單元[8]。齒面加載接觸分析有限元網(wǎng)格模型[9-10]如圖1所示。
圖1 輪齒有限元模型Fig.1 Finite element models of gear teeth
齒面網(wǎng)格柔度矩陣用于反映齒面抵抗外加載荷變形的能力。依次對(duì)齒面上的每個(gè)網(wǎng)格結(jié)點(diǎn)施加單位法向載荷,求得每次加載時(shí)齒面所有結(jié)點(diǎn)沿該處法線方向的位移,也就是各結(jié)點(diǎn)的柔度系數(shù),將它們按行排列,即可組成齒面網(wǎng)格柔度矩陣。圖1中的輪齒工作齒面上共有N個(gè)結(jié)點(diǎn),則齒面網(wǎng)格柔度矩陣為
式中,腳標(biāo)G為齒面網(wǎng)格;矩陣中的元素fij為在網(wǎng)格結(jié)點(diǎn)j作用單位載荷時(shí)網(wǎng)格結(jié)點(diǎn)i產(chǎn)生的位移。
齒面網(wǎng)格柔度矩陣反映了作用在輪齒齒面上各結(jié)點(diǎn)的載荷與變形關(guān)系,而在弧齒錐齒輪加載接觸分析過程中,需要用到反映接觸橢圓長(zhǎng)軸上離散點(diǎn)的載荷與變形關(guān)系,即接觸點(diǎn)柔度矩陣。
接觸點(diǎn)柔度矩陣的含義與齒面網(wǎng)格柔度矩陣相同,即依次對(duì)各接觸點(diǎn)加載時(shí),求得每次加載時(shí)所有接觸點(diǎn)沿法線方向的位移,也就是各接觸點(diǎn)的柔度系數(shù),將它們按行排列,即可組成齒面接觸點(diǎn)柔度矩陣。
設(shè)齒面瞬時(shí)接觸橢圓長(zhǎng)軸上的離散點(diǎn)有n個(gè),按照接觸點(diǎn)柔度矩陣的定義,可得齒面瞬時(shí)接觸橢圓長(zhǎng)軸上離散點(diǎn)的柔度矩陣為
式中,腳標(biāo)C為齒面接觸點(diǎn);矩陣中的元素fij為在接觸點(diǎn)j作用單位載荷時(shí)接觸點(diǎn)i產(chǎn)生的位移。
從[F]G與[F]C的定義可知,盡管兩者物理含義相同,但不能通過簡(jiǎn)單的網(wǎng)格點(diǎn)加密插值,實(shí)現(xiàn)從齒面網(wǎng)格柔度矩陣到接觸點(diǎn)柔度矩陣的轉(zhuǎn)換,而要依據(jù)柔度矩陣的力學(xué)含義,基于有限元數(shù)值計(jì)算原理,實(shí)現(xiàn)兩者的轉(zhuǎn)換。
為了在保證精度的前提下簡(jiǎn)化計(jì)算,將齒面網(wǎng)格和接觸橢圓長(zhǎng)軸上的離散點(diǎn),轉(zhuǎn)換到齒輪的旋轉(zhuǎn)軸截面中,如圖2所示。圖中○表示網(wǎng)格結(jié)點(diǎn),□表示橢圓長(zhǎng)軸上的接觸點(diǎn)。具體插值分兩步完成。
第一步,載荷等效即形成包圍接觸點(diǎn)的網(wǎng)格結(jié)點(diǎn)的等效結(jié)點(diǎn)載荷系數(shù)。以圖2中左邊第一個(gè)接觸點(diǎn)為例,設(shè)想在該點(diǎn)施加單位法向載荷,按照虛功原理,將該載荷等效轉(zhuǎn)換到包圍它的四個(gè)網(wǎng)格結(jié)點(diǎn)上,即
式中,F(xiàn)C為接觸點(diǎn)單位法向載荷;Ni(i=1,2,3,4)為四個(gè)網(wǎng)格結(jié)點(diǎn)的形函數(shù)在接觸點(diǎn)處的值;FGi為等效載荷。因?yàn)镕C為單位載荷,所以FGi也就是等效載荷系數(shù)。
圖2 齒面網(wǎng)格旋轉(zhuǎn)投影示意圖Fig.2 Schematic of rotary prejection of tooth surface grid
這里的形函數(shù)采用雙線性等參數(shù)單元,即
式中,(ξ,η)為局部坐標(biāo),如圖3所示。
圖3 單元坐標(biāo)系Fig.3 Element coordinate system
與網(wǎng)格結(jié)點(diǎn)相關(guān)的形函數(shù)定義于局部坐標(biāo)系中,而實(shí)際齒面網(wǎng)格結(jié)點(diǎn)和接觸點(diǎn)都定義于總體坐標(biāo)中,如圖2和圖3所示。為了計(jì)算形函數(shù)在接觸點(diǎn)處的值,需要建立如下坐標(biāo)變換關(guān)系:
式中,(xi,yi)為網(wǎng)格結(jié)點(diǎn)坐標(biāo);(x,y)為接觸點(diǎn)坐標(biāo)。將(xi,yi)和(x,y)代入上式,解之,得到與(x,y)對(duì)應(yīng)的(ξ,η),再帶入式(4),即可求得Ni(i=1,2,3,4)。
第二步,位移插值即從齒面網(wǎng)格柔度矩陣出發(fā),計(jì)算接觸點(diǎn)的柔度系數(shù),形成齒面柔度矩陣。在計(jì)算時(shí)假設(shè)輪齒為線彈性、小變形,因此可以采用疊加原理。
在接觸點(diǎn)施加單位載荷后,各接觸點(diǎn)均會(huì)產(chǎn)生變形。在得到包圍施加單位載荷的接觸點(diǎn)的網(wǎng)格結(jié)點(diǎn)的等效載荷系數(shù)之后,就可以應(yīng)用疊加原理,將網(wǎng)格結(jié)點(diǎn)的變形等效為接觸點(diǎn)處的變形,即
式中,對(duì)i的求和,表示在同一網(wǎng)格結(jié)點(diǎn)施加等效載荷時(shí),對(duì)接觸點(diǎn)的插值;對(duì)j的求和,表示包圍同一接觸點(diǎn)的四個(gè)網(wǎng)格結(jié)點(diǎn)在接觸點(diǎn)的插值之和。
依次對(duì)各接觸點(diǎn)施加單位法向載荷,采用上述載荷-位移混合插值方法,即可得到單齒接觸點(diǎn)的柔度矩陣。將相互嚙合的輪齒接觸點(diǎn)柔度矩陣相加,最終得到輪齒加載接觸分析(LTCA)所需的接觸點(diǎn)柔度矩陣。
得到齒面網(wǎng)格柔度矩陣以后,圖形顯示其計(jì)算結(jié)果。圖4和圖5為大輪凸面上靠近小端處的網(wǎng)格結(jié)點(diǎn)受到法向單位載荷作用時(shí),整個(gè)齒面及其局部的變形情況。
從圖5可以看出,沿著齒長(zhǎng)方向分布的5個(gè)結(jié)點(diǎn)1、2、3、4、5,加載點(diǎn)3的變形最大,結(jié)點(diǎn)2和4的變形比結(jié)點(diǎn)1和5的變形小,但它們卻比結(jié)點(diǎn)1和5更靠近加載點(diǎn)3,即齒面出現(xiàn)畸變。齒面畸變有可能導(dǎo)致接觸點(diǎn)柔度矩陣的非對(duì)角元出現(xiàn)負(fù)剛度,從而影響接觸點(diǎn)載荷的分布。這種現(xiàn)象與彈性力學(xué)中的布辛涅斯克的計(jì)算結(jié)果相左。為此,作者將有限元模型導(dǎo)入到Ansys中,計(jì)算齒面網(wǎng)格在某點(diǎn)施加法向單位載荷時(shí)的齒面變形,如圖6所示。
圖6 Ansys中分析大輪輪齒模型Fig.6 Ansys analysis model for bull wheel gear tooth
計(jì)算結(jié)果如下:加載點(diǎn)3與圖5中加載點(diǎn)對(duì)應(yīng)。結(jié)點(diǎn)1、2、3、4、5點(diǎn)的法向位移見表1。
表1 Ansys齒面變形計(jì)算結(jié)果Table 1 Results of tooth distortion calculated with Ansys ×10-6
由以上數(shù)據(jù)可以看出,Ansys有限元分析結(jié)果也存在齒面畸變的現(xiàn)象。
在齒面上其他網(wǎng)格結(jié)點(diǎn)加載時(shí),也發(fā)現(xiàn)類似齒面畸變現(xiàn)象。
當(dāng)網(wǎng)格結(jié)點(diǎn)間距過大時(shí),一個(gè)單元內(nèi)可能包含多個(gè)接觸點(diǎn)。如圖2所示,接觸點(diǎn)1、2、3同時(shí)落在一個(gè)網(wǎng)格單元內(nèi)。按照前述載荷-位移插值法,計(jì)算接觸點(diǎn)柔度矩陣時(shí),因載荷與位移插值分別進(jìn)行,有可能導(dǎo)致接觸點(diǎn)柔度矩陣主對(duì)角線元素不占優(yōu)。這不符合柔度矩陣應(yīng)有的物理意義。
以圖2中的接觸點(diǎn)1、2、3為例。設(shè)它們的局部坐標(biāo)分別為(-0.25,-0.25)、(0.25,025)和(0.75,0.75)。當(dāng)對(duì)第二個(gè)接觸點(diǎn)作用法向單位載荷時(shí),從已知的齒面網(wǎng)格結(jié)點(diǎn)柔度矩陣中提取該網(wǎng)格單元結(jié)點(diǎn)的柔度系數(shù)為
把三個(gè)接觸點(diǎn)的局部坐標(biāo)帶入式(4)分別得到三個(gè)接觸點(diǎn)的形函數(shù)為
把接觸點(diǎn)2的形函數(shù)帶入式(3),得到等效載荷系數(shù)為
FGi=[0.1406 0.2344 0.3906 0.2344]
將上述數(shù)據(jù)分別代入式(6)得到三個(gè)接觸點(diǎn)的位移值分別為
由計(jì)算結(jié)果可以看出,加載接觸點(diǎn)2的變形已經(jīng)小于同一網(wǎng)格單元內(nèi)的接觸點(diǎn)1的變形。說明柔度矩陣的主對(duì)角線元素已經(jīng)不再占優(yōu)。
要解決接觸點(diǎn)柔度矩陣主對(duì)角線不占優(yōu)的問題,就要密化齒面網(wǎng)格,使一個(gè)單元內(nèi)最多包含一個(gè)接觸點(diǎn),如圖7所示。從而避免出現(xiàn)柔度矩陣主對(duì)角線不占優(yōu)的現(xiàn)象。
為了保證計(jì)算符合物理意義,同時(shí)兼顧計(jì)算效率,密化網(wǎng)格的尺寸應(yīng)與接觸點(diǎn)距相協(xié)調(diào)。以一個(gè)網(wǎng)格單元正好包含一個(gè)接觸點(diǎn),且相鄰接觸點(diǎn)落在相鄰網(wǎng)格單元為宜。
圖7 密化后的齒面網(wǎng)格Fig.7 Tooth surface grid after densification
(1)本文提出了一種由齒面網(wǎng)格柔度矩陣計(jì)算齒面接觸點(diǎn)柔度矩陣的載荷-位移插值方法。
(2)在求解齒面網(wǎng)格柔度矩陣時(shí),發(fā)現(xiàn)了齒面網(wǎng)格結(jié)點(diǎn)位移畸變現(xiàn)象。
(3)采用載荷-位移插值時(shí),可能導(dǎo)致齒面接觸點(diǎn)柔度矩陣的主對(duì)角線不占優(yōu)。為此,應(yīng)使齒面網(wǎng)格單元的尺寸與接觸點(diǎn)距離相協(xié)調(diào),保證一個(gè)網(wǎng)格單元僅包含一個(gè)接觸點(diǎn)。
[1] BIBEL G D.Contact stress analysis of spiral bevel gears using finite element analysis[J].ASME Journal of Mechanical Design,1995,117(2):235-240.
[2] LITVIN F L,F(xiàn)UENTES A,HAYASAKA K.Design,manufacture,stress analysis,and experimental tests of low-noise high endurance spiral bevel gears[J].Mechanism and Machine Theory,2006,41(1):83-118.
[3] 唐進(jìn)元,蒲太平.弧齒錐齒輪動(dòng)態(tài)嚙合有限元數(shù)值分析[J].機(jī)械科學(xué)與技術(shù),2009,28(8):981-985.
[4] 鄭昌啟.螺旋錐齒輪輪齒加載接觸分析計(jì)算原理[J].機(jī)械工程學(xué)報(bào),1993,29(4):50-54.
[5] 李潤(rùn)方.運(yùn)轉(zhuǎn)中嚙合輪齒的三維應(yīng)力應(yīng)變數(shù)值分析及實(shí)驗(yàn)研究[J].機(jī)械工程學(xué)報(bào),1994,30(2):38-44.
[6] 方宗德.齒輪輪齒承載接觸分析(LTCA)的模型和方法[J].機(jī)械傳動(dòng),1998,22(2):1-3.
[7] 田行斌.弧齒錐齒輪的有摩擦承載接觸分析[J].西北工業(yè)大學(xué)學(xué)報(bào),2000,18(1):19-22.
[8] 劉更.結(jié)構(gòu)動(dòng)力學(xué)有限元程序設(shè)計(jì)[M].北京:國(guó)防科技大學(xué)出版社,1993.
[9] LITVIN F L.Computerized design,simulation of meshing,and contact and stress analysis of facemilled formate generated spiral bevel gears[J].Mechanism and Machine Theory,2002,37:441-459.
[10] 劉光磊,樊紅衛(wèi).含齒根過渡曲線的弧齒錐齒輪有限元精確建模[J].機(jī)械科學(xué)與技術(shù),2010,29(12):1595-1601.