系統(tǒng)參數(shù)對超空泡航行體非線性動力學(xué)特性的影響

呂一品， 熊天紅， 易文俊， 管 軍

(南京理工大學(xué) 瞬態(tài)物理國家重點(diǎn)實驗室，南京 210094)

當(dāng)航行體在水中運(yùn)動時，若水中某點(diǎn)的壓力降低到臨界壓力以下，水將發(fā)生汽化，先是微觀的，進(jìn)而發(fā)展成宏觀的小氣泡，然后在水中或水與航行體的交界面上，匯合形成較大的蒸汽與氣體空腔，稱為空泡，若整個航行體表面以及尾端附近的水中都出現(xiàn)空泡，猶如一個大汽袋，將整個航行體裝于其中，這種空泡稱為超空泡[1-2]。形成超空泡之后，水的阻力就變?yōu)橛伤魵鈽?gòu)成的空泡的阻力，大大地提高了運(yùn)動速度，實現(xiàn)了航行體在水中飛行。但由于復(fù)雜的水動力特性具有時滯性、不確定性和時變性等問題，超空泡航行體的運(yùn)動極容易受到周圍干擾因素的影響，系統(tǒng)參數(shù)如空化數(shù)、空化器直徑等能夠改變超空泡的尺寸，造成航行體尾部與空泡之間的碰撞，產(chǎn)生非線性滑行力，進(jìn)而造成航行體的振動[3-6]，使航行體的運(yùn)動具有強(qiáng)烈的非線性動力學(xué)特性。

在工程應(yīng)用中，非線性動力學(xué)特性是超空泡航行體穩(wěn)定控制的重要理論依據(jù)，Lin等[7]提出了超空泡航行體系統(tǒng)隨著空化數(shù)的變化具有復(fù)雜的非線性動力學(xué)行為；白濤等[8]分叉分析了當(dāng)尾翼偏轉(zhuǎn)角變化時航行體的運(yùn)動穩(wěn)定性；熊天紅等[9]探討了超空泡系統(tǒng)隨空化器偏轉(zhuǎn)角的變化而產(chǎn)生的非線性物理現(xiàn)象?，F(xiàn)有的研究多停留在單一參數(shù)的變化對超空泡航行體非線性動力學(xué)特性的影響。本文將從非線性的角度對超空泡航行體系統(tǒng)進(jìn)行分析，運(yùn)用Lyapunov穩(wěn)定性理論確定超空泡系統(tǒng)動力學(xué)行為與參數(shù)的對應(yīng)關(guān)系，探討系統(tǒng)隨多個參數(shù)變化產(chǎn)生的非線性物理現(xiàn)象，揭示其穩(wěn)定運(yùn)動的參數(shù)條件以及運(yùn)動特性，為超空泡航行體控制器的設(shè)計提供理論研究基礎(chǔ)。

1 水下超空泡航行體的動力學(xué)描述

1.1 超空泡航行體的外形

超空泡航行體的外形結(jié)構(gòu)如圖1所示，主要構(gòu)成有：頭部直徑為D的圓盤空化器，前端圓錐段和后端半徑為R的圓柱段，以及采用十字形布局的鰭舵，圓錐段長度是圓柱段的一半, 4個尾翼對稱地布置在航行體的尾端。

圖1 超空泡航行體示意圖Fig.1 Schematic diagram of a supercavitating vehicle

流場的外部條件和水下航行體幾何形狀等都會使空泡出現(xiàn)不同的狀態(tài)，通常用空化數(shù)σ來表征相似的空泡狀態(tài)，定義為σ=2(p∞-pc)/ΡV2，其中，P∞為無窮遠(yuǎn)處的壓力，Pc為空泡內(nèi)部的壓力，ρ為水的密度，V為航行體的合速度。

當(dāng)航行體高速運(yùn)動時，與水接觸的只有空化器和尾翼，空泡尺寸會隨著空化數(shù)的減小而增大，可以通過改變空化數(shù)來改變流場的空化狀態(tài)；合理調(diào)節(jié)尾翼偏轉(zhuǎn)角度不僅能夠控制航行體在超空泡內(nèi)的姿態(tài)，還能有效地抑制非線性滑行力的出現(xiàn)，提高航行體的穩(wěn)定性；另外，超空泡形態(tài)的優(yōu)良率與空化器直徑密切相關(guān)[10]。故通常將空化數(shù)、空化器直徑和尾翼偏轉(zhuǎn)角作為可變參數(shù)來分析其對超空泡航行體動力學(xué)特性的影響。

