Hamilton體系的辛算法

王顏超，鄧亞莉

（內蒙古大學，內蒙古呼和浩特010021）

Hamilton體系的辛算法

王顏超，鄧亞莉

（內蒙古大學，內蒙古呼和浩特010021）

本文介紹了辛空間的基本理論,通過構造辛算法得到Hamilton方程的數值解，并與分離變量法得到的解進行了比較.

Hamilton體系；辛算法；辛空間

在自然界中很多數學物理問題均可以用Lagrange方程或Hamilton方程表示,而Lagrange方程可以通過Legendre變換為Hamilton方程.經典力學具有三種表示形式,設有n個自由度的運動,位置向量記為q=(q1,q2，…,qn)T,勢能函數為V=V(q),則是運動方程的標準形式,它是在n維位形空間Rn中的2階微分方程組,通常稱為經典力學的標準形式即牛頓形式.Lagrange通過引進動能,勢能差數的作用量并利用變分原理將運動方程寫為它被稱為經典力學的變分形式即Lagrange形式.后來Hamilton利用動量p=Mq˙和總能量H=T+V將運動方程寫為稱為Hamilton正則方程，它是2n維相空間中相變量(p1,…,pn；q1,…,qn)的1階微分方程組.

作為動力系統的重要體系,Hamilton體系可以用來表示一切真實的,耗散可忽略不計的物理過程,而Hamilton體系的基礎是辛幾何.Hamilton將系統能量以廣義坐標及廣義動量來表示,構成一個2n維相空間,故經典力學可以被看成相空間中的一種幾何學即辛幾何學.經典力學基本定理用辛幾何的語言就表示為“一切哈氏體系的動力演化都使辛度量保持不變,即就是辛(正則)變換”.因此解哈氏方程的“正確”的離散算法就應是辛變換,這種算法稱為辛(正則)算法或哈密爾頓算法.

1 辛空間基本理論

通常我們比較熟悉的是以實數為坐標的歐式空間Rn,現將Rn推廣到以復數為坐標的酉空間,Cn為n維復向量空間,C表示全體復數,Z=(z1,z2,…,zn),zk∈C,k=1,2,…,n.在Cn中引進Hermite內積

具有Hermite內積的復線行空間稱為酉空間.在式(1.1)中,因zk∈C,故可以寫為zk=xk+iξk,xk,ξk∈R,因而Cn中的向量與Rn中的向量建立一一對應關系如下:

Cn?R2n

且對向量加法及用實數相乘仍保持,即:

可以看出Hermite內積(Z,W)由兩部分組成,實部是R2n中歐式內積,虛部可以看成是R2n中兩個向量x=(x1,x2,…,xn, ξ1,ξ2,…,ξn)與y=(y1,y2,…,yn,η1,η2,…,ηn)的一種新的“內積”,我們稱這種內積為辛內積[1].

定義1設W是實數域R上的一個2n維相空間,對W中的任意兩個向量α,β依一定法則對應著一個實數,這個數稱為辛內積,記作〈α,β〉,且辛內積〈α,β〉運算滿足下列4個條件:

（1）〈α,β〉=-〈β，α〉（1.2）

（2）〈kα,β〉=k〈α,β〉,?k∈R（1.3）

（3）〈α+γ,β〉=〈α,β〉+〈γ,β〉,γ為W中任意向量（1.4）

（4）若向量α對W中任意向量β均有〈α,β〉=0,則α=0（1.5）

稱定義有這種辛內積的相空間為辛空間,辛空間與研究長度等度量性質的歐幾里得空間不同,它研究的是面積,或說是研究做功的,利用它可以描述許多重要的物理現象,故近年來受到很到理論與實際工作者的重視.

R2n中兩個向量x=(x1,x2,…,xn,ξ1,ξ2,…,ξn)與y=(y1,y2,…,yn, η1,η2,…,ηn)定義辛內積,即上述Hermite內積的虛部〈x,y〉=

定義2若向量α,β的辛內積〈α,β〉=0,則稱α與β辛正交,記為α⊥β;否則稱α與β辛共軛.與歐式空間不同的是?α均有〈α，α〉,則α⊥α,即辛空間中任一向量必與自身正交.

若向量組(α1,α2,…,αr,β1,β2,…,βr)(r≤n)滿足:

定理1[2]共軛辛正交向量組是線性無關向量組.

辛矩陣具有以下性質:

(1)其逆矩陣為辛矩陣;(2)其轉置矩陣為辛矩陣;

(3)其行列式值為1或-1;(4)其乘積為辛矩陣.

定義4[3]一個2n階矩陣B稱為無窮小辛陣,如果滿足JB+BTJ=0.若B為無窮小辛陣,則eB為辛陣;若C為對稱陣,當且僅當B=JC時,則B為無窮小辛陣.

定義5若2n×2n矩陣H對任意2n維向量x,y有:＜x, Hy＞=＜y,Hx＞,則稱H為Hamilton矩陣.易證Hamilton矩陣可等價定義為:(JH)T=JH或JHJ=HT.

定理2若u為Hamilton矩陣H的本征值,重數為m,則-u也是其本征值,重數也是m;若Hamilton矩陣H具有零本征值,則其重數為偶數.稱本征值±u為Hamilton矩陣的互為辛共軛本征值,零本征值是特殊的辛本征值,與其自身互為辛共軛.

