帶支腿浮式結(jié)構(gòu)水動力系數(shù)快速算法研究

武海浪，陳徐均，黃亞新，沈海鵬，苗玉基

（解放軍理工大學(xué) 野戰(zhàn)工程學(xué)院，南京210007）

0 引 言

帶支腿浮式結(jié)構(gòu)是指在浮式平臺的基礎(chǔ)上設(shè)計加裝一定數(shù)量的可升降支腿而形成的海洋工程結(jié)構(gòu)。結(jié)構(gòu)在提樁和放樁過程中，樁腿并不觸及海底，主體平臺處于帶樁腿漂浮的狀態(tài)。對于一個確定的浮式結(jié)構(gòu)，按照船舶水彈性理論對其進(jìn)行水動力分析時，一般需要進(jìn)行濕面網(wǎng)格定義、格林函數(shù)分析和水動力響應(yīng)計算幾個步驟[1-5]。要研究水下樁腿長度變化以后的浮式結(jié)構(gòu)的水彈性性質(zhì)，就要針對變化的參數(shù)，重新進(jìn)行浮式結(jié)構(gòu)的干模態(tài)參數(shù)計算、重新定義新的濕面網(wǎng)格，然后進(jìn)行水動力計算。文獻(xiàn)[6-8]就是采用了根據(jù)水下樁腿長度分工況的方法進(jìn)行分別計算，認(rèn)為水下樁腿的長度對浮式結(jié)構(gòu)的水動力系數(shù)、波浪激勵力以及水動力響應(yīng)都有顯著的影響。

三維線性水彈性分析的核心是格林函數(shù)分析，是比較消耗CPU時間的，尤其是當(dāng)參與計算的面元較多的時候，計算效率會明顯地下降。前人對三維水彈性分析算法效率給予了關(guān)注，致力于提高計算效率：戴愚志等[9-10]用GMRES（Generalized Minimal RESidual algorithm）方法求解了大型離岸結(jié)構(gòu)水彈性分析所獲得的復(fù)系數(shù)線性方程組，并按照向后誤差分析方法給出了它的終止準(zhǔn)則。數(shù)值試驗(yàn)表明GMRES方法可提高效率，減少計算時間，優(yōu)于直接方法。Wang等[11]引入兩項(xiàng)技術(shù)到三維水彈性計算當(dāng)中，以提高超大型浮式結(jié)構(gòu)水彈性計算的效率。其一是使用新的截斷標(biāo)準(zhǔn)來降低格林函數(shù)及其導(dǎo)數(shù)的分析時間，其二是使用了交互式稀疏求解方法代替高斯方法求解水彈性方程，結(jié)果顯示有效地提高了計算效率。Wu[12]為了提高水彈性程序的計算效率在其程序中使用了對稱復(fù)合勢方法，這種分析方法也一直沿用至今，所以在后面的分析中，本文也只建立了一半的濕面單元。

浮式結(jié)構(gòu)提樁放樁時只是部分參數(shù)在確定的狀態(tài)下進(jìn)行了微調(diào)，并未發(fā)生本質(zhì)的變化。前面的研究方案盡管直觀、易行，但這種參數(shù)微調(diào)就不得不重新建模、反復(fù)分析和重新計算，這似乎并不夠高效，每一次新的計算包含了大量的重復(fù)工作。若采用這樣的計算方案，在這類浮式結(jié)構(gòu)的設(shè)計工作中，或許將降低工作效率。本文在研究者所熟知的三維線性水彈性理論體系下，利用格林函數(shù)矩陣的性質(zhì)，針對帶支腿浮式結(jié)構(gòu)變參數(shù)的水動力系數(shù)計算問題提出了一種相對較高效的分析方法，主要是在格林函數(shù)分析和線性方程組求解方面作了相關(guān)的研究。

1 算法原理

1.1 變物面條件問題

假設(shè)彈性浮體周圍為均勻不可壓縮、無粘的理想流體，流場的運(yùn)動是無旋的，自由表面波為微幅且?guī)е雀∈浇Y(jié)構(gòu)在波浪中平衡位置附近作微幅運(yùn)動，則流場的運(yùn)動可以在歐拉坐標(biāo)系中采用三維勢流理論來描述。

根據(jù)對流場的基本假定，在空間固定坐標(biāo)系中帶支腿浮式結(jié)構(gòu)周圍，如圖1所示，邊界SL+SF+S∞+SB所包圍的流場內(nèi)輻射勢定解問題可以描述為：

