基于解析法和有限元法融合的準(zhǔn)雙曲面齒輪接觸分析研究

劉少杭，聶少武，蔣闖，柴福博

(河南科技大學(xué)機(jī)電工程學(xué)院，河南洛陽 471003)

0 前言

準(zhǔn)雙曲面齒輪具有重合度大、傳動(dòng)平穩(wěn)、承載能力高等優(yōu)點(diǎn)，廣泛應(yīng)用于各種傳動(dòng)裝置，在汽車驅(qū)動(dòng)后橋中的應(yīng)用尤為廣泛。近些年來，為了對(duì)汽車驅(qū)動(dòng)橋準(zhǔn)雙曲面齒輪的嚙合性能進(jìn)行評(píng)價(jià)，國內(nèi)外學(xué)者對(duì)準(zhǔn)雙曲面齒輪的齒面接觸分析(Tooth Contact Analysis，TCA)進(jìn)行了一系列的研究。王星等人[1-2]建立HGT準(zhǔn)雙曲面齒輪切齒加工數(shù)學(xué)模型，對(duì)齒面接觸和邊緣接觸進(jìn)行了分析。杜進(jìn)輔等[3]通過分析克林貝格制擺線齒準(zhǔn)雙曲面齒輪，建立一種改進(jìn)的齒面接觸分析模型并對(duì)其進(jìn)行分析，避免了TCA結(jié)果出現(xiàn)幾何上的不準(zhǔn)確。杜進(jìn)輔等[4]通過建立奧利康制準(zhǔn)雙曲面齒輪的全齒面模型，對(duì)其進(jìn)行了齒面接觸分析。SHIH等[5]建立了采用端面銑削和端面滾切加工的準(zhǔn)雙曲面齒輪齒面接觸分析的統(tǒng)一數(shù)學(xué)模型。王崢嶸[6]對(duì)采用端面滾齒法加工準(zhǔn)雙曲面齒輪提出了兩種基于數(shù)值迭代的齒面接觸分析方法，避免了復(fù)雜的主曲率主方向的推導(dǎo)。卓耀彬等[7]建立了適合準(zhǔn)靜態(tài)齒面接觸分析的準(zhǔn)雙曲面齒輪傳動(dòng)系統(tǒng)，并對(duì)其進(jìn)行了準(zhǔn)靜態(tài)齒面接觸分析。蘇進(jìn)展等[8]對(duì)弧齒錐齒輪進(jìn)行了數(shù)值化齒面接觸分析，實(shí)現(xiàn)了數(shù)值化齒面的幾何特征修改。

目前，對(duì)準(zhǔn)雙曲面齒輪副嚙合性能的評(píng)價(jià)通常是采用輪齒接觸分析(TCA)方法。傳統(tǒng)TCA方法是通過解析法求解齒面接觸區(qū)，因涉及到主方向、主曲率等復(fù)雜曲率計(jì)算，推導(dǎo)過程十分繁瑣且計(jì)算量大，甚至有時(shí)出現(xiàn)齒面接觸橢圓求解異常的情況。因此，本文作者提出一種基于解析法和有限元法融合的齒面接觸分析方法，即利用解析法計(jì)算齒面?zhèn)鲃?dòng)誤差，通過有限元軟件分析加載接觸，得到齒面接觸區(qū)，通過與專業(yè)軟件進(jìn)行對(duì)比，驗(yàn)證該方法的可行性。

1 HFT數(shù)學(xué)模型的建立

根據(jù)砂輪結(jié)構(gòu)建立刀盤坐標(biāo)系[9]，如圖1所示。其中，坐標(biāo)系St(Xt,Yt,Zt)與刀盤固連；Gg和Gp是內(nèi)外刀刀刃上的一點(diǎn)；rg和rp分別是內(nèi)外刀刀尖半徑；αg和αp分別為內(nèi)外刀齒齒形角；αg取正值；αp取負(fù)值；ug和up為刀齒曲面參數(shù)。