1.2 超空泡航行體動力學(xué)模型

圖2所示，超空泡航行體動力學(xué)建模的體坐標(biāo)系原點(diǎn)位于航行體頭部空化器頂端面的圓心，X軸的方向與航行體的中心軸線重合指向前，Y軸的方向垂直于中心軸線指向下，把地面系當(dāng)作慣性系, 本文研究航行體在縱平面內(nèi)的運(yùn)動。建模采用航行體所處的深度z、垂直速度w、俯仰角θ、俯仰角速度q作為狀態(tài)變量來描述其動力學(xué)，垂直速度w的方向與Y軸方向一致，V為縱平面內(nèi)航行體的合速度，X軸方向的速度為Vx。另外，系統(tǒng)具有兩個控制輸入，分別為尾翼偏轉(zhuǎn)角δe和空化器偏轉(zhuǎn)角δc。

圖2 超空泡航行體縱向建模示意圖Fig.2 The longitudinal model of a supercavitating vehicle

根據(jù)航行體各部分所受的流體動力[11]，結(jié)合動量和動量距定理可以推出上述變量有如下關(guān)系[12]

(1)

(2)

(3)

超空泡航行體模型參數(shù)值如表1所示。

表1 超空泡航行體的參數(shù)值

當(dāng)空化器直徑D較小時，不能形成空泡或者形成部分空泡，而過大的空化器直徑在降低航行體摩擦阻力的同時大幅度增加了壓差阻力，導(dǎo)致減阻效果不明顯，即D過大過小都不利于超空泡航行體的航行。因此本文選取空化器直徑D的研究范圍為[0.010 0,0.030 0] m。

Dzielski 等提出了經(jīng)典控制律δe=0、δc=15z-30θ-0.3q，但由于尾翼偏轉(zhuǎn)角為零，往往會使航行體缺乏尾翼提供的支持力，當(dāng)空化器提供的升力不能平衡航行體的重力時，將導(dǎo)致航行體在重力的作用下浸入空泡而失穩(wěn)。因此，本文令控制律為δe1=kq、δc1=15z-30θ-0.3q，k為航行體俯仰角速度q的反饋控制增益，整理式(1)～式(3)式可得到以空化器直徑D、空化數(shù)σ、尾翼偏轉(zhuǎn)角反饋控制增益k為可變參數(shù)的動力學(xué)方程式

(4)

式中：C=158.914 7(1+σ)D2；a22=-4.207 9CV+0.425 0；a24=-6.154 0CV+14.710 1V；a42=1.515 4CV；a44=-1.304 6CV；b21=0.616 7CV2；b22= -2.805 4CV2；b41=-0.724 9CV2；b42=2.240 6CV2；d2=-1.226 6；d4=1.449 2。Fplaning為航行體尾部與空泡碰撞產(chǎn)生的非線性滑行力，該力的出現(xiàn)往往會導(dǎo)致航行體產(chǎn)生振動與沖擊?；辛Ρ磉_(dá)式

(5)

式中：R′=Rc-R，Rc為空泡半徑，R為航行體半徑；h為航行體尾部探出空泡的深度，稱作浸沒深度；α為航行體中心線與空泡中心線之間的幾何角為浸沒角，其表達(dá)式為

(6)

(7)

2 超空泡航行體的動力學(xué)特性

2.1 超空泡航行體的動力學(xué)地圖

為完整地展現(xiàn)超空泡航行體的動力學(xué)行為與參數(shù)的依賴關(guān)系，本文基于動力學(xué)模型(4)，隨機(jī)地選取初始條件，依照Lyapunov穩(wěn)定性理論[13]將模型的穩(wěn)定解、周期解、混沌解用不同的顏色在圖3中表示出來，繪制出超空泡航行體的動力學(xué)地圖。圖中淺灰色區(qū)域LG表示航行體穩(wěn)定運(yùn)動；深灰色區(qū)域DG表示航行體周期振蕩；白色區(qū)域W表示航行體會出現(xiàn)劇烈的振動與沖擊；黑色區(qū)域D表示系統(tǒng)發(fā)散。圖3 (a)為空化器直徑D=0.019 1 m時，分岔參數(shù)σ和k所定義的動力學(xué)地圖；圖3(b)為尾翼偏轉(zhuǎn)角反饋控制增益k=0時，分岔參數(shù)σ和D所定義的動力學(xué)地圖。