其中：▽2為拉普拉斯算子，φr表示輻射波速度勢的空間分量，當(dāng)浮體在平衡位置附近以入射波圓頻率ω，以第r階干模態(tài)振型作單位主坐標(biāo)幅值振蕩時所誘導(dǎo)的流場周圍的流體運(yùn)動速度勢。SL為物面邊界條件，表示帶支腿浮式結(jié)構(gòu)的濕面積，主要由浮式結(jié)構(gòu)平臺的底部、側(cè)舷和樁腿三個部分組成。由于樁腿是可以上下移動的，物面條件是可變的，物面條件的改變將導(dǎo)致方程解的改變。SF表示自由面邊界條件，S∞為無窮遠(yuǎn)邊界條件，SB為海底邊界條件。

1.2 計算方法

假定：（1）帶支腿浮式結(jié)構(gòu)中樁腿是平直的柱狀結(jié)構(gòu)；（2）忽略因?yàn)闃锻鹊纳仙蛘呦路哦鴮?dǎo)致的主體平臺水線面的變化；（3）樁腿底板的影響很小。選取帶支承腿浮式結(jié)構(gòu)水下樁腿最長L1時為初始狀態(tài)，即將樁腿盡可能向下放置但不觸及海底的狀態(tài)，物面條件記為SL1。采用格林函數(shù)法對結(jié)構(gòu)的初始狀態(tài)進(jìn)行水動力分析。

對初始狀態(tài)進(jìn)行水動力建模。圖2為任意一帶支腿浮式結(jié)構(gòu)，此時樁腿處于水下伸長狀態(tài)，濕面網(wǎng)格分為船身主體部分和樁腿部分。整個結(jié)構(gòu)一共劃分的單元數(shù)記為N，主體平臺單元編號記為1～Nhull，后面Nhull+1～N個編號給樁腿。對樁腿濕面元編號原則參照表1，樁腿濕面元編號的方法總體上遵循了以水深為標(biāo)尺，編號大小遵循“由上到下”的原則。即離水線面向下越遠(yuǎn)的面元，編號越靠后，豎直方向上，處于同一高度的面元編號越為接近。

表1 樁腿狀態(tài)編號及濕面元編號Tab.1 The numbering schedule of the wetted panels

樁腿處于初始狀態(tài)時，經(jīng)開拓后場內(nèi)任一點(diǎn)P的輻射勢可按源分布表示為：

其中：G（ P， ）Q 為格林函數(shù)，應(yīng)由初始物面條件SL1決定。使用文獻(xiàn)[12]給出的有限水深無航速問題的三維脈動源格林函數(shù)，其級數(shù)表達(dá)形式為：

積分方程離散為一組線性代數(shù)方程組，表示為：

上式可寫成矩陣形式，有：

其中：

當(dāng)樁腿同時緩慢、均勻地由初始狀態(tài)長度L1恰好抬升到Li時，從實(shí)際效果上看是樁腿的水下長度變短。若重新進(jìn)行水動力建模以及新的格林函數(shù)分析，可得到新的格林函數(shù)矩陣Gnew和新的影響矩陣Cnew。

其中：g′和c′分別代表樁腿底部的面元的影響。基于假定（3），可以認(rèn)為樁腿底板的作用很小，忽略矩陣中g(shù)′和 c′在矩陣中的影響，即：

這樣，當(dāng)樁腿由初始狀態(tài)連續(xù)地作抬升動作時，可以視為格林函數(shù)矩陣規(guī)模不斷縮小的過程，即連續(xù)地有矩陣中相應(yīng)的行列退出計算。如圖4（b）所示，根據(jù)編號法則，位于樁腿底部面元的行列編號靠后，當(dāng)樁腿由初始狀態(tài)作抬升動作時，等效于樁腿底部面元退出計算，相應(yīng)的就是將格林函數(shù)矩陣中靠右的部分一層一層地剝離矩陣。直到樁腿全部抬升起來，等同于編號為Nhull+1～Nn-1的網(wǎng)格退出計算。這樣有規(guī)律的按編號進(jìn)行矩陣縮減可以得到一系列格林函數(shù)矩陣，實(shí)際上是給出了流體域內(nèi)速度勢定解問題中物面由SL1逐步變化到SLn的方程組的近似格林函數(shù)解?？捎霉奖硎緸椋?/p>

格林函數(shù)及其影響矩陣的變化過程表示為

在濕面上進(jìn)行源強(qiáng)積分計算，可得到輻射勢：

以上便為帶支腿浮式結(jié)構(gòu)水動力系數(shù)快速計算的原理。根據(jù)算法原理，編寫程序，成為帶支腿浮式結(jié)構(gòu)水動力系數(shù)快速計算子模塊。