在HFT法加工過程中，根據(jù)刀具與工件的相對(duì)運(yùn)動(dòng)和位置變化，建立了小輪磨齒加工數(shù)學(xué)模型，如圖2所示。

坐標(biāo)系Sm(Xm,Ym,Zm)與加工機(jī)床固連稱為機(jī)床坐標(biāo)系，Om與機(jī)床中心重合，Sc(Xc,Yc,Zc)與搖臺(tái)固連稱為搖臺(tái)坐標(biāo)系，加工前坐標(biāo)Sm和Sc重合。j為刀傾角，i為刀轉(zhuǎn)角，q1為角向刀位，Sr1為徑向刀位，Em1為垂直輪位，Xb為床位，Xp1為水平輪位，γm為工件安裝角度，θ為刀具齒面參數(shù)，Sd與工件固連稱為工件坐標(biāo)系，通過展成運(yùn)動(dòng)關(guān)系可知φ=Rbφ，Rb為切齒滾比。

根據(jù)圖1刀盤坐標(biāo)系中的矢量關(guān)系，可以推導(dǎo)出刀具曲面方程以及單位法矢。

(1)

nt=[-cosαkcosθk-cosαksinθksinαk]

(2)

式中：k=g、p,分別代表內(nèi)刀和外刀。

由St坐標(biāo)系中的刀具曲面方程和單位法矢經(jīng)過一系列坐標(biāo)變換和矢量運(yùn)算得到小輪齒面方程和單位法矢。

r1(uk,θk,φ)=M1hMhsMsmMmcMcjMjtrt(uk,θk)

(3)

n1(uk,θk,φ)=LdhLhsLsmLmcLcjLjtnt(uk,θk)

(4)

式中:Ldh、Lhs、Lsm、Lmc、Lcj、Ljt分別為矩陣Mdh、Mhs、Msm、Mmc、Mcj、Mjt去掉最后一行及最后一列得到的3×3矩陣。

大輪采用成形法加工，成形法加工數(shù)學(xué)模型是刀傾法加工數(shù)學(xué)模型的一種特殊情況，因此只需令圖2中的j=0、Xb=0、Em1=0、Rb=0、φ=0、φ=0，即可推導(dǎo)出大輪齒面方程和單位法矢。

2 解析法求解傳動(dòng)誤差

齒面接觸分析是通過數(shù)學(xué)方法模擬出齒輪副的嚙合過程，依據(jù)準(zhǔn)雙曲面齒輪副在嚙合過程中的位置和運(yùn)動(dòng)關(guān)系[10]，建立準(zhǔn)雙曲面齒輪嚙合數(shù)學(xué)模型，如圖3所示。

圖3中小輪軸線和大輪軸線夾角為90°，軸線之間距離為ΔE，S2與大輪固連，S1與小輪固連，Sb和Sd為輔助坐標(biāo)系，嚙合過程中S2繞大輪軸線旋轉(zhuǎn)φ2，S1繞小輪軸線旋轉(zhuǎn)φ1。

通過上述建立的刀傾法加工數(shù)學(xué)模型即可推導(dǎo)出小輪齒面方程r1(u1,θ1)、小輪法矢n1(u1,θ1)和大輪齒面方程r2(u2,θ2)、大輪法矢n2(u2,θ2)。分別將它們進(jìn)行坐標(biāo)變換，變換到固定坐標(biāo)系Sa中，得到如下：

ra1(u1,θ1,φ1)=MadMd1r1(u1,θ1)

(5)

ra2(u2,θ2,φ2)=MabMb2r2(u2,θ2)

(6)

na1(u1,θ1,φ1)=LadLd1n1(u1,θ1)

(7)

na2(u2,θ2,φ2)=LabLb2n2(u2,θ2)

(8)

式中：Lad、Ld1、Lab、Lb2分別為矩陣Mad、Md1、Mab、Mb2去掉最后一行及最后一列得到的3×3矩陣。

在齒輪副嚙合過程中，齒面Σ1和齒面Σ2持續(xù)相切接觸，故在任一時(shí)刻兩齒面都有公共接觸點(diǎn)，在公共接觸點(diǎn)處兩齒面有相同的徑矢ram和相同法矢nam。

