多向不规则波群传播的数值模拟
刘思, 张永良     
清华大学 水沙科学与水利水电工程国家重点实验室, 北京 100084
摘要:波群的模拟是与波群有关的数值和物理模型试验的前提。实际的波浪是多向的不规则波, 为了在实验水池中模拟多向不规则波群, 该文首先根据给定的波要素、群性参数和方向分布参数, 给出了多向不规则波群的数值模拟方法, 在此基础上利用有限元法求解改进的Boussinesq方程的数值计算模型, 建立了数值水池在指定位置模拟所要求群性的多向不规则波的入射波浪边界条件的计算方法。典型算例表明: 按照该文建议的方法确定入射波浪边界条件而建立的数值模型, 可以在水池指定位置处产生满足给定群性要求的多向不规则波浪。此外, 对多向不规则波群的传播进行了数值研究。
关键词多向波    波群    群性参数    数值模拟    
Numerical simulation of multidirectional irregular wave group propagation
LIU Si, ZHANG Yongliang     
State Key Laboratory of Hydroscience and Engineering, Tsinghua University, Beijing 100084, China
Abstract: Irregular wave group motion was simulated for comparison with physical model tests. Real waves are multidirectional waves. A numerical method was developed to simulate multidirectional irregular wave groups for various wave parameters, grouping factors and directional spreading parameters for use in designing an experimental wave basin for multidirectional irregular wave groups. The simulations of the wave time series of the incident waves was based on the numerical model of modified Boussinesq equations solved by finite element. Examples show that multidirectional waves containing the desired wave groupiness can be generated at a specified position in the numerical wave basin. Transformations of the wave groupings in the basin are also numerically analyzed.
Key words: multidirectional waves    wave groups    grouping factor    numerical simulation    

海浪的群性是海浪场的一种重要特性,它对海堤[1]、 港湾共振[2, 3]、 系泊浮体运动[4, 5, 6]等有明显影响,因此,在海洋工程设计中需要进行波群对建筑物作用的试验,这就要求模拟波群。

关于如何模拟波群各国学者提出了很多不同的方法。Rye[7]认为仅需严格模拟频谱谱形来模拟波群; Funke和Mansard[8]认为不仅要模拟频谱,同时要模拟相位谱,并建议通过SIWEH(smooth instaneous wave energy history)来模拟; 徐德伦等[9]建议把波包谱和波浪频谱结合起来模拟波群。俞聿修和桂满海[10]基于由实测资料的分析而给出的波包谱经验公式并结合波浪频谱,模拟出了不同群性的不规则波列。柳淑学等[11]在俞聿修和桂满海[10]基础上,引入方向分布函数,给出了多向波群波浪波列的数值模拟方法。然而,目前关于不规则波波群模拟的研究一般基于给定波浪及其波群参数,探究数值模拟方法产生相应的波列,而对于不规则波波群在数值波浪水池或物理实验水池中的传播模拟方法的研究还很少,其问题的根本在于,欲在水池中产生所需要的波群波浪,如何确定相应的数值水池入射边界条件或物理实验水池的造波边界条件,尤其对于多向不规则波波群波浪的模拟,现仍未见报道。

考虑俞聿修和桂满海[10]依据的波包谱实测资料较少,刘思等[12]在收集了更多的实测波浪资料基础上,分析得到一个可反映天然海浪的波包能量分布的波包谱经验公式。本文基于该海浪波包谱[12],参考柳淑学等[11]关于模拟多向不规则波群的研究成果,首先给出多向不规则波群的数值模拟方法,在此基础上,利用有限元法求解改进的Boussinesq方程的数值计算模型,提出可在数值水池中产生所要求的多向不规则波波群波浪的入射边界条件的计算方法,通过不同计算示例验证所提方法的有效性,同时对多向不规则波群的传播特性进行数值模拟和分析。

1 多向不规则波群波列的数值模拟

参考柳淑学等[11]模拟多向不规则波群的方法,考虑到波浪的多向性,引进方向分布函数,把方向谱和波包-方向谱结合起来模拟多向波群,由波包面A(x,y,t)和相位函数φ(x,y,t),可构造出三维波面序列η(x,y,t),即

