0 引　言

在海洋资源的开发过程中，畸形波等极端波浪时有发生，其严重威胁海洋装备的安全性。目前有关海洋畸形波的实测数据非常稀缺，基于聚焦波模拟极端波浪已成为模拟极端波浪的重要手段^[1]。

为了深入理解聚焦波的演化特性及更好获得极端波浪的模拟结果，国内外相关学者对于聚焦波的数值模拟方法展开了相关研究。Osborne等^[2]发现，在特定的初始边界条件下，波浪的非线性聚焦，可导致波高异常大的聚焦波。Westphalen等^[3]通过STAR-CCM+软件，建立数值波浪水池，利用“新波”理论，模拟获得所预期的瞬态聚焦波。对于上述2种聚焦波的生成方法，前者利用的是波能的非线性聚焦，而后者运用的是波能的线性叠加。由于波能的非线性聚焦，只有在特定条件下才能发生，因而前者很难控制聚焦波的产生结果。裴玉国等^[4]基于Longuet-Higgins模型通过调整组成波分量的初相位叠加生成聚焦波，研究发现组成波分量初相位的变化范围越小，所产生聚焦波的畸形程度越大。基于相同的生成模型，崔成等^[5]通过傅里叶分解“新年波”的实测波面获得波浪分量，重现模拟“新年波”的产生过程，最大波面处的模拟结果与实测波面符合很好，这表明聚焦波可以用来模拟海洋中的畸形波。

在波浪的数值模拟过程中，减小反射波对于模拟结果的影响，是获得高质量波浪计算结果的关键，国内外学者对于数值水池中的消波方法进行了研究。Orlanski^[6]提出一种Sommerfeld线性辐射条件，将其作为数值水池尾端的边界条件，比较结果表明该方法对于抑制反射波具有很好效果。然而Israeli等^[7]对该辐射边界条件进行了考察，结果表明当波浪的入射方向与边界的法向不同时，该方法对于反射波的抑制效果会变差。其后相关学者在波浪的模拟中提出了消波区方法，如Kim等^[8]通过在流体的动量方程中施加额外源项，强迫CFD计算结果在某一段流体计算域内趋于指定的结果；Choi等^[9]在流场中设置消波区，在消波区内将流体的能量被额外添加的阻尼项吸收，从而减少波浪的反射。苏绍娟等^[10]基于STAR-CCM+探究了阻尼消波和力消波2种不同消波方法，对规则波模拟结果的影响，比较结果表明力消波方法更适合长时间对于规则波的模拟；廖芸蓉等^[11]基于DualSPHysics平台结合阻尼消波法建立了一种中尺度波浪数值水槽，通过对比不同工况下的平均波高与水质点速度，获得了阻尼消波区的最佳长度；韩朋等^[12]通过在数值水池的消波区设定具有不同线性阻尼系数的海绵层，分析获到了不同随机波工况下阻尼系数的最佳取值范围。此外，有一些学者提出在数值水池中建立多孔结构进行消波，徐瑶瑶等^[13]基于多孔材料消波，研究了透水潜堤对聚焦波的影响，并分析了波高、堤顶水深等多种参数对消波效果的影响，结果证实相较于不透水潜堤，透水浅堤在消波方面更具优越性。

由上述研究结果可得，目前国内外学者对聚焦波数值模拟的研究主要集中于聚焦波的形成机理和模拟方法，且对于数值水池中消波的研究，主要针对规则波工况，对于聚焦波数值模拟中的消波问题却鲜有探讨。本文基于CFD方法建立数值水池，采用波能聚焦原理生成聚焦波，然后与物理水池中聚焦波的测量结果进行比较，验证本文数值方法的可行性。在此基础上，考察2种不同消波方法对聚焦波生成中反射波的消波效果，为数值水池中模拟聚焦波，选取消波方法和设置相应消波参数提供参考。

1 数值模拟理论 1.1 数值水池及消波方法

本文基于STAR-CCM+建立数值水池，采用有限体积法求解雷诺平均的N-S方程，选用Realizable k-ε湍流模型使方程封闭，通过VOF法捕捉自由表面，流场的相关控制方程为^[14]：