ra1(u1,θ1,φ1)=ra2(u2,θ2,φ2)

(9)

na1(u1,θ1,φ1)=na2(u2,θ2,φ2)

(10)

由式(9)(10)可得到5個(gè)相互獨(dú)立的標(biāo)量方程，5個(gè)標(biāo)量方程中含有6個(gè)未知數(shù)，通過給定φ1值，求解出其他5個(gè)未知數(shù)(u1,u2,θ1,θ2,φ2)，得到第一個(gè)接觸點(diǎn)，將它設(shè)為參考點(diǎn)。以此為初值，通過以一定的步長改變?chǔ)?值求解出齒面上的所有點(diǎn)。

由傳動(dòng)誤差的定義，得到其表達(dá)式[11-12]如下:

(11)

3 齒面接觸區(qū)有限元仿真

通過加工數(shù)學(xué)模型求解出齒面方程，采用旋轉(zhuǎn)投影原理進(jìn)行網(wǎng)格劃分進(jìn)而求解出齒面點(diǎn)，將齒面點(diǎn)導(dǎo)入U(xiǎn)G中建立準(zhǔn)雙曲面齒輪模型[13]。為便于有限元軟件的仿真分析，將準(zhǔn)雙曲面齒輪模型分割為單齒模型。為便于六面體網(wǎng)格的劃分，建立輔助面將單齒模型進(jìn)行切分，最后將模型和輔助面導(dǎo)出為(.x_t)文件。

ABAQUS齒面接觸有限元仿真流程如圖4所示。仿真流程主要分為前處理、求解和后處理3個(gè)部分。在前處理過程中，通過導(dǎo)入上述(.x_t)文件，進(jìn)行裝配，完成幾何模型的構(gòu)建。通過設(shè)置材料密度、彈性模量、泊松比完成材料參數(shù)的定義。通過對(duì)大小輪各邊布點(diǎn)，劃分網(wǎng)格完成網(wǎng)格系統(tǒng)的構(gòu)建。

在求解過程中，求解器的設(shè)定采用靜力通用類型。連接關(guān)系采用面與面接觸。邊界條件的設(shè)定通過設(shè)置邊界約束以及施加負(fù)載扭矩實(shí)現(xiàn)。

后處理過程是當(dāng)上述所有操作設(shè)置完畢后，提交作業(yè)，得到場輸出數(shù)據(jù)中的接觸應(yīng)力。

接觸橢圓是齒輪副嚙合產(chǎn)生的瞬時(shí)接觸區(qū)，不同時(shí)刻瞬時(shí)接觸橢圓的集合稱為齒面接觸區(qū)[14-15]。因此，本文作者通過使用Python語言編寫腳本，對(duì)分析結(jié)果進(jìn)行提取。Python提取ABAQUS中輸出結(jié)果的流程如圖5所示。

Python提取流程圖主要分為讀取數(shù)據(jù)和輸出數(shù)據(jù)兩部分。在讀取數(shù)據(jù)階段，首先打開odb文件讀取場輸出結(jié)果中的齒面接觸應(yīng)力，然后定義2個(gè)數(shù)組，分別儲(chǔ)存單元編號(hào)和對(duì)應(yīng)的接觸應(yīng)力值，最后關(guān)閉結(jié)果文件。在輸出數(shù)據(jù)階段，首先打開新建的儲(chǔ)存文件，然后通過比較每個(gè)單元在不同時(shí)刻的接觸應(yīng)力，找到最大接觸應(yīng)力，并將每個(gè)單元的最大接觸應(yīng)力值寫入文件中，最后關(guān)閉文件。

ABAQUS軟件在分析結(jié)束后，關(guān)閉odb文件，運(yùn)行上述的腳本文件，即可得到主從動(dòng)輪接觸面上的每個(gè)單元在整個(gè)分析過程中不同時(shí)刻的最大接觸應(yīng)力。每一時(shí)刻顯示的接觸應(yīng)力云圖形狀等同于齒面瞬時(shí)接觸橢圓。通過提取不同時(shí)刻的最大接觸應(yīng)力，得到齒面接觸區(qū)。

