戴余良，王長湖，苗 海，劉祖源

（1海軍工程大學 船舶與動力學院，武漢 430033；2青島92196部隊裝備部，山東 青島 266011；3武漢理工大學 交通學院，武漢 430063）

潛艇水下運動穩(wěn)定性非線性分析研究

戴余良1，王長湖2，苗 海1，劉祖源3

（1海軍工程大學 船舶與動力學院，武漢 430033；2青島92196部隊裝備部，山東 青島 266011；3武漢理工大學 交通學院，武漢 430063）

應用非線性系統運動穩(wěn)定性理論和同倫延拓數值計算方法，分析了潛艇水下運動的穩(wěn)態(tài)響應及運動穩(wěn)定性，潛艇應急上浮運動穩(wěn)定性實例分析結果通過其大攻角六自由度運動方程數值積分的動態(tài)仿真進行了驗證，表明穩(wěn)定性分析結果是正確的。從而說明了采用同倫延拓法分析潛艇水下運動穩(wěn)定性是有效的。

潛艇；運動穩(wěn)定性；非線性動力學；同倫延拓法

1 引 言

潛艇是主要在水下活動的航行體。為了保證潛艇水下航行安全，特別是在艙室破損進水、舵卡等事故狀況下的安全，清楚地了解潛艇的動態(tài)響應特性和穩(wěn)定性是至關重要的。

上述穩(wěn)定性評價指標基于如下假設：在小擾動基準直線運動中，水平面運動與垂直面運動之間的耦合相對較弱，可以忽略。但是，對于高速航行的潛艇，在航行的極限狀態(tài)或應急狀態(tài)下，上述關于水平面與垂直面之間運動不耦合的假設不再成立，此時在所有六個自由度都可能發(fā)生大幅運動，各種運動模態(tài)之間的非線性相互作用變得非常顯著。因此，簡單地使用指標Gv和Gh可能導致不正確的結論。此時必須考慮水平面與垂直面運動之間相互耦合的運動特性。

潛艇的水下運動，尤其是應急情況下的操縱運動通常是非線性的。隨著非線性理論的日益成熟，運用非線性動力學理論對潛艇的運動微分方程進行分析，能夠更準確地對潛艇的運動規(guī)律作出描述，從而得到更準確的結論。本文運用非線性系統運動穩(wěn)定性理論，對潛艇水下運動穩(wěn)定性的分析方法進行了深入研究。

2 非線性系統運動穩(wěn)定性

運動穩(wěn)定性理論自Lyapunov于十九世紀九十年代創(chuàng)建以來，在物理科學和工程技術等各個領域都獲得了廣泛的應用。

Lyapunov意義下的運動穩(wěn)定性理論是研究干擾因素對物體運動的影響。眾所周知，微小的干擾因素對物體運動的影響，對于不同的運動是不一樣的。對于一些運動，這種影響可能并不顯著，因而受擾運動與未擾運動相差很??；反之，對于另外一些運動，干擾的影響可能很顯著，以至于干擾力無論多么小，受擾運動與未擾運動隨著時間的推移而可能相差得很大。簡單地說，屬于前者的運動稱為穩(wěn)定的，而屬于后一類型的運動則稱為不穩(wěn)定的。由于在現實世界中，干擾因素總是不可避免地存在著，所以運動穩(wěn)定性的問題有很重要的理論意義和實際意義，這也正是穩(wěn)定性理論蓬勃發(fā)展的重要原因。

2.1 運動穩(wěn)定性判定方法

在研究工程問題時，系統的漸近穩(wěn)定性具有特別重要的意義。那么如何判斷所考察的系統是穩(wěn)定的還是不穩(wěn)定的？Lyapunov直接法和Lyapunov一次近似法是兩種最基本的分析方法。由于Lyapunov直接法的成功與否取決于構造合適的Lyapunov函數。Lyapunov函數的構造需要一定的經驗和技巧，目前尚無一般規(guī)律。因此，Lyapunov一次近似法被廣泛應用于非線性系統的穩(wěn)定性分析。