$\eta \left( {x,y,t} \right) = A\left( {x,y,t} \right)\cos \left[{\varphi \left( {x,y,t} \right)} \right]$ (1)

模拟得到的多向不规则波η(x,y,t)的频谱和波包谱应与靶谱相符,其方向分布应与输入的方向分布相似。波包面A(x,y,t)和相位函数φ(x,y,t)的计算可参见柳淑学等[5]给出的方法。对于多向波,方向分布函数取光易型分布函数[13],即

$\begin{array}{l} G\left( {\omega ,\theta } \right) = {G_0}\left( s \right){\cos ^{2s}}\frac{{\theta - {\theta _0}}}{2},\\ {G_0}\left( s \right) = \left[{\int_{{\theta _{\min }}}^{{\theta _{\max }}} {{{\cos }^{2s}}\frac{{\theta - {\theta _0}}}{2}{\rm{d}}\theta } } \right] \end{array}$ (2)

其中: 参数s为方向分布函数集中度参数; [θmin,θmax]为波浪的方向分布范围。

考虑波群具有群高和群长这2个方面的特征,波包靶谱采用文[12]建议的经验公式,即

$\begin{array}{l} \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{S_A}\left( f \right) = \\ \left[{0.042 + 0.019\left( {f/{f_{PA}}} \right)} \right]\pi {m_0}{\left( {GFH} \right)^2}/{f_{PA}},\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;0 \le f \le {f_{PA}};\\ {S_A}\left( f \right) = 0.085\pi {e^{ - \frac{1}{{3.1}}\left( {f/{f_{PA}}} \right)}}{m_0}{\left( {GFH} \right)^2}/{f_{PA}},\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;f > {f_{PA}} \end{array}$ (3)

其中:

$\begin{array}{l} 当GFH \le 0.7,f/{f_{PA}} = 5 \sim 15;\\ 当GFH > 0.7,f/{f_{PA}} = 1028 \end{array}$ (4)

其中: m0为波谱的零阶矩; fPfPA分别为波谱的峰频与波包谱峰频; $GFH = \sqrt 2 {\sigma _A}/\bar A$和$GLF = {f_P}/{f_{PA}}$分别为表征波浪群性高度和长度的参数,称为群高和群长参数[12],式中σA和$\bar A$分别为波包线的方差和均值。可见,式(3)和式(4)涵括群高与群长这2方面的因素,模拟中可通过给定不同的波包靶谱参量GFH和GLF,来控制波浪的群高和群长特征。

考虑多向不规则波浪方向谱分析时阵列的布置原则[13],采用图 1所示的五角形浪高仪阵列模拟波浪,波浪的方向谱分析采用Bayes方法[14]。为了验证本文方法的有效性,按上述方法分别模拟了群性参数(GFH=0.4、 0.5,GLF=5、 10、 15)和(GFH=0.6、 0.7、 0.8、 0.9,GLF=10、 15、 20、 25) 的多向不规则波浪,给定的波要素为有效波高 ${H_{\frac{1}{3}}} = 0.07m$ ,有效周期${T_{{H_{\frac{1}{3}}}}} = 1.5s$,水深h= 1.0 m,采用合田改进的JONSWAP谱[15]作为目标频谱,取谱峰升高因子为3.3,方向集中度参数s分别取为10、 30、 70,方向分布范围为[-ππ]。

图 1 浪高仪布置形式图

作为示例,图 2给出了群高参数GFH=0.9、 群长参数GLF=15、 方向集中度s=30时模拟得到的1号测点处的波面序列及波包线序列,及其计算频谱、波包谱和对应峰频fp及相邻频率处的方向分布与给定谱及方向分布的比较。从图 2可以看出,波浪的频谱、波包谱和方向分布满足模拟要求。本文方法可以有效地数值模拟同时满足给定群高和群长要求的多向不规则波群波浪序列。

图 2 s=30、 GFH=0.9、 GLF=15 时数值模拟波浪分析结果
2 多向不规则波群传播的数值模拟计算和分析

孙忠滨等[16]建立了非结构化网格有限元法求解改进的Boussinesq方程的数值计算模型,能较好的模拟多向不规则波浪的传播。本文利用此数值模型,给出可以模拟多向波群波浪的入射波浪边界条件,从而可在水池指定位置处产生满足给定群性要求的多向不规则波浪,并对多向不规则波浪场中群性的变化情况进行分析。

