于龍飛，張亞群

（1. 廣州航海學(xué)院，海運(yùn)學(xué)院，廣州 510725；2. 中國(guó)科學(xué)院廣州能源研究所，中國(guó)科學(xué)院可再生能源重點(diǎn)實(shí)驗(yàn)室，廣州 510640）

海浪譜的二維仿真及驗(yàn)證*

于龍飛1?，張亞群2

（1. 廣州航海學(xué)院，海運(yùn)學(xué)院，廣州 510725；2. 中國(guó)科學(xué)院廣州能源研究所，中國(guó)科學(xué)院可再生能源重點(diǎn)實(shí)驗(yàn)室，廣州 510640）

海浪仿真是海洋工程的熱點(diǎn)話題。本文總結(jié)了一種可以適應(yīng)于各種隨機(jī)波浪譜來(lái)仿真海浪環(huán)境的方法，利用C語(yǔ)言將該方法轉(zhuǎn)換為過(guò)程簡(jiǎn)單、編寫方便的計(jì)算機(jī)程序得到大量仿真波浪數(shù)據(jù)，并作為造波機(jī)輸入數(shù)據(jù)，通過(guò)實(shí)測(cè)數(shù)據(jù)檢驗(yàn)仿真結(jié)果。試驗(yàn)證明，采用文中提出的數(shù)值方法，采集模擬波高數(shù)據(jù)利用相關(guān)函數(shù)法進(jìn)行頻譜估計(jì)，得到的模擬譜和靶譜誤差較小。該數(shù)值仿真方法在海浪仿真上具有實(shí)際應(yīng)用價(jià)值。

海浪譜；隨機(jī)數(shù)；靶譜；頻譜估計(jì)

0 引 言

海洋工程的研究基礎(chǔ)建立在對(duì)海浪的研究之上，海浪的仿真為在實(shí)驗(yàn)室展開(kāi)海洋工程研究提供了基本條件。海浪仿真的準(zhǔn)確度越高，實(shí)驗(yàn)室試驗(yàn)與實(shí)海況試驗(yàn)越接近，越有利于提高試驗(yàn)?zāi)M的精度、試驗(yàn)數(shù)據(jù)的可信度。海浪是不規(guī)則的、隨機(jī)的，通常以海浪譜的形式描述海浪的隨機(jī)過(guò)程，以達(dá)到統(tǒng)一的標(biāo)準(zhǔn)便于建模仿真等相關(guān)的研究分析。經(jīng)過(guò)長(zhǎng)期的實(shí)際觀測(cè)研究，研究人員已經(jīng)總結(jié)出各種海浪譜的經(jīng)驗(yàn)公式。這些海浪譜通常是以頻率為自變量的函數(shù)，無(wú)法直接應(yīng)用于實(shí)驗(yàn)室造波試驗(yàn)中。把這些海浪譜由頻域轉(zhuǎn)化為時(shí)域，再應(yīng)用于海浪仿真中，是目前常用的海浪仿真方法[1]。但此類方法大多比較繁雜且具有局限性，無(wú)法適應(yīng)于不同的海浪譜。本文概括了一種系統(tǒng)的海浪數(shù)值仿真方法，克服了傳統(tǒng)的波浪模擬中采用單一海浪譜為目標(biāo)、試驗(yàn)過(guò)程中多次輸入多次生成大量模擬數(shù)據(jù)所導(dǎo)致的局限性，改進(jìn)了循環(huán)嵌套程序造成的運(yùn)算量過(guò)大的問(wèn)題。該方法將波譜分析?波譜模擬?波浪生成集成為一個(gè)整體，適應(yīng)于各種類型波浪譜的模擬，加快了模擬的速度，提高了模擬精度。

本文雖然以PM譜為例，但各個(gè)參數(shù)的取值范圍較寬，使該方法不僅適用于成熟的波浪譜，而且適合一般海域的觀測(cè)譜。采用計(jì)算機(jī)C程序語(yǔ)言使仿真過(guò)程程序化之后，該方法將更加方便。采集模擬海浪的數(shù)據(jù)，采用相關(guān)函數(shù)法進(jìn)行頻譜估計(jì)，驗(yàn)證仿真結(jié)果與實(shí)際值的一致性。

1 海浪數(shù)值仿真方法[2-6]

1.1 海浪譜

經(jīng)分析，總結(jié)出國(guó)際上規(guī)定的海浪譜都符合以下公式[2]：