根據(jù)圖3(a)能夠確定航行體穩(wěn)定運(yùn)動的σ和k的取值范圍。當(dāng)σ確定時，只要在對應(yīng)LG區(qū)域內(nèi)調(diào)節(jié)k的取值，就能夠有效實現(xiàn)超空泡航行體的穩(wěn)定航行，對航行體的穩(wěn)定性控制具有指導(dǎo)意義。前文提到航行體空化器直徑D過小或者過大都不利于超空泡減阻效果，結(jié)合圖3(b)可以選擇合適的D，既能形成將航行體包裹住的完整空泡又不至于形成的空泡厚度過大而增加航行阻力，對超空泡航行體穩(wěn)定高速地運(yùn)動具有實際指導(dǎo)意義。

圖3 動力學(xué)地圖 Fig.3 Regions of different dynamical behaviors

分別在圖3中的LG區(qū)域、DG區(qū)域和W區(qū)域中各隨機(jī)地選取一點(diǎn)，圖4為三種狀態(tài)的相軌跡圖。LG區(qū)域中選取的系統(tǒng)參數(shù)為σ1=0.031 20，D1=0.019 1 m，k1=6，如圖4(a)所示，垂直速度w和俯仰角θ在反饋控制律的作用下被吸引到穩(wěn)定平衡點(diǎn)上，航行體在該參數(shù)的作用下穩(wěn)定航行。

DG區(qū)域中選取的參數(shù)為：σ2=0.028 10，D2=0.018 0 m，k2=0，相軌跡如圖4(b)所示，系統(tǒng)映射形成了閉合極限環(huán)，極限環(huán)與臨界切換分界線w=wth相交，其中wth=(Rc-R)V/L，Rc為空泡半徑，垂直速度w時而大于wth,時而小于wth，航行體尾部不時地擺動與空泡碰觸，引起滑行力周期性的變化，尾部時而穿過空泡插入水中，產(chǎn)生滑行力，時而被空泡包裹，無滑行力，航行體周期振蕩。

W區(qū)域中選取的參數(shù)為:σ3=0.032 60，D3=0.019 1 m，k3=-0.030，相軌如圖4(c)所示，混沌吸引子的出現(xiàn)表明航行體的運(yùn)動具有復(fù)雜的非線性動力學(xué)行為。

圖3(a)相對完整地反應(yīng)了航行體隨參數(shù)σ、k變化時所處的不同狀態(tài)，深灰色周期振蕩是主要的運(yùn)動狀態(tài)，周期振蕩對應(yīng)的極限環(huán)也是主要的動力學(xué)行為。當(dāng)k∈[1.041,5.122]時，航行體穩(wěn)定運(yùn)動的空化數(shù)σ取值范圍較大；當(dāng)k∈[5.121,27.553]時，淺灰色穩(wěn)定點(diǎn)星羅棋布地灑在深灰色周期區(qū)域中，σ的微小變化總是會引起運(yùn)動狀態(tài)在穩(wěn)定與周期之間轉(zhuǎn)變；研究該區(qū)域還發(fā)現(xiàn)，在某些參數(shù)處，通過仿真計算，隨機(jī)設(shè)置初始值，航行體表現(xiàn)出不同的動力學(xué)行為。

隨機(jī)選取(σ,k)=(0.030 68, 25.521)，當(dāng)初始值為α1(0.135 2, 0.515 2, 0.261 4, -0.941 5)時，相軌在w-θ平面上的投影如圖5(a)所示，形成了一個穩(wěn)定平衡點(diǎn)。圖5(b)為隨時間演化的 Lyapunov 指數(shù)譜，系統(tǒng)在反饋控制律的作用下迅速進(jìn)入收斂狀態(tài)，L4

圖4 系統(tǒng)相軌在w-θ平面上的映射 Fig.4 Phase track diagram for diffident parameters

圖5 σ=0.030 68, k=25.521時的運(yùn)動現(xiàn)象Fig.5 Motion phenomenon for σ=0.030 68, k=25.521