2.1 入射波浪边界条件

为了利用数值模型在水池指定位置处产生满足给定群性要求的多向不规则波浪,需要基于前述多向不规则波群的数值模拟波浪序列结果η(x,y,t)建立适当的入射波浪边界条件,包括入射波浪波面及其对应的沿水深平均速度。

对于单向规则波,当入射波面η0给定后,其对应的沿水深平均速度(u0,v0)可通过线性化改进的Boussinesq方程求得,即

${u_0} = \frac{\omega }{{k\left( {h + {\eta _0}} \right)}}{\eta _0}\cos \theta $ (5)
${v_0} = \frac{\omega }{{k\left( {h + {\eta _0}} \right)}}{\eta _0}\sin \theta $ (6)

其中: θ为相对于x轴的波浪传播方向角; h为水深; ω为频率; k为波数。对于多向不规则波,入射波浪边界条件可由线性叠加法来确定。首先将多向不规则波群的数值模拟结果η(x,y,t)进行Fourier展开,得到对应各频率ωi组成波的振幅ai和相位εi。如果各组成波传播方向θi已知,依据线性波浪理论,对于每个组成波,不同位置处的波面可用式(7)表示,即

${\eta _i}\left( t \right) = {a_i}\cos \left( {{w_i}t + {k_i}x\cos \theta + {k_i}y\sin \theta + {\varepsilon _i}} \right)$ (7)

假定(x, y)为水池欲产生具有给定波群特性波浪位置的坐标,(x0,y0)为入射边界的坐标,由线性叠加法,根据式(5)、 式(6)和式(7),入射边界条件可由式(8)—式(10)确定:

${\eta _0}\left( {{x_0},{y_0},t} \right) = \sum\limits_{i = 1}^N {{a_i}\cos \left[{{\omega _i}t + {k_i}\left( {x - {x_0}} \right)\cos {\theta _i} + {k_i}\left( {y - {y_0}} \right)\sin {\theta _i} + {\varepsilon _i}} \right]} $ (8)
${u_0}\left( {{x_0},{y_0},t} \right) = \sum\limits_{i = 1}^N {\frac{{{\omega _i}{a_i}\cos \left[{{\omega _i}t + {k_i}\left( {x - {x_0}} \right)\cos {\theta _i} + {k_i}\left( {y - {y_0}} \right)\sin {\theta _i} + {\varepsilon _i}} \right]}}{{{k_i}\left( {h + {a_i}\cos \left[{{\omega _i}t + {k_i}\left( {x - {x_0}} \right)\cos {\theta _i} + {k_i}\left( {y - {y_0}} \right)\sin {\theta _i} + {\varepsilon _i}} \right]} \right)}}} \cos {\theta _i}$ (9)
${v_0}\left( {{x_0},{y_0},t} \right) = \sum\limits_{i = 1}^N {\frac{{{\omega _i}{a_i}\cos \left[{{\omega _i}t + {k_i}\left( {x - {x_0}} \right)\cos {\theta _i} + {k_i}\left( {y - {y_0}} \right)\sin {\theta _i} + {\varepsilon _i}} \right]}}{{{k_i}\left( {h + {a_i}\cos \left[{{\omega _i}t + {k_i}\left( {x - {x_0}} \right)\cos {\theta _i} + {k_i}\left( {y - {y_0}} \right)\sin {\theta _i} + {\varepsilon _i}} \right]} \right)}}} \sin {\theta _i}$ (10)

其中: η(x0,y0,t)为入射边界处入射波浪波面; u(x0,y0,t)和v(x0,y0,t)为其对应的沿水深平均水质点速度; ki为组成波的波数。

然而如何确定各组成波的方向θi是利用式(8)—式(10)计算入射边界波浪条件的关键。本文借鉴多向不规则波数值模拟中的单叠加法[13],关于方向分布方向的分割采用等分能量的方法: 对方向分布函数G(ω,θ)积分可得到其累积分布函数p(θ)为

$p\left( \theta \right) = \int_{{\theta _{\min }}}^\theta {G\left( {\omega ,\theta } \right)} {\rm{d}}\theta ,0 \le p\left( \theta \right) \le 1$ (11)

