高階Maggi方程的Birkhoff化及其辛算法*

薛冰 解加芳 張可心

(北方工業(yè)大學(xué) 理學(xué)院, 北京 100144)

引言

非完整系統(tǒng)是一類(lèi)受到不可積微分約束的動(dòng)力學(xué)系統(tǒng)[1],廣泛應(yīng)用于場(chǎng)論、機(jī)電動(dòng)力系統(tǒng)、控制理論、工程科學(xué)等領(lǐng)域[2].20世紀(jì)80年代我國(guó)學(xué)者梅鳳翔研究了帶參數(shù)約束的一類(lèi)可控系統(tǒng)[3]、變質(zhì)量非完整約束系統(tǒng).Maggi在1896年推廣了拉格朗日第二類(lèi)方程[4],對(duì)線(xiàn)性非完整約束系統(tǒng)得到一類(lèi)動(dòng)力學(xué)方程,后人稱(chēng)為Maggi方程[5],這些方程后來(lái)被推廣到非線(xiàn)性非完整系統(tǒng)[6].Maggi方程是力學(xué)系統(tǒng)[7]各大運(yùn)動(dòng)方程的中間產(chǎn)物,對(duì)研究非完整系統(tǒng)的運(yùn)動(dòng)具有重要意義.

Birkhoff動(dòng)力學(xué)理論是Hamilton動(dòng)力學(xué)的自然推廣,它是包括齊次Hamilton系統(tǒng)和非齊次Hamilton系統(tǒng)的更一般動(dòng)力學(xué)理論,是最一般辛結(jié)構(gòu)的局部實(shí)現(xiàn),只有Birkhoff系統(tǒng)與一般辛幾何結(jié)構(gòu)之間才有一一對(duì)應(yīng)關(guān)系.因此 Birkhoff 系統(tǒng)動(dòng)力學(xué)[8]的研究對(duì)于完善和深化分析力學(xué)的理論體系具有重要意義, 尤其是對(duì)于非齊次 Hamilton 動(dòng)力學(xué)系統(tǒng)的幾何結(jié)構(gòu)分析具有重要應(yīng)用價(jià)值[9].本文針對(duì)非完整系統(tǒng)高階Maggi方程,在其滿(mǎn)足一定條件下,將其進(jìn)行Birkhoff化。并通過(guò)一個(gè)算例驗(yàn)證上述理論分析的正確性,再分別采用Runge-Kutta方法和Birkhoff辛算法對(duì)其進(jìn)行數(shù)值計(jì)算[10],并將數(shù)值結(jié)果進(jìn)行比較,給出Birkhoff辛算法在長(zhǎng)時(shí)計(jì)算時(shí)的優(yōu)越性.

1 高階非完整系統(tǒng)Maggi方程的Birkhoff化

設(shè)力學(xué)系統(tǒng)的位形由n個(gè)廣義坐標(biāo)qs(s=1,…,n)確定,系統(tǒng)受有g(shù)個(gè)理想m階非完整約束[11]

(1)

其中

根據(jù)d’Alembert-Lagrange原理可以導(dǎo)出Maggi形式為:

(σ=1,…,ε)

(2)

式中Qσ為廣義力,T為系統(tǒng)動(dòng)能.

令

(ν=1,…,ε;k=1,...,n)

(3)

則方程(2)有形式

(4)

(5)

約束對(duì)初始條件的限制為:

(6)

(7)

約束對(duì)初始條件的限制為:

(8)

將式(7)化為標(biāo)準(zhǔn)一階形式,令

(s=1,…,n)

(9)

則式(7)可以寫(xiě)成形式

(10)

其中

aν=xν,σs=xn+s,…,

σ(m-2)n+s=x(m-1)n+s,

σ(m-1)n+s=hs(s=1,…,n)

(11)

為將式(10)表為Birkhoff形式,其階必為偶數(shù)[9],即mn=2N,如果mn為奇數(shù)2N-1,可增加一個(gè)方程

(12)

使其成為偶階.從而,要使高階非完整系統(tǒng)的Maggi方程可表示成Birkhoff形式