$ \rho \frac{{\partial \overline {{{{u}}_{{i}}}} }}{{\partial {{t}}}} + \rho \frac{{\partial \left( {\overline {{{{u}}_{{i}}}}\; \overline {{{{u}}_{{j}}}} } \right)}}{{\partial {{{x}}_{{j}}}}} = \rho {{g}} - \frac{{\partial {{p}}}}{{\partial {{{x}}_{{j}}}}} + \frac{\partial }{{\partial {{{x}}_{{j}}}}}\left( {\mu \frac{{\partial \overline {{{{u}}_{{i}}}} }}{{\partial {{{x}}_{{j}}}}} - \rho \overline {{{u}}{{'}_{{i}}}{{u}}{{'}_{{j}}}} } \right)，$

(1)

$ \frac{{\partial {{u}}}}{{\partial {{x}}}} + \frac{{\partial {{v}}}}{{\partial {{y}}}} + \frac{{\partial {{w}}}}{{\partial {{z}}}} = 0，$

(2)

$ \frac{{\partial {{{a}}_{{q}}}}}{{\partial {{t}}}} + \mu \left( {\frac{{\partial {{{a}}_{{q}}}}}{{\partial {{x}}}}} \right) + {{v}}\left( {\frac{{\partial {{{a}}_{{q}}}}}{{\partial {{y}}}}} \right) + {{w}}\left( {\frac{{\partial {{{a}}_{{q}}}}}{{\partial {{z}}}}} \right) = 0 。$

(3)

式中：$\overline {{u}_{i}} $为x、y、z方向上的雷诺平均速度分量；ρ为流体密度；g为重力加速度；μ为常数表示动力粘度；u、v、w为流体速度矢量在x、y、z三个坐标轴上的分量；t为时间；a_q为第q项中流体在单元中的体积分数，在本模拟中只存在气液两相。在数值水池中进行聚焦波的模拟，拟分别采用力消波与阻尼消波2种消波方法。

1.1.1 力消波方法

力消波法^[8]是将流场分为2个空间区域，一部分为CFD计算域，另一部分为强迫消波区域，后者与前者存在部分重叠，在二者的重叠区内，在动量方程中施加额外源项，利用强迫参数γ，使前者CFD的计算结果逐渐退化为后者的指定波浪结果，该方法的消波效果主要受到消波区长度X_F和强迫参数γ的影响。强迫参数γ的表达式为：

$ \gamma = {\gamma _{{\text{max}}}}{\cos ^2}\frac{{{\text{π}} {{X}}}}{{2{{{X}}_{{F}}}}}。$

(4)

式中：X为重叠域内的空间x纵坐标；X_F为力消波区域设定的长度。

1.1.2 阻尼消波方法

阻尼消波^[9]为在消波区内，在流场的动量方程中施加阻尼源项$S_z^d $，从而达到消波的目的：

$ S_z^d = \rho ({f_1} + {f_2}\left| {{w}} \right|)\frac{{{e^k} - 1}}{{{e^1} - 1}}{{w}}，$

(5)

$ k = {\left( {\frac{{{\text{X}} - {X_{sd}}}}{{{X_{{\text{e}}d}} - {X_{sd}}}}} \right)^{{n_d}}}。$

(6)

式中：X_sd为流场中波浪阻尼作用的起点；X_ed为波浪阻尼区的结束点；f₁、f₂、n_d分别为阻尼消波法的参数；w为流场中水质点的垂直速度分量。

1.2 聚焦波浪模型

本文基于线性叠加原理，在数值水池中生成聚焦波，在数值水池中二维聚焦波的波面方程为：

$ \eta (x,t) = \sum\limits_{i = 1}^{{N}} {{{{a}}_i}} \cos \left( {{k_i}\left( {x - {x_{{f}}}} \right) - 2{\text{π}} {f_i}\left( {t - {t_{{f}}}} \right)} \right)。$

(7)

式中：N为组成波个数；a_i为波浪分量的波幅；f_i为波分量的频率；k_i为波数；x_f为聚焦位置；t_f为聚焦时刻。

本文数值水池如图1所示，水池左侧设置为速度入口，右端设置为压力出口且在水池右端设置消波区域。为探究不同消波方法中的相应参数对聚焦波的消波影响及满足后续数据分析的需要，在水池中设置监测点A，其所测量的时历曲线需满足，最高频率波分量反射返回监测点时，最低频率波浪分量在入口边界处，发生二次反射未到达监测点A，即满足下式：