对于任一频率ωi,其方向角θi可随机选取,但其分布应满足G(ω,θ)。计算中,对于每一个ωi,随机地选取p(θi),由式(11)通过数值迭代即可确定θi,这样选取的θi比较集中在主方向的附近。

2.2 单向不规则波群传播的模拟

依据上述确定入射波浪边界条件的方法,对于二维问题,由该数值模型模拟得到的数值水槽内指定位置处的波浪应与不规则波群的数值模拟结果 η(t)相同。刘思等[17]曾在二维数值模拟波群的基础上,将浪高仪置于水槽距造波机x=5 m处物理模拟了不同群性的波浪,同时也测得了这些不同群高和不同群长的波浪其群性参数的沿程变化情况。为了验证上述确定入射边界条件方法的有效性,首先应用改进的Boussinesq方程的数值计算模型,对单向不规则波群的传播进行了数值模拟,模拟结果与实验结果进行了对比。

采用上述方法由数值模拟波群波浪序列的结果确定入射波浪边界条件,此时θi取为0,并代入指定位置x=5 m。对应于二维物理模型实验,在数值计算模型中采用波浪要素为有效波高${H_{\frac{1}{3}}} = 0.07{\rm{m}}$,有效周期${T_{{H_{\frac{1}{3}}}}} = 1.5s$,水深h=0.5 m,时间间隔0.01s,波浪历时327.68 s,计算区域采用三角形单元。在计算过程中,分别计算了4个断面上的群高参数GFH和群长参数GLF。图 3图 4分别为数值模拟分析结果与实验分析结果,图 3图 4中横轴x为测点距造波板或造波边界的距离。从图 4可以看出,在x=5 m处能够模拟产生指定群性的波浪,且不同群性的波浪随波浪的传播其群性差异逐渐减小。与图 3对比,考虑到采用线性的边界条件数值计算与物理模拟波浪之间的差别,总体上来讲,数值模拟结果与实验结果基本一致。由此可见,按本文方法给定入射波浪边界而建立的数值模型是有效的。

图 3 实验测得群高和群长沿水槽的变化
图 4 数值模拟群高和群长沿水槽的变化
2.3 多向不规则波群传播的模拟和分析

对于多向不规则波波群,θi的选取则更为关键,为了验证能否通过本文的选取方法,由该数值模型在水池指定位置处产生多向不规则波群的波浪,取计算区域如图 5所示。数值模拟的水池长 18 m、 宽45 m,入射边界相邻的两侧边界及右侧边界均为吸收波浪边界。为了在有效试验区内点P3(6 m,22.5 m)处产生同时满足给定的波要素、方向分布和群性要求的波浪,按照2.1所述方法确定各组成波的方向、振幅和相位,代入坐标x=6 m、 y=22.5 m,由式(8)—式(10)式确定入射边界条件。数值模拟中采用的波浪要素为有效波高${H_{\frac{1}{3}}} = 0.07{\rm{m}}$,有效周期${T_{{H_{\frac{1}{3}}}}} = 1.5s$,时间间隔0.01 s,波浪历时327.68 s,计算区域采用三角形单元。在计算过程中,分别在B1B3P1P6点记录波浪,并在其中B1B3P3P4点布置按照方向谱分析要求设置阵列记录波浪[18],采用Bayes法进行方向谱的分析。

图 5 数值计算区域

图 6给出方向分布集中度s=30时,GFH=0.9、 GLF=15数值水池模拟波群波浪和直接由数值模型计算的波面分析所得频谱、波包谱及峰频波浪的方向分布与其对应的目标值的比较; 同时,对由数值模型计算得到的P3测点处的波面与数值水池波浪传播模拟波面进行了比较。从目标谱和模拟分析所得结果的比较可以看出,数值水池波浪传播模拟波群波浪和数值模型模拟波面计算所得的频谱、波包谱以及方向分布与目标值基本一致,同时数值水

图 6 s=30, GFH=0.9、 GLF=15时模拟分析结果

池模拟波面与数值计算所得波面η(x,y,t)很相近,结果是令人满意的。由此可见,按照前述各组成波方向θi的选取方法确定入射波浪边界条件而建立的数值模型,可以有效地在水池指定位置 (x y) 处模拟得到同时满足波要素、方向谱和群性要求的多向不规则波浪。