圖3(b)相對完整地反應(yīng)了航行體隨參數(shù)σ、D變化時所處的不同狀態(tài)，三種顏色代表的三種運(yùn)動狀態(tài)呈階梯式分布，在D取值較小的黑色發(fā)散區(qū)域中，存在個別的深灰色周期點(diǎn)和白色混沌點(diǎn)，這是由于D過小，不能產(chǎn)生空泡或者產(chǎn)生了部分空泡，航行體幾乎在整個空化數(shù)有效范圍內(nèi)嚴(yán)重失穩(wěn)；當(dāng)D>0.026 2 m時，運(yùn)動過程中會形成將整個航行體包裹在內(nèi)的大空泡，航行體全局穩(wěn)定運(yùn)動。當(dāng)D∈[0.010 8, 0.025 7]時，深灰色周期振蕩是主要的運(yùn)動區(qū)域，但在σ=0.031 48和σ=0.032 55附近，深灰色區(qū)域中穿插著白色混沌點(diǎn)，D的微小變化會引起運(yùn)動狀態(tài)在周期與混沌之間的轉(zhuǎn)變；通過仿真計算發(fā)現(xiàn)，當(dāng)該區(qū)域中某些系統(tǒng)參數(shù)一定時，隨機(jī)設(shè)置初始值的不同會導(dǎo)致航行體具有不同的動力學(xué)行為。

隨機(jī)選擇(σ,D)=(0.032 55, 0.019 1)，當(dāng)初始值為β1(1.436 7, -1.960 9, -0.197 7, -1.207 8)時，相軌跡在w-θ平面上的投影如圖6(a)所示，出現(xiàn)了極限環(huán)，形成了一個周期吸引子。圖6(b)為隨時間演化的 Lyapunov 指數(shù)譜，L1≈0，L40，L49.20 s時L1≈0，根據(jù)Lyapunov穩(wěn)定性理論，上述混沌現(xiàn)象是瞬態(tài)的，系統(tǒng)最終將收斂到周期狀態(tài)。由上述分析可知，在一些參數(shù)組合處，當(dāng)初始值不同時，航行體會出現(xiàn)穩(wěn)定周期航行和瞬態(tài)混沌振蕩兩種不同的運(yùn)動軌跡。

圖6 σ=0.032 55, D=0.019 1時的運(yùn)動現(xiàn)象 Fig.6 Motion phenomenon for σ=0.032 55, D=0.019 1

2.2 參數(shù)變化對系統(tǒng)動力學(xué)的影響分析

2.2.1 空化數(shù)對超空泡航行體動力學(xué)特性的影響

令D=0.019 1 m，k=0，σ∈[0.019 80,0.036 80]，即在經(jīng)典控制律的作用下，系統(tǒng)隨空化數(shù)σ變化的分岔情況如圖7所示。

當(dāng)系統(tǒng)處于小空化數(shù)σ∈[0.019 80, 0.024 20]時，圖7 (a)中的信號點(diǎn)為穩(wěn)定平衡點(diǎn)，航行體穩(wěn)定運(yùn)動。當(dāng)σ=0.024 25時，發(fā)生Hopf分岔，相軌由穩(wěn)定的有界點(diǎn)突變成極限環(huán)，極限環(huán)的出現(xiàn)使尾部滑行力開始周期性的變化，系統(tǒng)進(jìn)入周期軌道。圖7(b)是圖7(a)在σ∈[0.032 40,0.032 70]的放大，周期軌道經(jīng)過一系列倍周期分岔，于σ=0.032 46處進(jìn)入混沌狀態(tài)，混沌寬帶中穿插著多個周期窗，當(dāng)σ=0.032 70時，混沌寬帶突變回2周期軌道，之后，如圖7(a)所示，2周期軌道在σ=0.033 96處突變到1周期軌道，航行體開始劇烈的周期振蕩甚至傾覆。

圖7 系統(tǒng)隨σ變化的分岔圖 Fig.7 Bifurcation diagram of system for σ

總體來說，此時航行體穩(wěn)定運(yùn)動的空化數(shù)范圍是σ∈[0.019 80,0.024 10]，而在空化數(shù)的有效范圍σ∈ [0.024 20,0.036 80]下無法實現(xiàn)穩(wěn)定航行，很大程度上限制了航行體的可控性以及穩(wěn)定性。根據(jù)圖3(a)，在較大的σ范圍內(nèi)，調(diào)節(jié)k的取值，能夠有效實現(xiàn)超空泡航行體的穩(wěn)定航行，新的控制律有助于提高航行體的可控性和穩(wěn)定性。

2.2.2 尾翼偏轉(zhuǎn)角對超空泡航行體動力學(xué)特性的影響