(μ=1,…,2N)

(13)

即要求滿(mǎn)足

(μ=1,…,2N)

(14)

設(shè)式(10)的2n個(gè)第一積分Iμ(t,a)彼此無(wú)關(guān),即

(μ=1,…,2n)

(15)

(16)

根據(jù)Hojman方法,則式(14)的Birkhoff函數(shù)組Rμ由下式確定

(17)

Birkhoff量B為:

(18)

其中Gα需滿(mǎn)足條件

(19)

2 Maggi方程的廣義辛差分格式

(20)

其中

(21)

假定方程組(20)的一種離散可以記為:

[?tR(zi,ti)]i

(22)

式(22)稱(chēng)為式(20)的一個(gè)離散格式,如果它決定的格式Φ保持離散的K(z,t)辛格式,即

(23)

(24)

以及

(25)

此映射為:

(26)

式(26)的Jacobi矩陣為:

(27)

映射α要滿(mǎn)足

(28)

其中

當(dāng)K與z無(wú)關(guān)時(shí),根據(jù)式(27)和式(28)可以得到

(29)

且需滿(mǎn)足下面的截面條件

|CαM+Dα|≠0

設(shè)Bα=Cα=0,則式(29)可成

(30)

由式(30)可以得到

(31)

生成函數(shù)φ(ω,t,t0)=α(t,t0),可以構(gòu)造廣義Birkhoff辛差分格式.當(dāng)步長(zhǎng)τ>0足夠小的時(shí)候,取

(m=1,2,…)

(32)

α1(zk+1,zk,tk+1,tk)

(33)

3 算例

假設(shè)某一力學(xué)系統(tǒng)的位形由2個(gè)廣義坐標(biāo)q1,q2確定,其系統(tǒng)動(dòng)能為

(34)

且該系統(tǒng)受到1個(gè)2階非完整約束:

(35)

(36)

顯然式(36)有解:

q1=1-arctan(t)

q2=-3t·arctan(t)+2ln(t2+1)+t

由此可以構(gòu)造Birkhoff函數(shù)Rμ,B如下[11]:

R2=0

R4=0

從而將該系統(tǒng)的Maggi方程可以Birkhoff化:

結(jié)合式(24)和式(25),我們可以確定生成函數(shù)φ(ω,t,t0),之后采用上述的Birkhoff廣義辛算法式(32), 得到此 Birkhoff 系統(tǒng)的二階K(z,t)離散格式[11-13].

對(duì)該題采用二階K(z,t)算法和二階Runge-Kutta算法進(jìn)行計(jì)算.在計(jì)算過(guò)程中,先取如下初值:q1=1,C1=C2=1,取步長(zhǎng)τ=0.01.并通過(guò)比較兩種數(shù)值方法計(jì)算所得數(shù)值解和解析解q1=1-arctan(t)之間的相對(duì)誤差來(lái)說(shuō)明兩種數(shù)值方法的差別.

對(duì)比圖1和圖2可以看出,Runge-Kutta方法在長(zhǎng)期跟蹤后與解析解有著大幅度的相對(duì)誤差,而B(niǎo)irkhoff辛算法則相對(duì)誤差非常的小.因此,Birkhoff辛算法結(jié)果更加精確.

圖1 二階Runge-Kutta算法相對(duì)誤差Fig.1 Relative error of second order Runge-Kutta algorithm

圖2 Birkhoff辛算法相對(duì)誤差Fig.2 Relative error of Birkhoff symplectic algorithm

4 結(jié)論

本文對(duì)一定條件下的非完整系統(tǒng)的高階Maggi方程(2)先進(jìn)行了Birkhoff化,得到廣義Birkhoff方程(13),并針對(duì)該方程,應(yīng)用Birkhoff廣義辛差分格式與傳統(tǒng)Runge-Kutta算法分別進(jìn)行計(jì)算,比較兩種算法,最后得出Birkhoff廣義辛差分算法在求解非完整系統(tǒng)高階Maggi方程中更加優(yōu)越.