PMUCR 方法對浮筏系統(tǒng)的全局分析與應(yīng)用

張樹楨，陳 前

（1．南京航空航天大學(xué) 機(jī)械結(jié)構(gòu)力學(xué)及控制國家重點(diǎn)實(shí)驗(yàn)室，南京 210016；2．中國直升機(jī)設(shè)計(jì)研究所，景德鎮(zhèn) 333001）

在船舶隔振分析中，由于受隔振器非線性剛度的影響，浮筏隔振系統(tǒng)的動(dòng)力學(xué)特性較為復(fù)雜。同時(shí)由于設(shè)備多樣，擾動(dòng)源的變化，也會(huì)影響浮筏的動(dòng)力學(xué)特性，所以浮筏隔振系統(tǒng)的全局分析是一個(gè)亟待解決的問題。

目前，全局分析的有效方法為胞映射方法。Hsu［1－2］提出了簡單胞映射法（SCM）以及廣義胞映射法（GCM），對胞映射的分析過程進(jìn)行了詳盡的描述，并進(jìn)行了相應(yīng)的改良。Tongue等［3］提出了插值胞映射（ICM）法，該方法除了記錄積分過程中胞的軌跡外，還記錄了積分時(shí)點(diǎn)的具體坐標(biāo)，因此可通過插值計(jì)算下一步的映射。Levitas等［4］提出龐加萊型胞映射法，該方法基于龐加萊截面將系統(tǒng)降維，有效地減少了胞映射的計(jì)算時(shí)長和難度。該方法將系統(tǒng)降低一維，對四維以下的系統(tǒng)效果較好，但隨著維數(shù)的升高，降低一維也難以解決問題。

文成秀等［5］提出了胞映射－點(diǎn)映射法，該方法將胞映射得到的周期解用點(diǎn)映射加以驗(yàn)證，從而保證了計(jì)算的精度和有效性。江俊等［6－7］對胞映射方法在實(shí)際動(dòng)力學(xué)系統(tǒng)的應(yīng)用進(jìn)行了大量的研究，提出了胞化積分軌跡法，胞參考點(diǎn)映射法。將點(diǎn)映射和胞映射進(jìn)行了有效結(jié)合，兼顧了兩者的優(yōu)點(diǎn)，有效地提高了分析效率。

隨著系統(tǒng)維數(shù)的增高，所需處理的數(shù)據(jù)急劇增加，因此需要對全局分析方法加以改進(jìn)。胞參考點(diǎn)映射法（PMUCR）結(jié)合了胞映射的高效性和點(diǎn)映射的精確性，且易于推廣到高維空間。與其它胞映射法相比，不需要在相空間的每個(gè)胞中選取初始點(diǎn)，并可采用動(dòng)態(tài)存儲的方法［8］，減少對內(nèi)存的消耗。

本文試圖采用PMUCR方法，進(jìn)行對稱激勵(lì)條件下高維浮筏隔振系統(tǒng)的全局分析，并判斷剛性和柔性筏架條件下，全局分析結(jié)果的差異。討論P(yáng)MUCR方法在高維系統(tǒng)中的運(yùn)用準(zhǔn)則。最后通過打靶法對吸引子的穩(wěn)定性進(jìn)行判斷。

1 計(jì)算方法

1.1 PMUCR方法的改良

對于低維系統(tǒng)，PMUCR方法可以有效地進(jìn)行分析。在高維條件下，當(dāng)系統(tǒng)增加一維，則導(dǎo)致增加一個(gè)數(shù)量級的計(jì)算量，同時(shí)PMUCR中的系統(tǒng)特性保持標(biāo)準(zhǔn)也需要進(jìn)行調(diào)整，以適應(yīng)多變量判斷的需求。對PMUCR方法在高維系統(tǒng)中的改良措施為：

