敬 謙，劉宏昭，王庚祥

(1.西安理工大學(xué) 機(jī)械與精密儀器工程學(xué)院，西安 710048；2.隴東學(xué)院 機(jī)械工程學(xué)院，甘肅 慶陽 745000)

機(jī)構(gòu)系統(tǒng)關(guān)節(jié)反力是約束反力的一種，關(guān)節(jié)反力計算的準(zhǔn)確性對機(jī)構(gòu)動力學(xué)分析、工作特性分析、結(jié)構(gòu)設(shè)計、零部件選取以及機(jī)器工作壽命預(yù)測等意義重大。近年來，機(jī)構(gòu)系統(tǒng)運動中如何引入拉格朗日乘子并明確其物理意義仍然是國內(nèi)外研究的熱點[1]。關(guān)于拉格朗日乘子與關(guān)節(jié)反力的關(guān)系問題，丁光濤[2]從理論分析的角度討論了完整約束和非完整約束力學(xué)系統(tǒng)中引入待定乘子兩種不同的途徑，著重研究了變分原理條件極值中引入待定乘子修正系統(tǒng)的拉格朗日函數(shù)的方式，給出了拉格朗日乘子與理想約束反力之間的關(guān)系表達(dá)式。

關(guān)節(jié)反力常見的計算方法有牛頓歐拉方程、達(dá)朗貝爾原理以及拉格朗日方程。趙燕等[3]利用牛頓歐拉法列出所有構(gòu)件的力和力矩平衡方程，確定了驅(qū)動力和平臺受到的約束反力，但同時也指出當(dāng)機(jī)構(gòu)結(jié)構(gòu)復(fù)雜或需要求解的主動力和運動副反力的數(shù)目更多時，這種方法很不方便。Marek 等[4-5]討論了包含有冗余約束和奇異約束的剛體系統(tǒng)約束反力的求解方法，給出了在包含非獨立約束條件下能夠確定單一關(guān)節(jié)反力或一組關(guān)節(jié)反力的前提條件。通過整車系統(tǒng)動力學(xué)分析,劉剛同樣選用牛頓歐拉方程對工程機(jī)構(gòu)鉸接系統(tǒng)鉸點動態(tài)約束反力進(jìn)行了研究[6]。文獻(xiàn)[7]建立了多對均布傳力桿空間RCCR機(jī)構(gòu)模型，利用達(dá)朗貝爾原理計算出了在動力條件下機(jī)構(gòu)慣性力引起約束反力的理論表達(dá)式。拉格朗日增廣法將加速度約束方程與動力學(xué)方程耦合進(jìn)行積分，求解時必須考慮速度和位移約束的違約修正,如果修正不當(dāng),往往會使得求解發(fā)散。為了避免違約修正，原亮明等[8]將完整約束方程進(jìn)行泰勒展開后與動力學(xué)方程聯(lián)立求解并以一個七桿機(jī)構(gòu)為例，驗證了算法的正確性。彭慧蓮等[9]以多體系統(tǒng)中的樹形結(jié)構(gòu)為例，詳細(xì)介紹了選用局部坐標(biāo)建立約束方程的方法和優(yōu)點，指出對具有光滑固定面約束和光滑柱鉸鏈約束的平面多剛體系統(tǒng)的關(guān)節(jié)約束反力與Lagrange 乘子一一對應(yīng)，并以曲柄滑塊機(jī)構(gòu)為例利用矢量方程法加以驗證。上述文獻(xiàn)對拉格朗日增廣法中待定乘子的引入原理進(jìn)行了詳細(xì)介紹，并從理論分析的角度詳述了機(jī)構(gòu)系統(tǒng)在受完整約束或非完整約束下關(guān)節(jié)反力的數(shù)值計算方法。另外，當(dāng)前利用增廣法求解機(jī)構(gòu)系統(tǒng)動力學(xué)的目的主要集中在對機(jī)構(gòu)運動特性的分析[10-12]，而對關(guān)節(jié)反力本身的物理含義及求解方法有效性問題研究較少，這對研究桿件壽命、關(guān)節(jié)磨損特性等問題形成一定阻礙。

針對以上問題，本文在對拉格朗日方程中約束反力與拉格朗日乘子關(guān)系的基礎(chǔ)上，指出乘子項是廣義位形空間中的廣義約束力，其實質(zhì)是廣義約束力中各個約束所占的權(quán)重。利用拉格朗日增廣法計算機(jī)構(gòu)的關(guān)節(jié)反力，詳細(xì)說明每個關(guān)節(jié)對應(yīng)關(guān)節(jié)反力的計算方法和篩選步驟。為了證明計算方法的正確性，利用牛頓歐拉方程對所選模型關(guān)節(jié)反力再次求解，并應(yīng)用ADAMS對同一模型進(jìn)行動力學(xué)仿真，將所得的三種結(jié)果進(jìn)行了對比。本文雖然以雙搖桿機(jī)構(gòu)為例，所受約束為理想約束，但由于求解過程對約束方程與結(jié)構(gòu)特性沒有特殊要求，保證了該方法對空間任意機(jī)構(gòu)系統(tǒng)關(guān)節(jié)反力求解的適應(yīng)性和準(zhǔn)確性。