2 算例介紹

選取某型帶支腿浮式結(jié)構(gòu)為研究對象，結(jié)構(gòu)長度約為130 m，寬38 m，型深8 m，吃水4 m。將主體平臺簡化為一個浮箱，6個升降腿分布于全船兩側(cè)，每個樁腿高68 m，截面為4 m×4 m方形。使用格林函數(shù)分析方法，分兩種方案計算水下樁腿長度為10 m，20 m，30 m，40 m和50 m時該型帶支腿浮式結(jié)構(gòu)的水動力系數(shù)。

第一種方案采用分工況法計算。按樁腿長度建立5種工況分別進(jìn)行水動力建模。如圖5所示，最短為工況1，最長為工況5（工況2～4圖略）。所研究的浮式結(jié)構(gòu)是對稱結(jié)構(gòu)，只需要定義一半濕面單元，即可進(jìn)行水彈性分析。這樣，工況1就有265個面元；工況2時有289個面元；工況3為313個面元；工況4為337個面元；工況5為361個面元。使用水彈性程序?qū)ξ宸N工況分別進(jìn)行數(shù)值計算。

第二種方案采用以矩陣縮減為核心的快速算法計算。以浮式結(jié)構(gòu)的樁腿向水下伸出50 m，即工況5為初始狀態(tài)，快速計算當(dāng)水下樁腿長度為40 m，30 m，20 m，10 m時的水彈性響應(yīng)。

表2 樁腿狀態(tài)編號及濕面元編號Tab.2 The numbering schedule of the wet panels

整個浮式結(jié)構(gòu)處于初始狀態(tài)時共有361個濕面元，對已經(jīng)劃分好的濕面元網(wǎng)格進(jìn)行編號，浮式結(jié)構(gòu)主體部分編號為1～239，對樁腿部分面元主要按照由淺及深，由大及小的原則。表2說明了算例中對樁腿編號的規(guī)律。處于樁腿50 m處的濕面元組的編號為355～361，其狀態(tài)定為1；處于樁腿10 m處的面元組編號為239～262，其狀態(tài)編號為5。

3 結(jié)果檢驗(yàn)

3.1 計算精度分析

如圖6～11，在坐標(biāo)平面上，有一系列的散點(diǎn)，每一個點(diǎn)的橫坐標(biāo)表示快速算法結(jié)果，縱坐標(biāo)表示普通算法結(jié)果，通過衡量圖中點(diǎn)與直線y=x的位置來判斷快速算法結(jié)果和一般算法結(jié)果的相符程度，其中（a）圖為附加質(zhì)量，（b）圖為附加阻尼。每一幅圖中各有四種顏色的點(diǎn)，分別代表帶支腿浮式結(jié)構(gòu)的樁腿處于四種長度，在圖的右下方的方框內(nèi)進(jìn)行了說明。每一幅圖左上方標(biāo)有一個4×3的數(shù)據(jù)矩陣。其中一行代表一個固定的樁腿水下長度；一列代表兩組計算結(jié)果定量的誤差分析結(jié)果，共有三列，分別為平均相對誤差（EAVE），均方相對誤差（RMVE）和相關(guān)系數(shù)（CC）。其計算方法為：

對比分析通過快速算法得到的前六階剛體模態(tài)的水動力系數(shù)。如圖6，圖7和圖11所示，縱蕩、橫蕩和艏搖的水動力系數(shù)散點(diǎn)與直線y=x重合程度極高，這說明對這三種剛體模態(tài)來說，快速算法結(jié)果和普通算法結(jié)果相符程度極高。計算得到的平均相對誤差、均方相對誤差和相關(guān)系數(shù)上也證實(shí)了這樣的結(jié)論，這三種模態(tài)的平均相對誤差和均方相對誤差都在0.5%以下，而兩者之間的相關(guān)系數(shù)都在99%以上，可以認(rèn)為采用快速算法得到的結(jié)果與普通算法結(jié)果一致；

和前面這三個剛體模態(tài)相比，能夠看到垂蕩、橫搖和縱搖三個模態(tài)相對應(yīng)圖中的散點(diǎn)相對直線y=x有明顯的偏移，如圖8，圖9和圖10所示，這表明快速算法結(jié)果和普通算法結(jié)果有一定的誤差存在，從數(shù)據(jù)上看，這幾個模態(tài)的平均相對誤差、均方相對誤差幾乎都在1%以上，樁腿處于10 m時的橫搖附加阻尼系數(shù)誤差達(dá)到了7%，如圖9（b）所示，誤差明顯增高；注意到圖中的散點(diǎn)位置高于直線y=x，如圖10（a）所示，這說明，和普通算法相比較，快速算法在一定程度上低估了橫搖、縱搖和垂蕩的附加質(zhì)量和附加阻尼。