令σ=0.032 55、D=0.019 1 m，尾翼偏轉(zhuǎn)角δe1=kq，得到系統(tǒng)隨k變化的分岔圖8，當(dāng)k<-0.6時，系統(tǒng)發(fā)散，此處不做研究。圖8(b)是圖8(a)在k∈[-0.500,0.500]范圍的放大，當(dāng)-0.50034.782時，尾翼偏轉(zhuǎn)角過大，產(chǎn)生的流體動力會使航行體傾覆，系統(tǒng)發(fā)散。

圖8 系統(tǒng)隨k變化的分岔圖 Fig.8 Bifurcation diagram of system for k

2.2.3 空化器直徑對超空泡航行體動力學(xué)特性的影響

令k=0、σ=0.030 00，圖9為系統(tǒng)運(yùn)動隨空化器直徑D變化的分岔圖，如圖9(a)所示，當(dāng)D∈[0.015 8, 0.175 2] m時，1周期軌道表明航行體周期振蕩，在D=0.017 5 m處，跳變到2周期軌道，圖9(b)是圖9(a)在D∈ [0.017 7,0.017 9] m范圍的放大， 由2周期軌道進(jìn)入混沌寬帶，其間經(jīng)歷較窄的周期窗，然后經(jīng)過一系列的反倍周期分岔，在D=0.017 9 m附近進(jìn)入3周期軌道，如圖9(a)所示，當(dāng)D=0.018 3 m時，3周期軌道突變到1周期軌道，隨著D的增大，在D=0.022 5 m處發(fā)生Hopf分岔，不穩(wěn)定平衡點(diǎn)演變成穩(wěn)定平衡點(diǎn)，系統(tǒng)由周期狀態(tài)突變到穩(wěn)定狀態(tài)，之后，圖中穩(wěn)定的平衡點(diǎn)說明系統(tǒng)一直處于穩(wěn)定狀態(tài)。在實際工程應(yīng)用中，可以根據(jù)航行體運(yùn)動特性隨空化器直徑變化的分岔圖，選擇合適的空化器直徑，避免在非線性現(xiàn)象豐富的空化器直徑處取值，以保證航行體的穩(wěn)定航行。

圖9 系統(tǒng)隨D變化的分岔圖 Fig.9 Bifurcation diagram of system for D

依照Lyapunov穩(wěn)定性理論可以確定系統(tǒng)動力學(xué)行為與任意兩個可變參數(shù)的對應(yīng)關(guān)系，依照分岔圖可以確定系統(tǒng)隨任一參數(shù)變化而呈現(xiàn)的動力學(xué)現(xiàn)象，前一種側(cè)重于為航行體的穩(wěn)定運(yùn)動提供參數(shù)依據(jù)，后一種側(cè)重于研究系統(tǒng)隨參數(shù)變化的運(yùn)動規(guī)律以及產(chǎn)生的具體非線性物理現(xiàn)象。結(jié)合兩種理論可以相對完整地分析系統(tǒng)的動力學(xué)特性。

3 不同參數(shù)下航行體運(yùn)動狀態(tài)分析

表2 不同運(yùn)動區(qū)域的參數(shù)取值

按照表2中的A、B、C三點(diǎn)對應(yīng)的參數(shù)取值，隨機(jī)地選取初始值，系統(tǒng)時域仿真結(jié)果和頻譜如圖10所示，每一個變量的頻譜對應(yīng)在該圖的下方，其中虛線、實線和點(diǎn)線分別表示A、B、C點(diǎn)處的運(yùn)動情況。在A點(diǎn)處，即σA=0.031 20，DA=0.019 1 m，kA=6時，狀態(tài)變量z、w、θ、q均被吸引到平衡點(diǎn)PA上；浸沒深度h=0.014 m，滑行力Fplaning=12.94 N，空化器偏轉(zhuǎn)角δc1=0.054 4 rad，尾翼偏轉(zhuǎn)角δe1=0；頻譜能量集中在0頻率處。在B點(diǎn)處，即σB=0.028 10，PB=0.018 0 m，kB=0時，四個狀態(tài)變量均圍繞平衡點(diǎn)PB周期振蕩；h和Fplaning分別在[0,0.077] m和[0,63.40] N之間周期振蕩，δc1在[0.009 9,0.139 6] rad之間周期振蕩，航行體尾部不斷與空泡壁面發(fā)生碰撞；頻譜均為多個等間隔的離散尖峰。而在參數(shù)取值為σC=0.032 60，DC=0.019 1 m，kC=-0.030的C點(diǎn)處，z、w、θ、q、h、Fplaning、δe1、δc1隨著時間變化發(fā)生了劇烈的非周期振蕩，h最大值為0.380 m，該長度超過了航行體直徑，對應(yīng)的Fplaning在300 N左右，在實際中不可能出現(xiàn)該現(xiàn)象，航行體已傾覆，對應(yīng)的頻域均是寬帶連續(xù)譜。