其中，指數(shù)p通常取4～6，q為2～4，A和B中包含風(fēng)要素（風(fēng)速、風(fēng)時(shí)、風(fēng)距）或波要素（波高、周期）作為參量。海浪可以看作一個(gè)各態(tài)歷經(jīng)的平穩(wěn)隨機(jī)過(guò)程。通?？梢园岩唤M不規(guī)則波浪分解成多個(gè)不同周期、不同波長(zhǎng)、不同波幅和不同隨機(jī)相位的相同傳遞方向的余弦波疊加而成。式（2）為任意時(shí)刻波面相對(duì)靜水面高度的表達(dá)式：

其中：ai為t時(shí)刻第i個(gè)組成波的振幅；ki、ωi為t時(shí)刻第i個(gè)組成波的波數(shù)和圓頻率；x為t時(shí)刻波面坐標(biāo)軸的橫坐標(biāo)值，通常x軸的正方向沿著波浪傳播方向；εi為t時(shí)刻第i個(gè)組成波的初相位，在頻域與時(shí)域轉(zhuǎn)換時(shí)，一般推薦取分布于（0, 2π）范圍內(nèi)的隨機(jī)數(shù)。

1.2 波浪疊加法

用波浪疊加法進(jìn)行海浪仿真，主要是分解波浪的頻譜，應(yīng)具有以下幾個(gè)步驟：

（1）頻率范圍的選取。根據(jù)各方面資料顯示，海浪模擬的頻率范圍一般以表1為參考。表中列出的仿真頻段為欲模擬的對(duì)象譜（靶譜）S(ω) 能量集中頻段，其余部分可以忽略不計(jì)。根據(jù)試驗(yàn)要求的有義波高H1/3選定仿真頻率段，確定最低頻率（ωL）、最高頻率（ωH）的值。通常頻率選擇越寬，模擬的波浪所包含的能量越大，仿真結(jié)果越逼近實(shí)際值。但是頻率太寬會(huì)造成計(jì)算量過(guò)大，計(jì)算速度較慢，且在離譜峰頻率較遠(yuǎn)的頻率區(qū)域中的數(shù)據(jù)對(duì)結(jié)果影響較小。

表1 不同海況仿真頻段Table 1 Simulation frequency zones under different sea conditions

（2）劃分頻率區(qū)間。一般劃分頻率區(qū)間的方法分為等分頻率和等分能量?jī)煞N方法，這里采用等分能量法。定義累積譜為：

（一）預(yù)防措施 仔豬副傷寒的預(yù)防，應(yīng)特別注意初春時(shí)節(jié)溫度變化不定，做好仔豬圈舍的保暖，中午仔豬可適當(dāng)曬太陽(yáng)；圈舍保持清潔干燥，食槽要經(jīng)常洗刷，糞便發(fā)酵處理。加強(qiáng)飼養(yǎng)管理，初生仔豬應(yīng)爭(zhēng)取早吃初乳，并提前補(bǔ)料，以防亂吃臟物，斷奶分群時(shí)，不要突然改變環(huán)境，豬群盡量分小一些。仔豬斷奶前后（建議28日齡），口服仔豬副傷寒弱毒凍干苗進(jìn)行預(yù)防。較大規(guī)模養(yǎng)豬場(chǎng)發(fā)病后，應(yīng)將病豬隔離治療，豬舍徹底消毒。未發(fā)病的豬可用藥物預(yù)防，在每噸飼料中加入金霉素100 g或氟苯尼考（50～100 ppm），可起較好的預(yù)防作用。

等分能量法將所模擬的譜的能量等分為M份，即分為M個(gè)波組成，由i=1, 2, ··, M表示。分界頻率ωi可以由下式計(jì)算：

其中，m0為譜的零階矩，m0=A/(4B)。這里要求對(duì)象譜S(ω)為可積分，如P?M譜和B?M譜等，由式（1）和式（3），p=5，q=4，可得：

將E(∞)=A/(4B)和式（5）代入式（4）中得：

當(dāng)譜形無(wú)法積分時(shí)，通過(guò)數(shù)值積分計(jì)算累積譜E(ω)，得出一個(gè)關(guān)于 ωi的算式f(ωi)，與式（4）一起得到f(ωi)=im0/M，解出此式得ωi。

為了避免長(zhǎng)時(shí)間的模擬波形可能出現(xiàn)重復(fù)現(xiàn)象，第i區(qū)間的代表頻率可用下面算式計(jì)算：

式中，i=1, 2, ··, M，表示第i個(gè)頻率的波；q為分布在（0, 1）之間的隨機(jī)數(shù)。