1．1．1 PMUCR方法同簡單胞映射方法的一個(gè)重要區(qū)別是判斷某點(diǎn)落入某胞后的狀態(tài)。簡單胞映射方法認(rèn)定落入同一胞內(nèi)的點(diǎn)具有相同的屬性。當(dāng)胞的尺寸較小且各周期邊界規(guī)則時(shí)，該判斷基本合理，但當(dāng)其它情況下則可能出現(xiàn)偽周期胞。PMUCR方法則通過點(diǎn)－點(diǎn)的間距來判斷是否具有相同的屬性，即間距是否滿足系統(tǒng)特性維持標(biāo)準(zhǔn)。從運(yùn)動(dòng)軌跡分析，如一條軌跡接近另一條已知?jiǎng)恿W(xué)特性的軌跡，且兩者之間的距離小于某一定值時(shí)，則可認(rèn)定兩條軌跡具有相同的特性，該距離即為系統(tǒng)特性維持標(biāo)準(zhǔn)中的推薦距離。該推薦距離需要同時(shí)滿足多變量的需求，而不是直接求解點(diǎn)－點(diǎn)之間的間距，如果系統(tǒng)存在多個(gè)變量且變量處于不同的數(shù)量級時(shí)，由于數(shù)值計(jì)算的誤差，可能會(huì)忽略某些變量的貢獻(xiàn)。因此推薦距離應(yīng)對多個(gè)變量同時(shí)分析，即每個(gè)變量采用獨(dú)立的推薦距離。同時(shí)，在新周期胞的判斷中，當(dāng)使用計(jì)算機(jī)本身的舍入誤差作為評價(jià)標(biāo)準(zhǔn)，雖滿足精度的需求，但是以Matlab為例，其舍入誤差為小數(shù)點(diǎn)后15位，必然導(dǎo)致計(jì)算規(guī)模極大。因此將其歸結(jié)到系統(tǒng)特性維持標(biāo)準(zhǔn)中進(jìn)行分析，提高其舍入誤差，保證計(jì)算速度。

1．1．2 在判斷得到 Gr（B）＝－n時(shí)，即出現(xiàn)處理中的已處理胞，按照原方法，需要對新周期胞組和已確定周期胞組進(jìn)行判斷。出現(xiàn)處理中的已處理胞原因在于不滿足新周期胞組的判斷，尤其是已確定周期胞組的判斷。因此當(dāng)再次映射到該胞時(shí)，只需直接判斷同之前落入該胞中所有點(diǎn)的距離滿足推薦距離或舍入誤差，從而減少計(jì)算步驟。

1．1．3 對相空間中選定的研究區(qū)域邊界處初始點(diǎn)的積分計(jì)算可能映射出定義的研究區(qū)域，從而成為陷胞。而實(shí)際上在若干積分周期后，該點(diǎn)可以回到分析的研究區(qū)域并具有一定的周期性。因此，設(shè)定軌跡離開研究區(qū)域時(shí)，進(jìn)行多步映射并觀察該軌跡是否回到分析的相空間中選定的研究區(qū)域。如回到分析的研究區(qū)域，則繼續(xù)計(jì)算，否則認(rèn)定其進(jìn)入陷胞。此外，在映射步數(shù)方面，當(dāng)某點(diǎn)經(jīng)過一定的限制映射步數(shù)P lim后，仍未形成或達(dá)到某周期胞，則為了保證計(jì)算時(shí)長，可判斷其進(jìn)入陷胞。

1．1．4 為保證計(jì)算效率，可設(shè)定較大的系統(tǒng)特性維持標(biāo)準(zhǔn)，提高計(jì)算速度，但可能引入偽周期胞。因此在得到周期胞后，對周期胞內(nèi)的吸引子進(jìn)行一定的數(shù)值計(jì)算，再進(jìn)行周期胞以及打靶法的判斷。這一點(diǎn)類似于文獻(xiàn)［5］中方法，不同的是在打靶法的過程中，可進(jìn)行穩(wěn)定性的分析。

1.2 周期運(yùn)動(dòng)的穩(wěn)定性分析

對于耗散系統(tǒng)，其相軌線會(huì)隨著時(shí)間延續(xù)而趨于某一穩(wěn)定的吸引子，該吸引子可能為平衡點(diǎn)、周期運(yùn)動(dòng)、概周期運(yùn)動(dòng)乃至混沌運(yùn)動(dòng)，對于平衡點(diǎn)和周期運(yùn)動(dòng)的分析是工程分析的重要內(nèi)容。