4 算例分析

以一對(duì)HFT法加工的準(zhǔn)雙曲面齒輪為例，分別進(jìn)行建模以及有限元接觸分析。表1所示為準(zhǔn)雙曲面齒輪幾何參數(shù)，表2所示為加工參數(shù)。

表1 準(zhǔn)雙曲面齒輪幾何參數(shù)

表2 加工參數(shù)

利用解析法對(duì)輪齒接觸原理進(jìn)行研究，通過 MATLAB軟件編寫傳動(dòng)誤差程序。圖6所示為由MATLAB程序計(jì)算得出的傳動(dòng)誤差，圖7所示為由格里森CAGE軟件得到的傳動(dòng)誤差。

由圖6—圖7可以看出：解析法得到的傳動(dòng)誤差幅值(交點(diǎn)對(duì)應(yīng)的縱坐標(biāo))為28 μrad，格里森CAGE軟件得到的傳動(dòng)誤差幅值為28.5 μrad，兩者傳動(dòng)誤差幅值基本相同，驗(yàn)證了解析法求解傳動(dòng)誤差的正確性。

在UG中建立準(zhǔn)雙曲面齒輪模型，如圖8所示。

將建立好的模型導(dǎo)入有限元軟件ABAQUS中進(jìn)行負(fù)載接觸分析，得到分析結(jié)果后再進(jìn)行后處理，即對(duì)分析結(jié)果進(jìn)行提取，得到齒面接觸區(qū)，如圖9—圖11所示。

然后將準(zhǔn)雙曲面齒輪的幾何參數(shù)以及加工參數(shù)輸入到軟件MASTA中，建立齒輪副模型，通過施加不同載荷，得到齒輪副在不同負(fù)載下的齒面接觸區(qū)，如圖12所示。

通過對(duì)比ABAQUS仿真齒面接觸區(qū)與MASTA軟件仿真結(jié)果，發(fā)現(xiàn)接觸區(qū)都是位于齒面中部，成內(nèi)對(duì)角接觸，齒面接觸區(qū)形狀大小趨于一致，驗(yàn)證了文中所提出的齒面接觸分析方法的正確性。

為進(jìn)一步對(duì)有限元仿真結(jié)果進(jìn)行驗(yàn)證，對(duì)表1中準(zhǔn)雙曲面齒輪進(jìn)行磨齒加工和齒面誤差測量實(shí)驗(yàn)，如圖13和圖14所示。通過對(duì)磨齒后的齒面進(jìn)行測量反調(diào)，保證了實(shí)際加工齒面與設(shè)計(jì)齒面一致，最后在滾檢機(jī)T50上進(jìn)行了齒面滾檢實(shí)驗(yàn)(加載扭矩40 N·m)，得到了輕載下大輪凹面接觸區(qū)，如圖15所示。

通過對(duì)比滾檢結(jié)果圖15(b)與有限元仿真結(jié)果圖11(b)可以看出：仿真接觸區(qū)的形狀、位置、大小與滾檢結(jié)果趨于一致，進(jìn)一步證明了所提出的齒面接觸分析方法的正確性與可行性。

5 結(jié)論

基于解析法和有限元法對(duì)準(zhǔn)雙曲面齒輪齒面接觸分析進(jìn)行了研究，通過解析法求解出了傳動(dòng)誤差曲線，通過有限元仿真獲得了齒面三維接觸區(qū)，避免了解析法在求解接觸區(qū)橢圓時(shí)復(fù)雜的曲率計(jì)算，同時(shí)使得齒面接觸區(qū)顯示更加真實(shí)直觀。通過對(duì)比仿真結(jié)果與MASTA軟件結(jié)果，以及齒面加工滾檢實(shí)驗(yàn)，驗(yàn)證了所提齒面接觸分析方法的正確性與可行性。