$ \left(L-L_1\right)-\left(\frac{L_1}{C_{g2}}-\frac{L_1}{C_{g1}}\right)C_{g1}+L+L_1 > C_{g2}\cdot \left(T_L+T_1\right)，$

(8)

$ {{L}} - {{{L}}_{{1}}} = \frac{{{{{T}}_{{L}}} \cdot {{{C}}_{{{g}}2}}}}{2}。$

(9)

式中：T₁和C_g1分别为最低频率波浪分量的周期和波能速度；C_g2为最高频率波浪分量的波能速度；T_L为数据分析的有效时长。

2 计算模型建立与校核

本文数值水池的计算网格如图2所示，为了减小数值计算的时间及准确地捕捉波浪自由表面的变化，在基础网格尺寸选定后，在水线面附近采用了网格加密。

2.1 收敛性分析

为了便于与物理水池中的测量结果进行比较，本文所研究的聚焦波工况，参考了物理水池中可能的聚焦波参数。所选取的聚焦波工况参数为：A_f=0.1 m，T_f=2.2～2.9 s，X_f=30 m，t_f=20 s，a_i=A_f/N，N=8。对于所建立的数值水池，左侧为造波边界及波浪产生域，右侧为消波区域，水池长度L=94.6 m，L₁=38 m，H=5.5 m，d=3.5 m。为了减小反射波对收敛性分析结果的影响，采用力消波法，消波长度取为最长波分量波长的2倍，力消波参数设置为10，根据后文力消波法的比较结果可得，基本没有波浪发生反射。此外，对于后文消波方法的研究，本文设置了对照组无反射工况，水池长度L=250 m，其他参数不变。

1）空间网格大小

流场的垂向基础网格大小ΔZ，分别取聚焦波幅A_f的1/10、1/20和1/40，波浪传播方向的网格尺寸ΔX设定为4倍ΔZ，不同ΔZ下计算所得的聚焦波时历曲线比较见图3。

A_f/20与A_f/40网格下计算所得的聚焦波在波峰处的不同大约为0.37%，这表明当ΔZ<A_f/20时，计算所得的聚焦波结果已经收敛，因而对于后续聚焦波工况的计算，选取A_f/20作为数值水池基础网格的尺寸。

2）时间步长

以聚焦波中最高频波浪分量对应的周期T₁为参考，时间步长Δt分别取T₁/500、T₁/1000和T₁/2000对聚焦波进行模拟。从图4中不同时间步长的聚焦波的比较结果可得，Δt=T₁/1000与T₁/2000二者聚焦波峰的差别约为0.3%，这表明当时间步长小于T₁/1000时，聚焦波的数值计算结果已经收敛，在本文的后续聚焦波的计算中取Δt= T₁/1000。

2.2 聚焦波生成数值方法验证

为了校核本文线性叠加生成聚焦波的准确性，对文献[15]中的物理水池聚焦波实验工况进行模拟，聚焦波试验工况为A_f=0.1 m，T_f=2.2～5.9 s，X_f=30 m，t_f=36 s，a_i=A_f/N，N=31。在聚焦波发生位置处（x=23.505 m），本文计算的波面与试验测量的时历曲线比较于图5中，总体上二者具有相类似的波面变化且在最大波峰处较符合较好，二者波峰处的误差约为2%。由此可得本文所建立的数值水池对于聚焦波的模拟结果具有很好的可靠性，可用于本文后续对于聚焦波模拟中消波方法的考察研究。

3 聚焦波生成的不同消波方法研究

本节的聚焦波工况，与收敛性分析的工况相同，参考了物理水池中的实验工况参数，而非海洋中真实的聚焦波浪，聚焦波工况的的具体参数为：A_f=0.1 m，T_f=2.2～2.9 s，X_f=30 m，t_f=20 s，a_i=A_f/N，N=8。

3.1 力消波方法对聚焦波生成的消波研究 3.1.1 力消波方法的X_F考察

保持数值波浪水池的总长度不变，改变消波区长度X_F，分析力消波法对聚焦波中波浪分量的消波作用。以聚焦波中最高频率波浪分量对应的波长λ_N为参考，取X_F的变化范围为0～4λ_N。对指定监测点A处，所测得的波面时历，进行快速傅里叶变换（FFT），获得聚焦波的波浪分量的幅值谱，对不同X_F相应的波浪分量幅值进行比较，如图6所示。