Lyapunov一次近似法的基本思想[2-3]是：將非線性系統線性化，僅保留一次項，研究該線性系統的穩(wěn)定性。（1）如果一次線性近似系統的所有特征值都具有負實部，則原非線性系統的零解（平衡點）是漸近穩(wěn)定的；（2）如果一次線性近似系統的特征值中只要有一個具有正實部，則原非線性系統的零解是不穩(wěn)定的；（3）如果一次線性近似系統存在實部為零的特征值，而其余所有特征值都具有負實部，此時非線性系統零解的穩(wěn)定性取決于系統的非線性項，不能用Lyapunov一次近似法來研究，它的判定方法十分復雜和困難，本文試圖利用奇異性和分叉理論來分析。

2.2 運動穩(wěn)定度

物體的運動穩(wěn)定性通常均可通過分析系統方程的特征值λ的性質來確定，當特征值的實部Re(λ)＜0時，物體的擾動運動逐漸衰減，且Re(λ)越負，衰減越快，系統越穩(wěn)定；反之，當Re(λ )＞0時，物體的擾動運動越來越強，且Re(λ)越大，系統越不穩(wěn)定；如果存在一對共軛復特征值，則物體運動會發(fā)生振蕩。因此，系統特征值λ實部的正負決定了系統的運動穩(wěn)定性，而λ實部的大小決定了物體擾動運動衰減的快慢，如果λ存在一對共軛復數，則物體擾動運動周期性振蕩。

本文定義系統所有特征值的最大實部為系統的穩(wěn)定度 （degree of stability），記為e，即e=max Re λi()｛｝i=1,2,…,n，其中λi為系統的特征值，則e越負系統越穩(wěn)定，e越正系統越不穩(wěn)定。

3 同倫延拓非線性數值分析方法

同倫延拓方法是求解非線性問題的一種非常有效的方法，它克服了傳統迭代法局部收斂的弱點，對初值的選取沒有嚴格限制，能夠全局收斂，可以求出問題的全部解。它雖然于二十世紀70年代才開始發(fā)展起來，但在短短的30多年里，同倫延拓方法已經獲得了廣泛應用[4-9]。傳統的求解方法如迭代法、消去法等存在的最大問題是依賴于初值的選取，初值選擇不當常常導致迭代過程不收斂，而且很難求出全部解。因此，同倫延拓方法被譽為二十世紀數學研究中一項突破性的新成果[4]。

3.1 同倫延拓法的基本思想

同倫延拓法的基本思想[4-9]是借助同倫函數H，從輔助映射g的零點集出發(fā)，跟蹤同倫方程組H(x, λ )=0的解曲線，走到目標映射f的零點集。

為了求解非線性方程組f(x)=0，其中光滑映射 f：Rn→Rn，x=｛x1,x2,…,xn｝T，選取具有已知零點的光滑輔助映射 g：Rn→Rn，將其與映射f構成一線性同倫函數H：Rn×Rn→Rn，使得H( x, 0 )=g(x)和H( x, 1 )=f(x)，典型的同倫函數為：

對于同倫方程組H( x, λ )=0，當λ=0時，其解與輔助方程組g(x)=0的解相同；讓同倫參數λ從0逐漸變化到1，并跟蹤同倫方程組的解，則當λ=1時，同倫方程組的解便是待解方程組f(x)=0的解。

3.2 同倫延拓算法

根據同倫延拓法的基本思想，利用同倫延拓法求解非線性方程組f(x)=0首先應該選擇合適的輔助函數g(x)。g(x)的選取方法很多，本文選擇g(x)=f(x)-f( x( 0 ))。那么同倫函數為：

因此，H( x, 0 )=f(x)-f( x( 0 ))和H( x,1)=f(x)，其中x(0)是系統f(x)=0對應于λ=0時的解，x(1)是系統f(x)=0對應于λ=1時的解，x(λ)是同倫方程組H( x( λ ),λ )=0的唯一解。

利用同倫延拓法，非線性方程組f(x)=0的求解問題轉化成了求解微分方程：

其中x(0)是系統f(x)=0的初始已知解。

微分方程（1）的求解可以采用任一種數值計算方法，本文采用4階Runge-Kutta法求解，則4階Runge-Kutta法的迭代公式為：

圖1 同倫延拓算法流程圖Fig.1 Flow of homotopy and continuation algorithm

式中：h為計算步長。同倫參數λ從0逐漸變化到1，則當λ=1時的解x(1)即是系統f(x)=0的近似解。