同時，這三幅圖中誤差統(tǒng)計量表明，快速算法結(jié)果的誤差隨著樁腿長度變短而增大，如圖9（b）中的數(shù)據(jù)，當(dāng)樁腿為40 m時，平均相對誤差為1.4597%和0.3421%，當(dāng)樁腿長度處于10 m狀態(tài)時，這兩個數(shù)據(jù)增加到7.1231%和1.4454%。原因在于普通算法中針對樁腿不同長度時的浮式結(jié)構(gòu)分別進(jìn)行水動力建模，樁腿的底部加了底板，使結(jié)構(gòu)物封閉。而在快速算法中，在進(jìn)行矩陣拆解時，這部分濕面元沒有參與計算，不包含樁腿底板的信息。當(dāng)樁腿長度較長時，浮式結(jié)構(gòu)的整體濕面積相對較大，底板的影響偏?。欢?dāng)樁腿長度為10 m時，盡管樁腿底板的面積并沒有變化，但是浮式結(jié)構(gòu)的整體濕面積變小了，樁腿底板在濕面積中占有的比例變大，因此底板影響就顯現(xiàn)出來了。

3.2 計算效率分析

表3給出了分工況方法和快速算法耗費(fèi)的CPU時間統(tǒng)計。在算例中一共分析了200個頻率，20個模態(tài)（包含14個彈性體模態(tài)）。

表3 CPU計算時間統(tǒng)計表Tab.3 Statistics of CPU time

當(dāng)使用分工況法進(jìn)行分析時，計算工況5耗時67 min；從工況1到工況4的計算時間和為208 min；當(dāng)使用快速算法進(jìn)行分析時，工況5被選為初始狀態(tài)，分析時間與傳統(tǒng)方法一致，但工況4到工況1的結(jié)果為一次性算出，消耗CPU時間為16 min，和前一種分析方法相比較，一共節(jié)約192 min，效率提高了92.3%。這還只是統(tǒng)計的在計算機(jī)上的絕對計算時間，尚不包括網(wǎng)格劃分、參數(shù)定義等前期手工準(zhǔn)備時間。計算效率有顯著的提升。效率提高的主要原因在于后續(xù)的分析主要是在進(jìn)行矩陣的縮減，實(shí)際上是在對電腦中已有的數(shù)據(jù)進(jìn)行行變換和列變換，而之前方法的計算主要是進(jìn)行格林函數(shù)分析和格林函數(shù)偏導(dǎo)數(shù)的分析，這樣節(jié)約出了大量的機(jī)時。若濕面網(wǎng)格劃分更多，格林函數(shù)矩陣的規(guī)模更大，計算效率或許將進(jìn)一步提高，算法的優(yōu)勢將進(jìn)一步體現(xiàn)。

4 結(jié) 論

本文在三維勢流理論體系下，利用格林函數(shù)矩陣的性質(zhì)，針對帶支腿浮式結(jié)構(gòu)變參數(shù)的水動力系數(shù)計算問題提出了一種相對較高效的方法。將求解格林函數(shù)的過程視為格林函數(shù)矩陣G和影響矩陣C矩陣規(guī)模減小的過程，即水下樁腿任意長度為Li時的格林函數(shù)矩陣和影響矩陣中包含著L＜Li時的信息。這樣在進(jìn)行下一步計算的時候，可以通過上述已經(jīng)確定下來的編號找到退出計算的網(wǎng)格，對其進(jìn)行消除，而將已有的信息保留下來，而不必重復(fù)計算，從而達(dá)到簡化計算的目的。本文還給出了通過快速算法計算結(jié)果的誤差和快速算法的效率。

結(jié)果表明：快速算法所計算的剛體的縱蕩、橫蕩和艏搖三種水動力系數(shù)與普通算法結(jié)果一致，可完全反應(yīng)這些模態(tài)的水動力特性；快速算法所計算的垂蕩、橫搖和縱搖三種剛體模態(tài)的結(jié)果與普通算法有較小的誤差，但基本可以反應(yīng)這些模態(tài)的水動力特性；當(dāng)使用格林函數(shù)矩陣縮減法進(jìn)行分析以后，算例的計算效率可以提高92%以上，非?？捎^。