图 7图 8分别为数值水池中模拟不同群性波浪传播时波浪场内不同测点处的频谱和波包谱的变化情况。从图 7a图 8a可以看出,P1P4点处频谱的形状基本相似,P1点处的频谱的谱密度最大,由P1P4,谱密度略呈减小趋势; 同时,从图 7c图 8c可以看出,P5P6P3点处频谱的形状基本不变。此外,从图 7b图 8b可以看出对于群性大的波浪(GFH=0.9、 GLF=15)、 P3波浪组次(GFH=0.5、 GLF=6),P3点处的波包谱比周围测点处的谱密度都小。

图 7 s=30、 GFH=0.9、 GLF=15时测点P1P6处的频谱和波包谱
图 8 s=30、 GFH=0.5、 GLF=6时测点P1P6处的频谱和波包谱

前面二维实验和数值模拟计算结果表明,指定位置处模拟得到的群性大的波浪随波浪传播群性逐渐减小,群性小的波浪随传播其群性逐渐增大,群高GFH渐渐都趋于0.7左右。已有分析结果表明,实测波浪资料的群高GFH也大多集中在0.7左右[12]。此外,仅模拟频谱的实验室波浪往往具有中等强度的群性[19]。与二维实验结果一致,本文模拟得到的波浪在达到指定点之前群性差别逐渐增大,在指定点处表现出明显不同的群性差异,满足模拟要求。随波浪继续传播,其群性差异逐渐减小,趋于中等强度。因此数值模拟得到的、在水池指定位置处群性大的波浪相比周围测点处波浪的群性都大; 同样的频谱条件下,表现在波浪的波包谱比周围测点处的谱密度都大; 同理,在指定点处模拟得到的群性小的波浪相比周围测点处波浪的群性都小,表现在波浪的波包谱比周围测点处的谱密度都小。

图 9图 10分别为方向分布集中度s=30时,数值水池内模拟的7组波浪分析所得群高和群长参数及有效波高沿水池纵向和横向的变化情况。正如所预期的,一方面,群高参数和群长参数在水池指定位置(6 m,22.5 m)处基本满足模拟要求,随波浪向前继续传播群性逐渐趋于一致。另一方面,模拟得到的各组次波浪的有效波高值随波浪传播略呈减小趋势,且在水池指定点(6 m,22.5 m)处波浪的有效波高值随群高的增大而增大。这是由于同样的频谱条件下,各组次波浪的波面平均能量基本相同,而群高大的波浪其波包线的波动幅度较大,正是这种波浪能量分布的变化使得其有效波高值相对较大,这与文[11]给出的结果是一致的。

图 9 数值模拟波浪群高和群长及有效波高沿水池纵向的变化
图 10 数值模拟波浪群高和群长及有效波高沿水池横向的变化
3 结 论

本文首先通过模拟实例表明,以方向谱和波包-方向谱结合起来作为靶谱可以有效地数值模拟多向不规则波群。然后利用有限元法求解改进的Boussinesq方程的数值计算模型,建立了数值水池在指定位置模拟所要求群性的多向不规则波的入射波浪边界条件的计算方法,分别对单向和多向不规则波群的传播进行了数值模拟和分析。典型算例表明,按照本文建议的方法确定入射波浪边界条件而建立的数值模型,可以在水池指定位置处产生满足给定群性要求的多向不规则波浪。 此外,数值模型计算给出了关于频谱、波包谱、方向分布及群性参数和波浪统计特征的分析结果。本文的研究提供了实验水池模拟多向不规则波群的造波信号的计算方法,为在实验水池模拟多向不规则波群奠定了理论基础。