（3）初相位的確定。隨機(jī)初相位εi應(yīng)在（0，2π）之間隨機(jī)分布。計(jì)算式為：

式中f也為（0, 1）之間的隨機(jī)數(shù)，初相位的隨機(jī)性決定了波浪重復(fù)性，故隨機(jī)數(shù)f的選取要注意質(zhì)量。

由上式可知，即當(dāng)M=50、n=1時(shí)，可得50個(gè)不同的ai、50個(gè)不同的εi，由此得到50個(gè)不同的ηi(Δt)。累加這50個(gè)ηi(Δt)，即可得到t=Δt時(shí)的波面高度η(Δt)。依此類推可以計(jì)算出不同時(shí)刻的波高值，輸出給造波機(jī)就可以模擬出靶譜。

1.3 編程實(shí)現(xiàn)

上述對(duì)海浪譜仿真過(guò)程的數(shù)值分析，可以通過(guò)調(diào)用C語(yǔ)言函數(shù)庫(kù)中現(xiàn)成的函數(shù)編程得到基本實(shí)現(xiàn)。某些情況下，只需做一些數(shù)學(xué)處理。例如，調(diào)用函數(shù)int random（int num），可以得到0～num之間的隨機(jī)整數(shù)。如要得到上述要求的0～1之間的隨機(jī)數(shù)，程序可以寫成：

在對(duì)不同海浪譜的海浪仿真時(shí)，由于對(duì)頻譜的處理過(guò)程都用波浪疊加法處理頻率計(jì)算，因此在編寫一個(gè)外部函數(shù)供不同頻譜仿真時(shí)調(diào)用，將大大減少計(jì)算機(jī)的計(jì)算量。并且，數(shù)值分析得到的波浪時(shí)域上的仿真數(shù)據(jù)（即模擬時(shí)輸入波高數(shù)據(jù)）可以存儲(chǔ)在同一個(gè)動(dòng)態(tài)數(shù)組當(dāng)中，節(jié)省了程序的內(nèi)部存儲(chǔ)空間。以PM譜為例介紹程序結(jié)構(gòu)[7]。

2 仿真結(jié)果的驗(yàn)證方法

仿真結(jié)果的驗(yàn)證是對(duì)仿真結(jié)果的定量和定性分析，確定其與理論結(jié)果一致性的程度[8]，同樣也是對(duì)仿真數(shù)值分析方法正確與否的一種檢驗(yàn)。本文使用上述仿真方法造波，用浪高儀在定點(diǎn)位置實(shí)時(shí)采集波高數(shù)據(jù)，再進(jìn)行波譜分析得到模擬的波譜，通過(guò)與靶譜進(jìn)行對(duì)比分析得到驗(yàn)證結(jié)論。

通常采用的譜分析有多種方法，如相關(guān)函數(shù)法、快速傅里葉變換法（FFT）、最大熵法等譜估計(jì)的方法。由于FFT需要采集的數(shù)據(jù)實(shí)際是一個(gè)完整的周期，這在實(shí)際操作中很難達(dá)到。而最大商法是建立在相關(guān)函數(shù)法的基礎(chǔ)上，雖然相關(guān)函數(shù)法需要在較大區(qū)間內(nèi)作積分求解，由于計(jì)算機(jī)的計(jì)算能力較強(qiáng)，在研究初期本文采用相對(duì)簡(jiǎn)單的相關(guān)函數(shù)法。該方法得到的頻譜光滑度不高，利用加窗平滑處理，提高譜的質(zhì)量。一般加窗處理有兩種選擇，一種為哈明（Hamming）窗，另一種為哈寧（Hanning）窗。

3 結(jié)果驗(yàn)證及分析

圖1 隨機(jī)海浪仿真結(jié)果Fig. 1 Simulation results of random wave

圖2為實(shí)際采集的模擬波高數(shù)據(jù)，采集頻率為50 Hz，采集時(shí)間為200 s。根據(jù)圖中的實(shí)測(cè)波浪數(shù)據(jù)統(tǒng)計(jì)，H1/3=0.219 9 m，誤差為0.019 9 m。

圖2 采集海浪波高數(shù)據(jù)Fig. 2 Gathering data of wave height

采集數(shù)據(jù)后進(jìn)行波譜分析得到波譜曲線（藍(lán)色曲線）與靶譜曲線（紅色曲線）比較如圖3所示。圖中橫坐標(biāo)表示圓頻率（w=2πf），縱坐標(biāo)表示譜密度函數(shù) [S(ω)]。