1 拉格朗日乘子的含義

第一類拉格朗日方程是應(yīng)用數(shù)學(xué)分析中的乘子法建立的動力學(xué)方程。對于機(jī)構(gòu)系統(tǒng)中任意構(gòu)件其完整的表達(dá)式為[13]

(1)

(2)

式中：F為廣義外力；Qc為約束反力。比較式(1)和式(2)明顯可得約束力與拉格朗日乘子之間的關(guān)系

(3)

2 關(guān)節(jié)反力計算與篩選

約束方程是在建立系統(tǒng)模型時，系統(tǒng)的狀態(tài)變量必須滿足一定條件所構(gòu)成的方程。機(jī)械系統(tǒng)關(guān)節(jié)約束方程數(shù)等于該處關(guān)節(jié)反力個數(shù)且與該關(guān)節(jié)自由度數(shù)之和為3(平面)或為6(空間)[14]。設(shè)q為機(jī)構(gòu)系統(tǒng)廣義坐標(biāo)，C為機(jī)構(gòu)約束方程，對C關(guān)于廣義坐標(biāo)求一階和二階偏導(dǎo)可得機(jī)構(gòu)的速度和加速度約束方程

(4)

(5)

(6)

當(dāng)約束方程中包含驅(qū)動約束方程時，雅克比矩陣為方陣，且通常為非奇異矩陣，即m=n，當(dāng)約束方程中不包含驅(qū)動約束方程時，雅克比矩陣不為方陣，即m≠n。將式(6)與動力學(xué)方程式(1)聯(lián)立可得到系統(tǒng)完整的動力學(xué)模型

(7)

將剛體質(zhì)量矩陣M、雅克比矩陣Cq、合外力矢量Qe以及Qd代入式(7)即可求得拉格朗日乘子λm×1，其標(biāo)量形式可以寫為

λ=[λ1，λ2,...,λm]T

(8)

假設(shè)機(jī)構(gòu)總運動副個數(shù)為A，構(gòu)件數(shù)為b，第a個運動副自由度為Fa(a的取值規(guī)則與關(guān)節(jié)反力計算結(jié)果無關(guān))，則機(jī)構(gòu)總的約束方程個數(shù)m和機(jī)構(gòu)廣義坐標(biāo)個數(shù)n分別為

(9)

圖1 機(jī)構(gòu)系統(tǒng)運動副Fig.1 Kinematic pairs in mechanism system

同樣，與該運動副相關(guān)的拉格朗日乘子(λ)k為機(jī)構(gòu)系統(tǒng)乘子項λ中的第l+1行至第l+3-Fk行，即

(10)

將式(9)轉(zhuǎn)置與式(10)代入式(3)得

(11)

(12)

再次取所有廣義關(guān)節(jié)反力項N中與該運動副相關(guān)的關(guān)節(jié)反力，即N中的第l+1行至第l+3-Fk行

(13)

上述計算方法和篩選步驟的整個求解過程對系統(tǒng)的結(jié)構(gòu)特性并無特殊要求，適用于空間機(jī)構(gòu)系統(tǒng)，因此可以得出選用拉格朗日增廣法求解任意機(jī)構(gòu)系統(tǒng)關(guān)節(jié)反力的一般步驟：

步驟1根據(jù)機(jī)構(gòu)運動學(xué)已知條件求解約束方程C，并對約束方程關(guān)于廣義坐標(biāo)求偏導(dǎo)得到雅克比矩陣Cq，提取雅克比矩陣Cq中與待求關(guān)節(jié)相關(guān)的分量并轉(zhuǎn)置；

步驟2同樣方法，提取拉格朗日乘子λ中與代求關(guān)節(jié)相關(guān)的拉格朗日乘子分量；

圖2 關(guān)節(jié)反力計算流程圖Fig.2 Calculated flow chart of reaction force

3 數(shù)值算例