1.2.1.2 成立在職培訓管理組:由區護士長任組長,指定護理組長一對一帶教;根據CSSD制定的崗前培訓計劃,專人負責,落實到位。

定理4設±u≠0為Hamilton矩陣H的本征值,重數為m，則必存在共軛辛正交的向量組(ψ(0),ψ(1),…,ψ(m-1),φ(m-1),…, φ(1),φ(0)),即:

其中(ψ(0),ψ(1),…,ψ(m-1))和(φ(0),φ(1),…,φ(m-1))分別為u和-u對應的基本本征向量和Jordan型本征向量.

2 辛算法的構造

下面我們考慮一般的線性Hamilton體系,p=[p1,p2,…,pn]T, q=[q1,q2,…,qn]T,則,設Hamilton函數是H(z1,其中Z=[p,q]T=[z1,z2,…,z2n],C為

2n階對稱矩陣,即CT=C,則

其中B=J-1C為無窮小辛陣,用分離變量法求得方程Z˙ =BZ的解為Z(0),求數值解時需要進行離散化,對變量t取步長為τ＞0,將連續依賴t的Z(t)離散化為{Z(kτ)},k=0,1,2,…則eBτZ(kτ),k=0,1,2,…下面我們利用有理Padé逼近構造eBτ保持典則變換的近似表達式,即用有理分式來逼近eBτ,其中Pm(x),Qn(x)分別為m,n次多項式.

定理5設B為2n階矩陣,滿足JB+BTJ=0,即B為無窮小辛陣,若的Padé逼近,則是算子etB保持為辛算法的近似表達式.

則2階精度的辛差分格式(隱式)為:Z((k+1)τ)=Z(kτ)+Bτ 2 (Z((k+1)τ)+Z(kτ)).

下面討論用分離變量法求解:Z˙=BZ,設Z(t)=ξ(t)ψ,其中ψ= (ψ1,ψ2,…,ψ2n)T與t無關,ξ(t)是t的函數,與向量ψ中的任何分量無關.則:Z˙=ξ˙(t)ψ=Bξ（t）ψ,故:得:Bψ=μψ,ξ(t) =eμt,而Bψ=μψ為Hamilton矩陣的本征問題.

于是:(B-μI)ψ=0,det(B-μI)=f(μ)=0,

故:f(H)=μ2nH2n+…+μ0I2n=0得:μi:ψ1,…,ψni

3 一個例子

易證:JHJ=HT,H為Hamilton矩陣.

且通過特征多項式求得其特征值為1重根,由分離變量法理論得知僅有上(2.1)式存在.

下面是辛算法(右圖)和分離變量法(左圖)求解x方向的結果:

分離變量X·方向相圖

平面簡諧振子X·方向相圖

我們模擬了y方向步長為1的相圖.此時的辛算法求得的結果與精確結果相差比較大,由于篇幅的原因這里就不再贅述了(有興趣的讀者可以參考我們的matlab程序).

function[ZZ]=Z(k,t,beta,z)

m=length(beta);

A=eye(m)-(t.*beta./2);

ZZ(:,1)=z;

fori=1:k

z=A(eye(m)+(t.*beta./2))*z;

ZZ(:,i+1)=z;

end

function[QQ]=FLBL(G,D,V)

fort=0:0.1:100

k=t*10+1;k=int64(k);

QQ(k,:)=G*[exp(D(1,1)*t)000;

0exp(D(2,2)*t)00;

00exp(D(3,3)*t)0;

000exp(D(4,4)*t)]*V;

end

下圖為分離變量法和辛算法求解y方向的結果(此圖的步長為0.1):

分離變量與辛算法組合圖

分離變量與辛算法比較圖

從圖上可以得出:用辛算法得到的解和用分離變量法得到的解幾乎一致.

辛算法保持了Hamilton體系的兩個守恒[4]:

（1）相空間體積的不變——Liouville-Poincaré守恒律?

（2）運動不變量:如能量,動量,角動量的守恒

辛算法能夠在數值計算中保持辛變換的結構,于是就會得到高的穩定性,辛算法的差分方法被認為是目前最穩定,高效的計算方法,適合用于經典力學體系.

我們知道,辛算法之所以逼近真實解的程度這么高是因為它的每一步迭代過程都保持了辛變換即典則變換,也就是說每次迭代都是能量守恒的.但是,19世紀龐加萊等人指出三體問題不可積,并意識到許多Hamilton系統是不可積的.其實并不是所有系統都是能量守恒的,耗散結構以及其他的動力系統都是實實在在存在的.但即便如此,Hamilton系統依舊是研究其他問題的重要體系.上世紀30年代,偉大的量子力學創始人之一Schr?dinger先生就曾經說過“Hamilton原理已經成為現代物理學的基石”.所以,研究辛空間以及辛算法都是非常有價值的.

〔1〕李世雄.波動方程的高頻近似與辛幾何[M].北京:科學出版社，2001.3.

〔2〕姚偉岸,鐘萬勰.辛彈性力學[M].北京:高等教育出版社. 2002.4.

〔3〕秦孟兆.辛幾何及計算哈密頓力學[J].力學與實踐，1990（2）.

〔4〕馮康,秦孟兆.哈密兒頓系統的辛幾何算法[M].杭州：浙江科學技術出版社，2003.

O241

A

1673-260X（2013）06-0006-03