(a) z、w、θ、q的時域、頻域圖

(b) h、Fplaning、δc、δe的時域、頻域圖圖10 A、B、C點(diǎn)處的時域、頻域圖Fig.10 Time response and frequency domain response for A , B , C point

總的來說，系統(tǒng)參數(shù)在A點(diǎn)取值時，由航行體空化器的偏轉(zhuǎn)產(chǎn)生的流體動力以及尾部伸出空泡與流體接觸所產(chǎn)生的滑行力來共同維持航行體的平衡，使其處于小攻角斜向運(yùn)動的穩(wěn)定狀態(tài)；參數(shù)在B點(diǎn)處取值時，航行體處于周期振蕩的運(yùn)動狀態(tài)；在C點(diǎn)處，航行體將失穩(wěn)進(jìn)而傾覆。由此可知，系統(tǒng)參數(shù)對超空泡航行體的運(yùn)動狀態(tài)有著至關(guān)重要的作用。

系統(tǒng)參數(shù)在D點(diǎn)處取值，即σD=0.030 68，DD=0.019 1，kD=25.521時，系統(tǒng)時域仿真結(jié)果和頻譜如圖11所示，每一個變量的頻譜對應(yīng)在該圖的下方，虛線表示初始值α1(z0,w0,θ0,q0)=(0.135 2, 0.515 2, 0.261 4, -0.941 5)時，系統(tǒng)四個狀態(tài)變量、浸沒深度，滑行力以及空化器和尾翼偏轉(zhuǎn)角隨時間變化的運(yùn)動情況以及對應(yīng)的頻譜，實線則表示初始值α2(z0,w0,θ0,q0)=(-0.714 5, 1.351 4, -0.224 8, -0.589 0)時，系統(tǒng)各個變量的時域和頻域響應(yīng)情況。

初始值為α1時，在控制律的作用下，狀態(tài)變量z,w,θ,q迅速穩(wěn)定在平衡點(diǎn)PD上，h=0.013 m，F(xiàn)planing=13.20 N，δc1=0.042 6 rad，δe1=0.092 9 rad，航行體的尾部穿過空泡浸入水中，產(chǎn)生的滑行力與空化器和尾翼提供的升力共同平衡航行體重力，使航行體穩(wěn)定運(yùn)動；對應(yīng)的頻譜能量集中在0頻處，說明系統(tǒng)很快進(jìn)入穩(wěn)定狀態(tài)，微弱幅度的周期振蕩不影響實際的穩(wěn)定性。

初始值為α2時，隨著時間的增大，系統(tǒng)的四個狀態(tài)變量均減幅振蕩，并逐漸向平衡點(diǎn)PD靠攏，直至t=11.02 s附近，z,w,θ,q穩(wěn)定在PD上；h和Fplaning一開始均為零，在t=10.88 s時突然開始周期振蕩，并于t=11.89 s穩(wěn)定到0.013 m和13.20 N上；δc1減幅振蕩并于t=11.22 s穩(wěn)定在0.042 6 rad上，δe1在t=10.85 s時發(fā)生狀態(tài)轉(zhuǎn)移于t=11.95 s穩(wěn)定在0.092 9 rad上?？梢姾叫畜w發(fā)射之后，在初始階段處于小攻角斜向的振蕩狀態(tài)，尾部一直處于超空泡內(nèi)，隨著時間的增加，各個變量均穩(wěn)定在固定值上，尾部刺破空泡與水接觸，航行體的運(yùn)動由周期狀態(tài)轉(zhuǎn)移到穩(wěn)定狀態(tài)。對應(yīng)的頻譜除在0頻外有一對顯著的峰值，系統(tǒng)存在較長時間的顯著周期振蕩。

(a) z、w、θ、q的時域、頻域圖

(b) h、Fplaning、δc、δe的時域、頻域圖圖11 D點(diǎn)處的時域、頻域圖Fig.11 Time response and frequency domain response for D point