以平面四桿機(jī)構(gòu)OABC為例，簡圖如圖3所示。表1給出各均質(zhì)桿件結(jié)構(gòu)參數(shù)，轉(zhuǎn)動慣量相對質(zhì)心。由于各桿長關(guān)系滿足l1+l4≤l2+l3，因此該機(jī)構(gòu)為雙曲柄機(jī)構(gòu)，即驅(qū)動桿和從動桿均可做整周回轉(zhuǎn)。桿l2為驅(qū)動桿，逆時針勻速轉(zhuǎn)動，轉(zhuǎn)速50π rad/s，計算時間t從0～0.4 s即一個周期內(nèi)各構(gòu)件的位置、速度和加速度。取桿件質(zhì)心建立局部坐標(biāo)系，系統(tǒng)坐標(biāo)系O1XY，原點位于O，X軸與桿l1重合，方向指向C點。

圖3 四桿機(jī)構(gòu)簡圖Fig.3 Diagram of planar four-bar mechanism

表1 四桿機(jī)構(gòu)參數(shù)Tab.1 Parameters of four bar mechanism

為了證明上述拉格朗日增廣法求解各個關(guān)節(jié)對應(yīng)關(guān)節(jié)反力選擇步驟和方法的正確性，利用牛頓歐拉方程對所選模型各關(guān)節(jié)反力重新計算。設(shè)桿件質(zhì)量為m，選用質(zhì)心坐標(biāo)為廣義坐標(biāo)q=(Rx,Ry,θ)，根據(jù)桿件受到與其相鄰構(gòu)件i和j的力和力矩平衡條件，可列出如下矢量表達(dá)式

(14)

同樣選用四桿機(jī)構(gòu)OABC為對象，利用ADAMS軟件建立三維仿真模型如圖4所示。并進(jìn)行動力學(xué)仿真分析，求解各個關(guān)節(jié)接觸反力。

ADAMS簡化模型在時間t從0～0.4 s即一個周期內(nèi)進(jìn)行仿真，得到各個關(guān)節(jié)的關(guān)節(jié)反力及x向和y向的分力，將結(jié)果與前兩節(jié)的數(shù)值結(jié)果進(jìn)行對比如圖5～圖8所示。

圖4 ADAMS仿真模型Fig.4 ADAMS simulation model

圖5 O關(guān)節(jié)反力及分量Fig.5 Reaction force and components at point O

圖7 B點關(guān)節(jié)反力與分量Fig.7 Reaction force and components at point B

分析圖5～圖8中數(shù)值計算與ADAMS仿真結(jié)果曲線可得以下結(jié)論：

(1)拉格朗日增廣法對所選雙搖桿機(jī)構(gòu)各關(guān)節(jié)處關(guān)節(jié)反力計算結(jié)果與其他兩種方法計算結(jié)果相比略有誤差，但不影響本次對關(guān)節(jié)反力計算方法與篩選步驟的討論。

(2)拉格朗日增廣法計算關(guān)節(jié)反力產(chǎn)生誤差的原因主要是由于初始坐標(biāo)位置精度不同、給定計算結(jié)果誤差等級不夠以及動力學(xué)方程組求解過程中出現(xiàn)的約束違約等原因造成。

(3)三種計算方法所得各關(guān)節(jié)處關(guān)節(jié)反力變化曲線基本重合，符合四桿機(jī)構(gòu)實際的運動規(guī)律和動力學(xué)特性，可以證明文中拉格朗日增廣法求解關(guān)節(jié)反力計算方法和篩選步驟的正確性。

4 結(jié) 論

(1)針對長期以來對拉格朗日乘子法中涉及拉格朗日乘子是否為系統(tǒng)真正關(guān)節(jié)反力的問題，明確了拉格朗日乘子中的乘子項是廣義位形空間中的廣義約束力而非真正的關(guān)節(jié)約束反力。

(2)利用拉格朗日增廣法計算機(jī)構(gòu)的關(guān)節(jié)反力，詳細(xì)說明每個關(guān)節(jié)對應(yīng)關(guān)節(jié)反力的計算步驟和方法，并以雙搖桿機(jī)構(gòu)為例對一個周期內(nèi)的所有關(guān)節(jié)的關(guān)節(jié)反力進(jìn)行數(shù)值求解。

圖8 C點關(guān)節(jié)反力及分量Fig.8 Raction force and components at point C

(3)為了驗證計算方法和選擇步驟的正確性，利用牛頓歐拉方程與ADAMS軟件對同一模型進(jìn)行動力學(xué)數(shù)值仿真，通過使用不同方法求得的各關(guān)節(jié)反力相同，且誤差在允許范圍之內(nèi)。

(4)關(guān)節(jié)反力計算方法與篩選步驟雖然以雙搖桿機(jī)構(gòu)運動一周各個關(guān)節(jié)所受關(guān)節(jié)反力的計算與篩選為例，但通過第二部分論述，該方法可使用于空間任意機(jī)構(gòu)系統(tǒng)的關(guān)節(jié)反力計算。