受周期激勵(lì)的n維非自治系統(tǒng)，其周期運(yùn)動(dòng)的存在性和數(shù)目較平衡點(diǎn)的求解更為困難，試圖采用打靶法對浮筏隔振系統(tǒng)的周期運(yùn)動(dòng)穩(wěn)定性進(jìn)行分析。打靶法的計(jì)算方法為Newton-Raphson迭代法，其通過構(gòu)造龐加萊截面，分析相軌線在經(jīng)過該截面時(shí)點(diǎn)的分布規(guī)律，從而得到周期運(yùn)動(dòng)以及對該周期運(yùn)動(dòng)穩(wěn)定性的分析。

2 算例分析

2.1 PMUCR方法對高維動(dòng)力學(xué)方程的分析

圖1 浮筏隔振系統(tǒng)空間結(jié)構(gòu)示意圖Fig．1 The sketch of floating raft

引入船舶中的浮筏隔振系統(tǒng)，其空間構(gòu)型以及力和速度條件分別如圖1，2所示。采用子結(jié)構(gòu)方法對該系統(tǒng)進(jìn)行建模分析。系統(tǒng)分為A機(jī)組、B上層隔振器、C筏架、D下層隔振器以及E基礎(chǔ)板五大子結(jié)構(gòu)。其中A子結(jié)構(gòu)為兩臺機(jī)組，是浮筏系統(tǒng)的本質(zhì)要素。由線性條件下的分析可知［11］，在低頻段，筏架的剛?cè)岢潭葘ο到y(tǒng)整體的動(dòng)力學(xué)響應(yīng)影響較小，甚至可忽略。本算例試圖分析在非線性條件下，是否具有類似的結(jié)論。

圖2 浮筏隔振系統(tǒng)力和速度條件簡圖Fig．2 The force and velocity condition of floating raft

設(shè)定隔振器的剛度為立方非線性即彈性恢復(fù)力的表達(dá)式為kx＋Ux3。在豎直方向上，柔性筏架柔性基礎(chǔ)的非線性浮筏系統(tǒng)動(dòng)力學(xué)方程可簡化如下式所示（因位置所限，暫不列出阻尼項(xiàng)）

上式表示系統(tǒng)處于十維的狀態(tài)空間。方程以隔振器原長為坐標(biāo)原點(diǎn)，其中mz分為mz1和mz2，代表兩臺機(jī)組的負(fù)載質(zhì)量，其位移分別表示為x1和x2。作用于機(jī)組的力FAe分別表示為F1和F2，其大小相等。上層隔振器的力為和表示為e臺機(jī)組的第i個(gè)隔振器的上部t和下部b的作用力，表達(dá)式為k1x＋Ux3，下層隔振器的力表達(dá)式為k2x＋Vx3。筏架的質(zhì)量為mzc，其密度、厚度和抗彎剛度分別表示為ρ1、h1和D1*，與此對應(yīng)的基礎(chǔ)板的密度、厚度和抗彎剛度為ρ、h和D。筏架的剛體位移和模態(tài)位移以及基礎(chǔ)板的模態(tài)位移分別為xcg、xcr和xer，表達(dá)式中的大寫字母R均為模態(tài)變換矩陣。

改良后PMUCR方法的流程圖如圖3所示。其計(jì)算流程較簡單胞映射復(fù)雜，分析判斷過程更為細(xì)致。

將剛性筏架和柔性筏架的系統(tǒng)分別進(jìn)行全局分析，并對PMUCR方法在高維系統(tǒng)中的運(yùn)用進(jìn)行討論。系統(tǒng)的物理參數(shù)為：負(fù)載機(jī)組1，2的幾何尺寸為0．8 m×0．8 m×0．4 m，質(zhì)量為1 996．8 kg。為分析多種工況，設(shè)定上層隔振器剛度為0．5×107N／m和1×107N／m。平置板式筏架的幾何參數(shù)為2 m×1 m×0．06 m，質(zhì)量為942 kg，損耗因子為0．01，自由邊界。設(shè)定下層隔振器剛度為1×107N／m。基礎(chǔ)板的幾何參數(shù)為2 m×1 m×0．02 m，損耗因子為0．01，簡支邊界。