在E點(diǎn)處，即σE=0.032 55,DE=0.019 1，kE=0，系統(tǒng)時域仿真結(jié)果和頻域響應(yīng)情況如圖12所示，為了清晰顯示航行體運(yùn)動狀態(tài)的轉(zhuǎn)移，取t∈[8.50,9.50] s為圖示的時間范圍。點(diǎn)線表示初始值β1(z0,w0,θ0,q0)=(1.436 7, -1.960 9, -0.197 7, -1.207 8)時，系統(tǒng)四個狀態(tài)變量、浸沒深度，滑行力以及空化器和尾翼偏轉(zhuǎn)角隨時間變化的運(yùn)動情況以及對應(yīng)的頻譜；實線為初始值β2(z0,w0,θ0,q0)=(0.318 8, -1.307 7, -0.433 6, 0.342 6)時系統(tǒng)各個變量的時域、頻域響應(yīng)情況。

當(dāng)初始值為β1時，z、w、θ、q分別以平衡點(diǎn)PE為中心周期振蕩，h的周期性變化為0.126 m、0 m、0.051 m，產(chǎn)生Fplaning對應(yīng)為138.30 N、0 N、69.22 N，航行體尾部在空泡內(nèi)外周期振蕩；δc1也在平衡點(diǎn)附近周期振蕩，k=0，δe1=0；對應(yīng)的頻譜均為多個等間隔的分離尖峰。

當(dāng)初始值為β2時，在t=9.20 s之前，z、w、θ、q以及h、Fplaning、δc1隨時間的變化均發(fā)生了劇烈的非周期振蕩，t=9.20 s附近，各個變量從非周期振蕩突變?yōu)橹芷谡袷帲笠恢狈€(wěn)定在周期狀態(tài)，與初始值為β1時的運(yùn)動狀態(tài)重合，航行體的運(yùn)動發(fā)生了狀態(tài)轉(zhuǎn)移，在有限時間尺度內(nèi)出現(xiàn)的混沌稱之為瞬態(tài)混沌[14]，該現(xiàn)象被稱為瞬態(tài)混沌穩(wěn)態(tài)周期；對應(yīng)的頻域是寬帶連續(xù)譜。

通過系統(tǒng)參數(shù)在D、E兩點(diǎn)取值時超空泡航行體運(yùn)動狀態(tài)的分析，不難發(fā)現(xiàn)，航行體的運(yùn)動狀態(tài)除了隨系統(tǒng)參數(shù)變化外，還依賴于系統(tǒng)的初始條件，在不同的發(fā)射初始條件下，航行體運(yùn)動狀態(tài)存在瞬態(tài)振蕩和運(yùn)動狀態(tài)轉(zhuǎn)移的現(xiàn)象。在實際工程應(yīng)用中，由于水下環(huán)境復(fù)雜，外界干擾因素多，系統(tǒng)參數(shù)一定時，發(fā)射初始條件一旦受到外界的微小擾動，航行體就會出現(xiàn)不同的運(yùn)動狀態(tài)。考慮到超空泡航行體在水下航行時間較短，針對以上航行體瞬態(tài)振蕩和運(yùn)動狀態(tài)轉(zhuǎn)移的穩(wěn)定控制，更具有實際工程意義。

(a) z、w、θ、q的時域、頻域圖

(b) h、Fplaning、δc、δe的時域、頻域圖圖12 E點(diǎn)處時域、頻域圖Fig.12 Time response and frequency domain response for E point

4 結(jié) 論

基于動力學(xué)地圖、分岔圖以及數(shù)值仿真等方法，本文研究了多個系統(tǒng)參數(shù)對超空泡航行體非線性動力學(xué)特性的影響，得到以下結(jié)論：

(1)隨著尾翼偏轉(zhuǎn)角、空化器直徑等參數(shù)的變化，超空泡航行體系統(tǒng)均會出現(xiàn)分岔、混沌等非線性物理現(xiàn)象。

(2) 超空泡航行體的運(yùn)動狀態(tài)不但受系統(tǒng)參數(shù)影響，還極端依賴系統(tǒng)的初始條件，在不同的初始條件下發(fā)射，航行體會出現(xiàn)不同的運(yùn)動軌跡，存在瞬態(tài)振蕩和運(yùn)動狀態(tài)轉(zhuǎn)移等非線性運(yùn)動特性。