為了避免由于Jacobi矩陣J(x)不可逆時上述Runge-Kutta算法失效，將計算ki的等式（2）改寫成求解如下線性方程：

具體算法流程見圖1。

4 潛艇水下運動非線性方程穩(wěn)態(tài)解

4.1 潛艇水下運動穩(wěn)態(tài)方程

為了研究潛艇運動的非線性特性，必須采用大攻角六自由度非線性運動方程，本文采用文獻[10-11]的潛艇運動方程，可化成一般形式：

式中狀態(tài)變量 x=[u,v,w,p,q,r,φ,θ,ψ ]，控制變量 u= [δr,δb,δs,δB,xGB,yGB,zGB,n ]，其中 δB=B-W，xGB=xG-xB，yGB=yG-yB，zGB=zG-zB，n 為螺旋槳轉速。

4.1.1 穩(wěn)態(tài)條件

潛艇水下運動的穩(wěn)態(tài)條件是：處于定常運動狀態(tài)（定常直航或定?；剞D狀態(tài)），即

4.1.2 穩(wěn)態(tài)運動方程

將穩(wěn)態(tài)條件（6）式代入方程（5）便可得到潛艇水下穩(wěn)態(tài)運動方程：

4.2 潛艇水下運動方程穩(wěn)態(tài)解

由于穩(wěn)態(tài)方程（7）是高度非線性的，其解可能出現分叉和多解現象，因此本文使用同倫延拓法求解穩(wěn)態(tài)方程（7）的可能解。

潛艇在水下運動過程中，潛浮速度Vζ、橫傾角φ、縱傾角θ和攻角α是影響潛艇安全的主要參數，所以重點分析穩(wěn)態(tài)的潛浮速度Vζ、橫傾角φ、縱傾角θ和攻角α隨控制參數u= [δr,δb,δs,δB,xGB,yGB,zGB,n ]的變化情況。

潛艇水下運動穩(wěn)態(tài)方程求解及運動穩(wěn)定性分析過程：

（1）通過動態(tài)方程（5）采用數值積分法求得某一開始穩(wěn)定狀態(tài)的狀態(tài)變量值，以此作為同倫延拓算法的初值；

（2）以某一控制參數作為同倫參數，其它控制參數的值保持不變，計算出各狀態(tài)變量隨同倫參數變化的穩(wěn)態(tài)解曲線及其Jacobi矩陣的特征值；

（3）改變除同倫參數以外的其它控制參數（作為擾動參數）的值，重復上面的計算，得出另一組穩(wěn)態(tài)解曲線。

（4）根據特征值的性質確定平衡解的穩(wěn)定性。

4.3 潛艇水下運動穩(wěn)定性

潛艇水下運動方程可寫成如下形式：

式中狀態(tài)變量 x=[u,v,w,p,q,r,φ,θ,ψ ]。

為了分析潛艇水下運動在某個平衡狀態(tài)x0的穩(wěn)定性，根據Lyapunov線性化原理，將水下運動方程（8）在平衡點 x0附近線性化，即：

其中A為f(x)在x0處的Jacobi矩陣，即

由運動穩(wěn)定性理論可知，若矩陣A的所有特征值有負實部，則平衡狀態(tài)x0是漸近穩(wěn)定的；只要矩陣A有一個特征值有正實部，則平衡狀態(tài)x0是不穩(wěn)定的。因此，當系統的穩(wěn)定度e＜0時，則平衡狀態(tài)x0是漸近穩(wěn)定的；當系統的穩(wěn)定度e＞0時，則平衡狀態(tài)x0是不穩(wěn)定的。

5 潛艇水下運動穩(wěn)定性分析實例

應急上浮是潛艇在艙室破損進水、舵卡等事故狀況下所采取的緊急操縱措施，是潛艇非戰(zhàn)時緊急脫險的最有效方法。為了確保采取的操縱方法正確有效，保證潛艇應急上浮安全，清楚地了解潛艇在應急狀況下的動態(tài)響應特性和穩(wěn)定性至關重要。潛艇應急上浮通常是一種大攻角強機動，描述潛艇的運動模型是強非線性的。下面針對潛艇應急上浮運動，以xGB為同倫參數分析潛艇運動穩(wěn)定性為例，說明本文提出的潛艇水下運動穩(wěn)定性非線性分析方法。潛艇模型參數及水動力系數采用文獻[10-11]提供的數據，其中艇重53.4kN，艇長5.3m。