圖3 PMUCR方法流程圖Fig．3 Flow chart for the PMUCR algorithm

圖4 激勵(lì)條件下的吸引子Fig．4 Attractor under conditon

依據(jù)研究區(qū)域覆蓋盡可能大的原則，在剛性筏架模型中，以狀態(tài)空間中子空間X＝［x1x′1x2x′2xcg40），x2∈ （－0．16，0．16），x′2∈ （－40，40），x3∈（－0．16，0．16），x′3∈（－40，40），x4∈（－0．16，0．16），x′4∈（－40，40）｝為研究區(qū)域，并均勻分布508個(gè)胞為胞參考空間。為分析最重要的機(jī)組狀態(tài)，選取該研究區(qū)域?qū)?yīng)的50×50個(gè)初始點(diǎn)。與此對應(yīng)的柔性筏架模型以及不同激勵(lì)條件下的模型，均根據(jù)剛性筏架模型得到其子空間。在激勵(lì)1條件和激勵(lì)2條件兩種工況。經(jīng)PMUCR方法計(jì)算得到的吸引子在初值平面內(nèi)的投影如圖4所示。

圖4（a）說明在可能的初始條件下，隔振系統(tǒng)最后可能呈現(xiàn)出兩個(gè)不同的周期運(yùn)動(dòng)，即周期1和周期2運(yùn)動(dòng)狀態(tài)。在兩個(gè)運(yùn)動(dòng)狀態(tài)下，系統(tǒng)的隔振性能不同。例如，當(dāng)任一初始條件使得系統(tǒng)處于周期1運(yùn)動(dòng)，并以功率流落差為評價(jià)標(biāo)準(zhǔn)，則功率流落差為23．61 dB，而由其他初始條件，使系統(tǒng)處于周期2運(yùn)動(dòng)狀態(tài)時(shí)，其功率流落差為14．82 dB。因此系統(tǒng)在周期1運(yùn)動(dòng)狀態(tài)下的隔振效果優(yōu)于周期2運(yùn)動(dòng)，所以當(dāng)船舶隔振系統(tǒng)受到一定條件的外部擾動(dòng)時(shí)，可能發(fā)生運(yùn)動(dòng)狀態(tài)的突變，導(dǎo)致隔振性能發(fā)生變化。在進(jìn)行非線性隔振系統(tǒng)的設(shè)計(jì)時(shí)需要考慮該因素。同時(shí)在激勵(lì)2條件下進(jìn)行功率流落差的分析，其功率流落差在周期1運(yùn)動(dòng)為18．33 dB，其余為15．37 dB左右。同時(shí)由隔振器不同剛度的分析表明，立方非線性剛度的隔振器表現(xiàn)出的周期運(yùn)動(dòng)較為一致。

圖5 吸引子的吸引域分布Fig．5 The domain of attraction

圖6 P－1和P－2周期運(yùn)動(dòng)的時(shí)間歷程及龐加萊截面Fig．6 Displacement and phase trajectory of P－1 and P－2

圖5中“·”區(qū)域?yàn)镻－1周期運(yùn)動(dòng)，“＋”區(qū)域?yàn)镻－2周期運(yùn)動(dòng)，空白區(qū)域?yàn)橄莅Ｓ蓤D4（b）和圖5可見，在相同條件下剛性和柔性筏架的系統(tǒng)，其周期運(yùn)動(dòng)的區(qū)域不同。

由PMUCR方法所得到的P－1和P－2周期運(yùn)動(dòng)處的點(diǎn)為積分的初始點(diǎn)，計(jì)算得到如圖6所示的時(shí)間歷程及龐加萊截面圖。由龐加萊截面分析，其P－1及P－2周期運(yùn)動(dòng)是正確的，而且經(jīng)過PMUCR方法的映射計(jì)算，已經(jīng)消除了系統(tǒng)瞬態(tài)響應(yīng)的影響。

2.2 PMUCR的參數(shù)設(shè)置分析