(3) 利用動力學(xué)地圖一方面可以確定航行體穩(wěn)定運(yùn)動的參數(shù)范圍，另一方面也為避開瞬態(tài)振蕩的參數(shù)范圍提供了理論依據(jù)，從而設(shè)置合理的系統(tǒng)參數(shù)，抑制航行體的非線性振動與沖擊，實現(xiàn)其穩(wěn)定運(yùn)動。

[ 1 ] 曹偉，魏英杰，王聰，等． 超空泡技術(shù)現(xiàn)狀問題與應(yīng)用[J]．力學(xué)進(jìn)展，2006，36(4)： 571-579．

CAO Wei，WEI Yingjie，WANG Cong，et al．Current status problems and applications of supercavitating technology [J]．Advances in Mechanics，2006，36(4)： 571-579．

[ 2 ] SINGHAL A K，LI H Y，AHAVALE M M，et al．Mathematical basis and validation of the full cavitating high-speed torpedo[J]．Journal of Fluids Engineering，2003，125： 459-468．

[ 3 ] VANEK B，BOKOR J，BALAS G J，et al．Longitudinal motion control of a high-speed supercavitation vehicle[J]．Journal of Vibration and Control，2006，13(2)： 159-184．

[ 4 ] WANG G，OSTOJA-STARZEWSKI M．Large eddy simulation of a sheet/cloud cavitation on a NACA0015 hydrofoil [J]．Applied Mathematical Modeling，2007，31(3)： 417-447．

[ 5 ] LI Qitao，HE Yousheng，XUE Leiping．A numerical simulation of pitching motion of the ventilated supercaviting vehicle around its nose[J]． Chinese Journal of Hydrodrnamics，2011， 26 (6)： 589-685．

[ 6 ] 魏英杰，王京華，張嘉鐘，等．水下超空泡航行體非線性動力學(xué)與控制[J]．振動與沖擊，2009，28(6)： 179-204．

WEI Yingjie，WANG Jinghua，ZHANG Jiazhong，et al．Nonlinear dynamics and control of underwater supercavitating vehicle[J]．Journal of Vibration and Shock， 2009，28(6)： 179-204．

[ 7 ] LIN G J，BALACHANDRAN B，ABED E H. Bifurcation behavior of a supercavitating vehicle[C]∥ 2006 ASME international Mechanical Engineering Congress and Exposition. Chicago: ASME, 2006： 5-10．

[ 8 ] 白濤，孫堯，莫宏偉．分叉分析在水下高速運(yùn)動體穩(wěn)定控制中的應(yīng)用[J]．哈爾濱工程大學(xué)學(xué)報，2008，29(10)： 1067-1075．

BAI Tao，SUN Yao，MO Hongwei．Application of bifurcation analysis to the stability control of underwater high-speed vehicles[J]．Journal of Harbin Engineering University，2008， 29(10)： 1067-1075．

[ 9 ] 熊天紅，包伯成．超空泡航行體閉環(huán)控制動力學(xué)特性研究[J]．振動與沖擊, 2015，34(17)： 168-173．

XIONG Tianhong，BAO Bocheng．Closed-loop control dynamics for supercavitating vehicles[J]．Journal of Vibration and Shock， 2015，34(17)： 168-173．

[10] 賈力平，王聰，于開平，等．空化器參數(shù)對通氣超空泡形態(tài)影響的實驗研究[J]．工程力學(xué)，2007，24(3)： 159-164．

JIA Liping，WANG Cong，YU Kaiping，et al．Experimental investigation of cavitator parameters effecting on ventilated supercavity shape[J]．Engineering Mechanics，2007，24(3)： 159-164．

[11] DZIELSKI J，KURDILA A．A benchmark control problem for supercavitating vehicles and an initial inves-tigation of solutions[J]．Journal of Vibration and Control，2003，9(7)：791-804．

[12] LIN G J，BALACHANDRAN B，ABED E H．Dynamics and control of supercavitating vehicles[J]. Journal of Dynamic Control Systems， Measurement， and Control，2008，130(2)： 281-287．

[13] 包伯成．混沌電路導(dǎo)論[M]．北京：科學(xué)出版社，2013： 5-38.

[14] 包伯成，胡文，許建平，等．憶阻混沌電路的分析與實現(xiàn)[J]．物理學(xué)報，2011，60(12)： 58-65．

BAO Bocheng，HU Wen，XU Jianping，et al．Analysis and implementation of memristor chaotic circuit[J]．Acta Physica Sinica，2011，60(12)： 58-65．