5.1 穩(wěn)態(tài)解及特征值

計算條件：δr=δs=δb=0，n=1 000rpm，yGB=0，zGB=0.061m，以 δB為擾動參數，分別取 δB＝1%W，5%W 和10%W，其中W為艇重，xGB的變化范圍為：-10%L～10%L，其中L為艇長。

計算結果如圖2～10所示。

5.2 結果與討論

由圖2～10可以看出：

圖2 對于不同δB航速U隨xGB的變化 Fig.2 Speed versus xGBfor variations in δB

圖3 對于不同δB橫傾角φ隨xGB的變化Fig.3 Roll angle versus xGBfor variations in δB

圖4 對于不同δB縱傾角θ隨xGB的變化Fig.4 Pitch angle versus xGBfor variations in δB

圖5 對于不同δB攻角α隨xGB的變化Fig.5 Attack angle versus xGBfor variations in δB

圖6 對于不同δB深度速率ζ˙隨xGB的變化 Fig.6 Speed depth versus xGBfor variations in δB

圖7 對于不同δB穩(wěn)定度隨xGB的變化Fig.7 Degree of stability versus xGBfor variations in δB

圖8 δB=1%W時系統特征值隨xGB的變化Fig.8 Eigenvalue versus xGBfor δB=1%W

圖9 δB=5%W時系統特征值隨xGB的變化Fig.9 Eigenvalue versus xGBfor δB=5%W

圖10 δB=10%W時系統特征值隨xGB的變化Fig.10 Eigenvalue versus xGBfor δB=10%W

（1）潛艇的航速、上浮速度和縱、橫傾角和攻角等的穩(wěn)態(tài)響應在正浮力的作用下，在xGB=2%L附近均有一個突變值，由此看到了當δr=0時（即不操方向舵），潛艇在正浮力作用下出現了突傾現象（見圖3），即出現了平面外運動的情況。在xGB的其它取值范圍內潛艇運動的變化相對平穩(wěn)。

（2）由圖7可知，當xGB=2%L左右時潛艇在正浮力作用下上浮速度、縱、橫傾發(fā)生突跳現象，主要是由于此時系統的穩(wěn)定度是一個較大的正數值。同時由圖8～10可以看出，系統特征值中出現了幅值較大的共軛虛數，此時運動將發(fā)生大幅振蕩。

（3）圖7系統穩(wěn)定度表明，潛艇在正浮力作用下當xGB≤0時，穩(wěn)定度均為負值，其運動是穩(wěn)定的，同時由圖8～10可以看出，系統特征值存在一對幅值較小的共軛復數，此時運動將發(fā)生輕微振蕩。xGB＞0時穩(wěn)定度基本上為正值，其運動是不穩(wěn)定的。

因此，當xGB≤0時潛艇在正浮力作用下，穩(wěn)定度均為負值，其運動是穩(wěn)定的，當xGB＞0時，穩(wěn)定度過零點，由負值變成正值，運動由穩(wěn)定變?yōu)椴环€(wěn)定，出現分叉點，同時由圖8～10可以看出，在xGB=0前后系統特征值均為實數，此時為靜態(tài)分叉。

6 數值仿真

為了驗證前面穩(wěn)定性分析實例結果的正確性，本節(jié)采用穩(wěn)態(tài)方程（7）所對應的文獻[10-11]提供的潛艇大攻角六自由度運動方程，潛艇模型參數及水動力系數與前面穩(wěn)定性分析實例相同，采用文獻[10-11]提供的數據，其中艇重53.4kN，艇長5.3m。利用Matlab通過龍格庫塔（Runge－Kutta）數值積分進行了仿真，將仿真結果與前述分析結果進行對比。

6.1 仿真結果

仿真條件：δr=δs=δb=0，n=1 000rpm，yGB=0，zGB=0.061m，δB＝1%W， 其中 W 為艇重，xGB=-5%L 和 2%L，其中L為艇長。

仿真結果如圖11所示。

6.2 結果討論