不同力消波区域长度X_F下相应频率波浪分量的幅值比较结果表明，消波区的长度对于反射波的消波效果具有显著影响，消波区如果设定的比较小，将使所生成聚焦波的低频波浪分量的幅值产生波动，随着消波区长度的变大，各波浪分量的幅值趋于稳定，且与相应不发生反射工况波分量的幅值大体相当，这表明当消波区长度超过一定值（如X_F=λ_N）时，反射波可被很好抑制。

为了量化不同消波区长度对于聚焦波中波浪分量的消波效果，将特定X_F下监测的各波分量幅值$a_i^1 $与相应无反射波分量幅值$a_i^0 $进行比较，计算它们相对误差的绝对值，然后求和得Er如式(10)所示，不同消波区长度X_F下的消波误差Er比较见图7。

$ Er = \sum\limits_{i = 1}^n {\frac{{\left| {a_i^1 - a_i^0} \right|}}{{a_i^0}} \times 100}{\text{%}} $

(10)

误差比较结果显示，当λ_N≤X_F≤3λ_N，监测所得的波分量的误差小于5%，且趋于稳定；当X_F>3λ_N，Er显示出了一定的增大，因而对于聚焦波的生成，在采用力消波法消除反射波时，建议设定的消波区域长度为λ_N≤X_F≤3λ_N，可很好抑制反射波对于聚焦波生成结果的影响。

3.1.2 力消波方法的消波参数γ考察

基于上文对于消波区长度的分析结果，指定消波区域长度为2λ_N，为了研究消波参数γ对于聚焦波工况消波效果的影响，选取参数γ的变化范围为2.5～60。

图8中不同γ参数下，相应频率波浪分量的比较结果表明，γ参数的取值过大（如γ=60），与无反射波浪分量幅值相比，可导致波浪分量的幅值不稳定，可能会明显增大（如ω₁对应波分量）或者显著降低（如ω₈对应波浪分量）波浪分量的幅值。出现上述结果，可能是由于大的消波系数γ，导致消波区内强迫作用增强，从而影响了所观测计算域内波浪的求解。

图9中不同参数γ的消波误差比较结果表明，γ参数越大，对于数值水池中反射波的消波效果越弱，当γ小于10时，力消波法的消波效果趋于稳定，消波误差小于4.2%，因而在数值水池中，采用力消波法消除反射波时，建议γ系数的取值小于10。

3.2 阻尼消波方法对聚焦波生成的消波研究 3.2.1 阻尼消波方法的X_D考察

保持数值水池总长度不变，取阻尼区的长度X_D=0.25λ_N～4λ_N，图10中相应频率波浪分量幅值的比较结果表明，如果阻尼区的长度取值较小（如X_D=0.25λ_N），可使监测点处波浪分量的幅值不稳定，可能会显著高于其无反射波分量的幅值，或明显偏小。随着阻尼区长度的增大，各频率波浪分量幅值的波动变弱，当阻尼区长度大于λ_N，监测点处所得的波浪分量的幅值趋于稳定且与无反射工况的幅值相当，这表明阻尼区对于反射波浪的抑制效果已达到最优。

图11中不同阻尼区长度下的消波误差比较结果显示，随着阻尼区长度的增大，数值水池中由于反射波产生的误差逐渐降低，当X_D≥2λ_N时，水池中所生成的聚焦波的波分量误差趋于稳定（约为6%）。因而，基于阻尼法消波时，阻尼区长度不应小于2λ_N，可以很好抑制反射波对于聚焦波生成结果的影响。

3.2.2 阻尼消波方法的消波参数f₁、f₂、n_d考察

在阻尼消波方法中，有3个阻尼消波参数f₁、f₂、n_d，为探究上述3个参数对消波结果的影响，采用单一变量法逐一进行分析。

1）消波参数f₁对消波效果的影响