参考文献
[1] Johnson R R, Mansard E P D, Ploeg J. Effects of wave grouping on breakwater stability [C]// Proceedings of the 16th International Conference on Coastal Engineering. Hamburg, Germany: The American Society of Civil Engineers, 1978: 2228-2243.
[2] Bowers E C. Harbour resonance due to set down beneath wave groups [J]. Journal of Fluid Mechanics, 1977, 79: 71-92.
[3] Wu J K, Liu P L F. Harbour excitations by incident wave groups [J]. Journal of Fluid Mechanics, 1990, 217: 595-613.
[4] Mansard E P D, Pratte B D. Moored ship response in irregular waves [C]// Proceedings of the 17th International Conference on Coastal Engineering. New York, USA: The American Society of Civil Engineers, 1982.
[5] Balaji R, Sannasiraj S A, Sundar V. Detection of wave groups from the motion behaviour of a discus buoy [J]. Journal of Hydro-Environoment Research, 2008, 1(3): 195-205.
[6] MA Xiaojian, SUN Zhaochen, ZHANG Zhiming, et al. The effect of wave groupiness on a moored ship studied by numerical simulations [J]. Journal of Hydrodynamics, 2011, 23(2): 145-153.
[7] Rye H. Ocean Wave Groups, UR-82-18 [R]. Trondheim, Norway: Department of Marine Technology, University of Trondheim, 1982.
[8] Funke E R, Mansard E P D. On the synthesis of realistic sea states in laboratory flume [C]// Proceedings of the 17th International Conference on Coastal Engineering. New York, USA: The American Society of Civil Engineers, 1980: 2974-2991.
[9] 徐德伦, 候伟. 随机波群的模拟 [J]. 海洋学报, 1993, 15(5): 27-35.XU Delun, HOU Wei. Simulation of random wave groups [J]. Acta Oceanologica Sinica, 1993, 15(5): 27-35. (in Chinese)
[10] 俞聿修, 桂满海. 波群的数值模拟和物理模拟 [J]. 大连理工大学学报, 1998, 38(1): 86-91.YU Yuxiu, GUI Manhai. Numerical simulation and physical simulation of sea wave groups [J]. Journal of Dalian University of Technology, 1998, 38(1): 86-91. (in Chinese)
[11] 柳淑学, 包艳, 俞聿修. 多向不规则波群的数值模拟 [J]. 水道港口, 2006, 27(5): 273-278.LIU Shuxue, BAO Yan, YU Yuxiu. Numerical simulation of multidirectional irregular wave groups [J]. Journal of Waterway and Harbour, 2006, 27(5): 273-278. (in Chinese)
[12] 刘思, 柳淑学, 俞聿修. 改进的海浪波包谱 [J]. 水道港口, 2010, 31(4): 229-235.LIU Si, LIU Shuxue, YU Yuxiu. Improved wave envelope spectrum [J]. Journal of Waterway and Harbour, 2010, 31(4): 229-235. (in Chinese)
[13] 俞聿修. 随机波浪及其工程应用 [M].大连: 大连理工大学出版社, 2003.YU Yuxiu. Random Wave and Its Applications to Engineering [M]. Dalian: Dalian University of Technology Press, 2003. (in Chinese)
[14] Goda Y A. A comparative review on the functional forms of directional wave spectrum [J]. Coastal Engineering, 1999, 41(1): 1-20.
[15] Goda Y A. Random Seas and Design of Maritime Structures [M]. Tokyo, Japan: World Scientific Publishing Company, 2000.
[16] 孙忠滨, 柳淑学, 李金宣. 基于改进Boussinesq方程的三角形网格有限元模型 [J]. 海洋工程, 2010, 28: 50-58. SUN Zhongbin, LIU Shuxue, LI Jinxuan. Triangular element FEM model based on modified Boussinesq equations [J]. Ocean Engineering, 2010,28: 50-58. (in Chinese)
[17] 刘思, 柳淑学, 李金宣, 等. 单向不规则波群的实验室模拟和分析 [J]. 水道港口, 2011, 32(5): 305-312.LIU Si, LIU Shuxue, LI Jinxuan, et al. Experimental simulation and analysis of irregular sea wave groups [J]. Journal of Waterway and Harbour, 2011, 32(5): 305-312. (in Chinese)
[18] YU Yuxiu, LIU Shuxue. Study on measurement methods of multi-directional waves by gauge array method [J]. China Ocean Engineering, 1992, 6(2): 233-242.
[19] 杨正已, 左其华, 苏卫东. 波群的实验室模拟和分析 [J]. 海洋工程, 1991, 9(1): 61-67.YANG Zhengyi, ZUN Qihua, SU Weidong. Experimental simulation and analysis of wave groups [J]. Ocean Engineering, 1991, 9(1): 61-67. (in Chinese)