如果PMUCR方法使用的系統(tǒng)特性保持標(biāo)準(zhǔn)不夠合理，則會(huì)引起一定的問題，若該標(biāo)準(zhǔn)過小且限制映射步數(shù)P lim過少，則導(dǎo)致計(jì)算過程難以收斂，從而產(chǎn)生陷胞及其它問題；當(dāng)限制映射步數(shù)過大，則增加了計(jì)算量。以激勵(lì)1條件下的剛性筏架為例，進(jìn)行討論分析。設(shè)定分析工況如表1所示，設(shè)定不同的舍入誤差以及映射步數(shù)，對PMUCR方法在全局分析中的參數(shù)設(shè)置進(jìn)行研究。

表1 PMUCR分析工況條件Tab．1 The analysis condition of PMUCR

在圖7中的左側(cè)圖為吸引子，右側(cè)圖中以“·”表示陷胞的分布。自上而下分為四列，分別對應(yīng)表1中自上而下四個(gè)條件的計(jì)算結(jié)果。

圖7 吸引子及陷胞的分布Fig．7 The attractor and the sink cell under different condition

由圖7（a），7（b）所示，當(dāng)系統(tǒng)特性保持標(biāo)準(zhǔn)和限制映射步數(shù)較小，導(dǎo)致大量的胞映射進(jìn)入陷胞，而且也未發(fā)現(xiàn)P－1周期胞。圖7（c），7（d）增大系統(tǒng)特性保持標(biāo)準(zhǔn)和限制映射步數(shù)，導(dǎo)致識別出大量的周期胞，改良方法4的計(jì)算表明，除圖4（a）中的周期胞外，其余均為偽周期胞。圖7（e），7（f）減小系統(tǒng)特性保持標(biāo)準(zhǔn)及增大限制映射步數(shù)，偽周期胞數(shù)量減少，陷胞數(shù)量減少。圖7（g），7（h）縮小系統(tǒng)特性保持標(biāo)準(zhǔn)及再次增大限制映射步數(shù)，偽周期胞消失，陷胞數(shù)量較圖7（f）多，并散布于分析平面內(nèi)。原理上，增大系統(tǒng)特性維持標(biāo)準(zhǔn)、延長計(jì)算時(shí)長及增多限制映射步數(shù)可以有效地減少陷胞的發(fā)生，但會(huì)導(dǎo)致結(jié)果精度的降低，尤其大幅增加計(jì)算時(shí)長。由圖7可見，改良方法4較為適用于高維系統(tǒng)的分析。

2．3 周期解的穩(wěn)定性分析

以激勵(lì)2條件下，柔性筏架的周期運(yùn)動(dòng)穩(wěn)定性分析為例，其P－1運(yùn)動(dòng)的特征值最大模為0．989 1，P－2周期運(yùn)動(dòng)的特征值最大模分別為0．990 8、0．990 6、0．991 4及0．982 7，均未超過 1，所以其周期運(yùn)動(dòng)穩(wěn)定。對其它情況分析，也具有類似的性質(zhì)。實(shí)際上，非自治系統(tǒng)周期運(yùn)動(dòng)的穩(wěn)定性受激勵(lì)頻率影響極大，尤其接近于主共振頻率時(shí)，存在不穩(wěn)定區(qū)域使其特征值大于1，而當(dāng)激勵(lì)頻率處于正常工作范圍，其穩(wěn)定性較好。

3 結(jié) 論

通過PMUCR方法對剛性和柔性筏架條件下的浮筏隔振系統(tǒng)進(jìn)行全局分析，研究其筏架對系統(tǒng)全局特性的影響。并對該方法在高維系統(tǒng)中的應(yīng)用進(jìn)行相應(yīng)的改良和討論。得到的主要結(jié)論為：

（1）PMUCR適用于高維系統(tǒng)的全局特性分析，并具有較好的精度。其計(jì)算量較胞映射大，計(jì)算結(jié)果的精度同時(shí)受到系統(tǒng)特性保持標(biāo)準(zhǔn)以及計(jì)算步數(shù)的影響，因此在使用的過程中需要同時(shí)兼顧計(jì)算量、特性保持標(biāo)準(zhǔn)以及映射步數(shù)三者的關(guān)系。