圖3 模擬譜與靶譜比較Fig. 3 comparison of simulate spectrum and target spectrum

對(duì)比圖3中的兩條曲線，模擬的海浪譜與靶譜比較相近。靶譜的譜峰密度為0.001 108 m2s、譜峰頻率為2.355 rad/s，實(shí)測(cè)譜譜峰密度為0.001 255 m2s、譜峰頻率為2.355 rad/s。兩者譜峰頻率沒(méi)有誤差，但是譜峰密度之間差距稍大。其中，造波機(jī)的造波精度、采樣點(diǎn)數(shù)有限、浪高儀采樣精度等都是造成誤差存在的主要原因。利用迭代的數(shù)值處理方法，可以減小兩者之間的誤差，并隨著迭代次數(shù)的增多，誤差將越小。

4 結(jié) 論

本文總結(jié)了一套條理清晰、步驟明確的方法仿真海浪場(chǎng)環(huán)境，可以應(yīng)用于絕大部分海浪譜，具有通用性，同時(shí)利用C語(yǔ)言將該方法編寫成通用的程序，又使海浪譜的頻域至?xí)r域轉(zhuǎn)化可以利用計(jì)算機(jī)來(lái)完成，使整個(gè)仿真過(guò)程更加便捷精確，為了進(jìn)一步驗(yàn)證仿真結(jié)果，采集模擬波高數(shù)據(jù)后通過(guò)頻譜分析，對(duì)比模擬譜和靶譜，驗(yàn)證了該方法的準(zhǔn)確性。

[1] 鄭茂琦, 馬春翔, 王志波, 等. 基于海浪譜的海浪模擬的改進(jìn)[J]. 系統(tǒng)仿真學(xué)報(bào), 2014, 36(2): 369-375.

[2] 俞聿修. 隨機(jī)海浪及其工程應(yīng)用[M]. 大連: 大連理工大學(xué)出版社. 2000.

[3] 王立軍. 船舶動(dòng)力定位海浪環(huán)境的實(shí)時(shí)仿真與海浪譜分析[J]. 華東船舶工業(yè)學(xué)院學(xué)報(bào), 2001, 15(1): 48-51.

[4] 李遠(yuǎn)林. 近海結(jié)構(gòu)水動(dòng)力學(xué)[M]. 肇慶: 廣東理工大學(xué)出版社. 1999.

[5] 文圣常. 海浪理論與計(jì)算原理[M]. 北京: 科學(xué)出版社. 1984.

[6] 徐德倫, 于定勇. 隨機(jī)海浪理論[M]. 北京: 高等教育出版社. 2001.

[7] 張鳳均. Labwindows/CVI開(kāi)發(fā)入門和進(jìn)階[M]. 北京:北京航空航天大學(xué)出版社. 2002.

[8] 楊惠珍, 康鳳舉, 褚彥軍, 等. 基于海浪譜的隨機(jī)海浪仿真及驗(yàn)證[J]. 系統(tǒng)仿真學(xué)報(bào), 2005, 17(10): 2324-2326.

2D Wave Simulation and Validation Based on Ocean Wave Spectrums

YU Long-fei1, ZHANG Ya-qun2
(1. Guangzhou Maritime Institute, Maritime, Guangzhou 510725, China； 2. CAS Key Laboratory of Renewable Energy, Guangzhou Institute of Energy Conversion, Chinese Academy of Sciences, Guangzhou 510640, China)

Wave simulation is a hot topic of Ocean Engineering. This article summarizes a method that can adapt to a wide range of random wave spectrums to emulate the waves environmentally. The method is converted to computer programs which are simple and convenient by language C and compute a large number of data used as data source to wave maker. The simulation result is verified by the experiment data. The test shows that the error between the analog spectrum and the target spectrum is lowest by using the numerical method presented in this paper. The acquisition of analog wave height data are used to estimate the spectrum by correlation function method. The numerical simulation method has a practical application value in the wave simulation.

wave spectrum； random number； target spectrum； evaluation of frequency spectrum

TK79；TP391.9

A

10.3969/j.issn.2095-560X.2015.03.009

2095-560X（2015）03-0214-04

于龍飛（1981-），男，博士，講師，主要從事海洋工程、船舶及新能源方面的研究。

2015-04-21

2015-05-08

廣州航海學(xué)院強(qiáng)校項(xiàng)目；中國(guó)科學(xué)院重點(diǎn)實(shí)驗(yàn)室基金（y407j31001）

? 通信作者：于龍飛，E-mail：yulongfei118@163.com