由圖11的仿真結果可以看出：

（1）當xGB=-5%L時潛艇縱傾角θ動態(tài)響應的穩(wěn)態(tài)解，由圖7可知該穩(wěn)態(tài)解是穩(wěn)定的。

（2）當xGB=2%L時潛艇縱傾角θ動態(tài)響應的穩(wěn)態(tài)解，由圖7、8可知該穩(wěn)態(tài)解是不穩(wěn)定的，此時潛艇的縱傾響應是一種周期性的振蕩運動。

因此，由上面的分析可知，通過數值積分前面穩(wěn)定性分析結果的正確性得到了初步驗證，它與數值仿真的結果完全一致，說明了前述穩(wěn)定性分析結果是正確的。

圖11 xGB=-5%L和2%L，δB=1%W時潛艇縱傾的動態(tài)響應Fig.11 Pitch response for xGB=-5%L and 2%L， δB=1%W

7 結 論

本文對非線性系統運動穩(wěn)定性理論及其數值分析方法進行了研究，提出了將同倫延拓非線性數值計算用于潛艇水下運動的穩(wěn)態(tài)響應及運動穩(wěn)定性分析的新方法，潛艇應急上浮運動穩(wěn)定性實例分析結果通過其大攻角六自由度運動方程數值積分的動態(tài)仿真進行了驗證，表明穩(wěn)定性分析結果是正確的。從而說明了本文提出的潛艇水下運動穩(wěn)定性非線性分析方法是正確有效的。

[1]施生達.潛艇操縱性[M].北京:國防工業(yè)出版社,1995.

[2]阮 炯,顧凡及,蔡志杰.神經動力學模型方法和應用[M].北京:科學出版社,2002.

[3]陸啟韶.常微分方程的定性方法和分叉[M].北京:北京航空航天大學出版社,1989.

[4]王則柯,高堂安.同倫方法引論[M].重慶:重慶出版社,1990.

[5]Lee Jaewook,Chiang Hsiao-Dong.A singular fixed-point homotopy method to locate the closest unstable equilibrium point for transient stability region estimate[J].IEEE Transactions Circuits and Systems-II:Express Briefs,2004,51(4):185-189.

[6]Burden R L,Douglas Faires J.Numerical Analysis[M].Seventh Edition.Thomson Learning,2000.

[7]陳 永,李 立等.一種基于同倫函數的迭代法─同倫迭代法[J].數值計算與計算機應用,2001(2):149-155.

[8]周 洲,祝小平,張肇熾等.一種新的同倫解曲線跟蹤算法[J].西北工業(yè)大學學報,1998,16(2):246-250.

[9]黃象鼎,曾鐘鋼,馬亞南.非線性數值分析[M].武漢:武漢大學出版社,2000.

[10]Ibrahim Aydin.Out of plane solutions of submarines in free positive buoyant ascent[R].AD-A271544,1993.

[11]Healey A J,Lienard D.Multivariable sliding mode control for autonomous diving and steering of unmanned underwater vehicles[J].IEEE Journal of Ocean Engineering,1993,18(3):327-339.

Study on nonlinear analysis of motion stability of submarines under water

DAI Yu-liang1,WANG Chang-hu2,MIAO Hai1,LIU Zu-yuan3
(1 Naval Architecture and Power College,Naval University of Engineering,Wuhan 430033,China;2 Equipment Department,Army 92196,Qingdao 266011,China;3 Transportation College,Wuhan University of Technology,Wuhan 430063,China)

The motion stability theory of a nonlinear system and the principles of continuation based on homotopy techniques were applied to stable state response and motion stability analysis of submarines during emergency ascent.The stability analysis results were verified by simulations using direct numerical integrations of the fully nonlinear,coupled six degrees-of-freedom equations of motion for the submarine.It indicates that homotopy and continuation theory techniques could be an effective analysis method for motion stability of submarines under water.

submarine;motion stability;nonlinear dynamics;homotopy and continuation theory techniques

U661.3

A

1007-7294（2011）08-0844-09

2010-10-11 修改日期：2011-03-28

國防預研項目資助（1010501020301）

戴余良（1966-），男，博士，副教授，主要從事船舶操縱運動建模、智能控制與實時仿真研究，E-mail:yuliang_dai@163.com。