（2）剛性和柔性筏架浮筏隔振系統(tǒng)的運(yùn)動(dòng)現(xiàn)象存在差異。PMUCR方法的計(jì)算表明，剛性和柔性筏架隔振系統(tǒng)所發(fā)生周期運(yùn)動(dòng)的類型較為相似，但其相對應(yīng)區(qū)域是不同的。該計(jì)算切實(shí)地證明了在非線性條件下的所有頻段內(nèi)，筏架的特性對系統(tǒng)的動(dòng)力學(xué)特性影響顯著，而不僅僅局限于線性系統(tǒng)的中、高頻段。

（3）非線性浮筏隔振系統(tǒng)的全局特性在不同的激勵(lì)條件以及初始條件下呈現(xiàn)出復(fù)雜的周期現(xiàn)象。在分析平面內(nèi)發(fā)現(xiàn)了較為典型的單周期運(yùn)動(dòng)和二倍周期運(yùn)動(dòng)。當(dāng)處于不同的周期運(yùn)動(dòng)時(shí)，系統(tǒng)的隔振效果存在差異。

［1］Hsu C S．A Theory of cell-to-cell mapping for dynamical system［J］．ASME，Journal of Applied Mechanics，1980，47：931－936．

［2］Hsu C S．Cel1-to-Cell mapping：A method of global analysis for nonlinear systems［M］．New York，Berlin，Heidelberg，Vienna，Tokyo：Springer-Verlag，1987：1－352．

［3］Tongue B H，Gu K．Interpolated cell mapping of dynamic system［J］．ASME，Journal of Applied Mechanics，1988，55（6）：461－466．

［4］Levitas J，Weller T，Singer J．Poincare-like simple cell mapping for non-linear dynamical systems［J］．Journal of Sound and Vibration，1994，176（5）：641－662．

［5］文成秀，姚玉璽，聞邦椿．動(dòng)力系統(tǒng)的胞映射－點(diǎn)映射法及應(yīng)用［J］．東北大學(xué)學(xué)報(bào)，1997，18（2）：119－122．WEN Cheng-xiu， YAO Yu-xi， WEN Bang-chun． Cellmapping point-mapping method and its application for dynamic system［J］．Journal of Northeastern University，1997，118（2）：119－122．

［6］江俊，徐健學(xué)．動(dòng)力系統(tǒng)的胞化積分軌跡法［J］．振動(dòng)工程學(xué)報(bào)，1993，6（2）：194－198．JIANG Jun，XU Jian-xue．Cellurated integration method of dynamical systems［J］．Journal of Vibration Engineering，1993，6（2）：194－198．

［7］Jiang Jun，Xu Jiang-xue．A method of point mapping under cell reference for global analysis of nonlinear dynamical systems［J］．Physics Letters A，1994，188：137－145．

［8］俞翔，朱石堅(jiān)，劉樹勇．PMUCR方法在高維非線性動(dòng)力學(xué)系統(tǒng)中的應(yīng)用［J］．應(yīng)用力學(xué)學(xué)報(bào)，2006，23（6）：232－236．YU Xiang，ZHU Shi-jian，LIU Shu-yong．PMUCR method with application to high-dimensional nonlinear system［J］．Chinese Journal of Applied Mechanics，2006，23（6）：232－236．

［9］Moresco P，Dawson S P．The PIM-simplex method：an extension of the PIM-triple method to saddles with an arbitrary number of expanding directions．Physica D，1999，126：38－48．

［10］高文輝，郭旭，江俊．基于胞參考點(diǎn)映射法對轉(zhuǎn)子／定子碰摩響應(yīng)的全局分析［J］．動(dòng)力學(xué)與控制學(xué)報(bào)，2011，9（4）：363－367．GAO Wen-hui，GUO Xu，JIANG Jun．Study on the global response characteristics of a rotor／stator contact system through PMUCR［J］．Journal of Dynamics and Control，2011，9（4）：363－367．

［11］張樹楨，陳 前．柔性浮筏隔振系統(tǒng)的理論建模和仿真研究［J］．船舶力學(xué)，2012，16（10）：1187－1198．ZHANG Shu-zhen，CHEN Qian．Modeling and simulation research on flexible floating raft isolation system［J］．Journal of Ship Mechanics，2012，16（10）：1187－1198．