取f₁的变化范围为2.5～40，f₂和n_d分别取默认参数10和2。图12中不同f₁下的波分量比较结果表明，f₁的取值偏小（如f₁=2.5）可使监测点处波浪分量的幅值不稳定，由阻尼消波法的式(5)可得f₁越小，阻尼消波模型中阻力项$S_z^d $越小，从而导致消波效果变差，导致波浪反射严重。当f₁≥10，计算所得的各频率波分量幅值基本趋于稳定，表明水池中的反射波得到了很好抑制。图13中不同f₁的误差比较结果亦表明，随着参数f₁的增大，数值水池中所测得反射波的误差逐渐降低，当f₁=10时，聚焦波的波浪分量由于反射波而产生的误差不足6%。

2）消波参数f₂对消波效果的影响

选取阻尼消波参数f₂的变化范围为2.5～40，f₁和n_d分别设定为10和2。从图14中不同f₂下波浪分量的幅值及图15中误差的比较结果可以发现，相应频率下的各波浪分量幅值的基本相同，此外不同f₂下的误差结果亦无显著差别，这可能是由于阻尼参数f₂与流场中水质点的垂向速度有关（由式(5)可得），本文所计算的聚焦波波幅为0.1 m，垂向速度分量相对较小，因而f₂参数对于数值水池中反射波的结果影响有限，这也暗示相对于f₁，反射波对于参数f₂的依赖性弱。

3）消波参数n_d对消波效果的影响

选取n_d的变化范围为1.0～3.0，f₁和f₂均分别设置为10。图16中不同n_d下波浪分量的幅值比较结果表明，当n_d=1.0时，与无反射波浪分量相比，ω₁对应的波浪分量的幅值明显偏大，然而当ω≥ω₆，各波浪分量的幅值低于其相应参照工况的幅值；对于n_d=3.0的波分量幅值比较结果与n_d=1.0的分析结果相类似。此外，对于其他n_d值，相应频率波分量的幅值基本与其参照工况的波幅大体相当。

图17中不同消波参数n_d的误差比较结果显示，当n_d=1时数值水池中指定点所测得的聚焦波分量的误差较大，随着n_d的增大，误差逐渐减小，然而当n_d>2.5时，误差指标Er亦开始表现出增大变化趋势，这表明n_d在1.0～3.0之间存在最佳值。根据本文的计算结果，在数值水池中进行聚焦波模拟时，建议n_d取值为2.0～2.5。

此外，还对聚焦波工况（A_f=0.1 m，T_f=2.2～5.9 s，X_f=30 m，t_f=36 s，a_i=A_f/N，N=31）进行了消波考察，所的研究结果与上述分析结果基本一致，这进一步验证本文分析所得的参数建议范围，可为数值水池中相应消波方法参数的设定提供参考。

4 结　语

本文基于STAR－CCM+软件建立数值水池，针对聚焦波模拟中的２种消波方法进行考察，以期获得聚焦波模拟时相应消波方法的最佳参数设置，根据本文聚焦波工况的计算结果和所定义的误差指标Er，获得的主要结论如下：

1）对于力消波法，数值水池中监测点的聚焦波的分析结果表明，消波区长度X_F对于反射波的消波效果具有显著影响，消波区的长度如果设定的偏小，将使所生成聚焦波的低频波浪分量的幅值产生波动，随着消波区长度的变大，各波浪分量的幅值趋于稳定。其γ参数的取值过大，可导致波浪分量的幅值不稳定，可明显增大或者显著降低波浪分量的幅值。由误差评价指标Er的比较结果可得，力消波区长度设定为λ_N≤X_F≤3λ_N，力强迫系数选取为γ≤10，对于数值水池中聚焦波的模拟，可以很好抑制反射波（Er≈4%）。

2）对于阻尼消波法，数值水池中所监测聚焦波的比较结果显示，阻尼区的长度取值偏小，可使监测点处波浪分量的幅值不稳定，可显著高于或明显降于其无反射波分量的幅值。随着阻尼区长度的增大，各波浪分量幅值的波动变弱。较小的f₁的取值，使监测点处波浪分量的幅值不稳定，随着参数f₁的增大，数值水池中所测得的反射波的误差逐渐降低。本文所定义的反射波的误差评价指标表明，对于阻尼消波法，其阻尼区长度X_D≥2λ_N，阻尼参数f₁≥10，n_d取值为2.0～2.5，可有效抑制反射波对于聚焦波生成结果的影响。此外，对于参数f₂，该参数对于反射波的影响不显著，在数值水池中进行聚焦波模拟时，f₂可取默认值10。