波浪在局部可渗透水平海床上传播的解析解

倪云林，滕斌

(1.大连理工大学 海岸和近海工程国家重点实验室，辽宁 大连 116024；2.浙江海洋大学 海洋工程装备学院，浙江 舟山 316022)

1 引言

多孔介质海床具有渗透性，当波浪在多孔介质海床上传播时会发生波高衰减和波能损失。但是大多数的水波理论假设海底是不可渗透的，这势必会对波浪长距离传播的推算造成误差，对海岸及近海工程中设计波浪要素的计算带来影响。

关于波浪在渗透海床上传播变形的问题，许多学者采用不同的方法开展了相关研究。其中，在实验研究方面，Savage和Fairchild[1]开展了波浪在粗砂和细砂海床上传播能量衰减的实验；Özhan和Shi-Igai[2]、Gu和 Wang[3]、Sawaragi和 Deguchi[4]先后开展了不同类型的波浪在多孔介质海床上波高衰减的实验研究。在理论研究方面，Putnam[5]基于线性波浪理论和达西定律，推导了因波浪运动产生渗流而导致波能损失的理论；Reid和 Kajiura[6]、Murray[7]、Dean和 Dalrymple[8]同样基于线性波浪理论，深入研究了刚性、无限深多孔介质海床上波浪衰减的边界值问题，并得到了复波数色散关系。在数值模型研究方面，主要是利用缓坡类方程[9–11]、Boussinesq 型方程[12]、Navier-Stokes方程[13–14]研究波浪与渗透潜堤的相互作用。

现有的理论研究主要针对波浪在全渗透海床上的传播问题，且大多数只考虑了传播模态。本文推导波浪在局部可渗透海床上传播，包含非传播模态的解析解，探讨渗透系数、相对水深和渗透海床长度对波高衰减的影响以及局部渗透海床对波浪的反射和透射。

2 数学公式

2.1 理论推导

考虑二维笛卡尔坐标系统下波浪在局部可渗透水平海床上的传播问题。如图1所示，整个流域水深为常数h，海床水平且局部可渗透，其中，区域II(−L/2≤x≤L/2)海底可渗透，长度为L，其下区域IV为不可压缩饱和多孔介质组成的无限深可渗透海床，其左侧区域I(x<−L/2)和右侧区域III(x>L/2)为海底不透水的半无限长区域。选取坐标系原点为渗透海床垂直平分线与静水位的交点，x轴水平向右为正，z轴竖直向上为正。

图1 波浪在局部可渗透海床上传播示意图Fig.1 Definition sketch of wave propagation over a local porous seabed

假设流体为无黏、不可压缩的均质液体，波浪运动无旋且波面升高和运动速度足够小，则可以应用线性波浪理论。对于不可渗透海床上波浪沿x轴正方向传播的情况，速度势 φ (x,z)为

式中，A为波浪振幅；g为重力加速度；为波浪圆频率；k0为波数。垂向特征函数Z0(k0z)=coshk0(h+z)/coshk0h，圆频率和波数满足色散关系

对于渗透海床上波浪沿x轴正方向传播的情况，根据Dean和Dalrymple[8]的推导，速度势 φp(x,z)为

式中，κ0=κ0r+iκ0i，其实部κ0r为波数；虚部κ0i是空间衰减系数。复波数可采用Mendez和Losada[15]提出的摄动法计算。

根据Darcy定律，可渗透海床内部的渗流速度为

式中，K为海床多孔介质的渗透特性；η为水的动力黏性系数；Ps为海床中的动水压强，满足Laplace方程

这样，区域I、II、III中速度势和区域IV中的动水压强可分别展开如下。

在区域I中

式中，Am(m=0,1,···,N)为待定系数，N为选取的项数；非传播模态的垂向特征函数Zm(kmz)=coskm(h+z)/coskmh，km是下述色散方程的正实根

在区域Ⅱ中

式中，Bm和Cm(m=0,1,···,N)是待定系数；非传播模态的垂向特征函数是下述色散方程的复数根

同样可采用摄动法求解。

在区域III中

式中，Dm(m=0,1,···,N) 是待定系数。

在区域IV中

式中，Em(m=0,1,···,M) 是待定系数，M为选取的项数，μm=2mπ/L；ρ为水体密度。

速度势和速度在区域I、II和区域II、III两个交界面上连续的匹配条件为

将式(13)、式(14)左右两边同乘以Zm(kmz)(m=0,1,···,N)，并对z沿水深积分

压强和渗流速度在区域II、IV两个交界面上连续的匹配条件为

将式(19)直接对x沿渗透海床长度积分，式(20)左右两边同乘以 co sµmx(m=1,···,M)后对x沿渗透海床长度积分

综上可得4×N+M+5个线性方程，进而确定4×N+M+5个展开系数。

2.2 收敛性分析

本章节分析非传播模态截断项N和渗透海床内部动水压强截断项M的取值对结果收敛性的影响。关于非传播模态截断项N，Kirby和Dalrymple[16]指出只有N的值取足够大，才能保证解的收敛性。他们通过研究波浪过海沟地形的问题，发现对于绝大多数的入射波而言，取N=16时就能保证收敛性，得到足够准确的结果。这一结论被应用于Bender和Dean[17]的研究工作。本文在后续的计算中，同样取N=16。

图2 动水压强截断项M 取值分析Fig.2 Analysis of the value of the pore pressure truncation itemM

关于可渗透海床内部动水压强截断项M的取值问题，本文通过增加M的值，计算z=−h处动水压强的变化情况，进而确定M的取值。计算时，水深h=15 m，入射波波长l=40 m（相应地波浪圆频率ω=1.229 6 Hz），渗透海床长度分别取L=6l和L=1 000l，海床渗透系数ks分别取 0.05 m/s和 0.2 m/s，则z=−h处无量纲化的动水压强Ps/ (ρgh)在渗透海床长度方向上（x / L）的变化情况如图2所示。可以看出，当M=20时，z=−h处动水压强与M=100时基本一致，本文认为取M=20就能保证可渗透海床内部动水压的收敛性。

3 计算结果和讨论

本节应用上述解析方法讨论海床渗透系数ks、 相对水深h/l、渗透海床相对长度L/l这3个因素对波高衰减和波浪反射、透射的影响。

3.1 渗透系数ks的影响

在本算例中，水深h=15 m，入射波波长l=40 m（相应地波浪圆频率ω=1.229 6 Hz），渗透海床长度L=6l。考虑到粗砂、细砂、珊瑚礁等渗透海床的渗透系数数量级介于 10−1～10−3m/s[18–19]，本算例首先计算渗透海床的渗透系数ks=0.5 m/s、0.2 m/s、0.08 m/s、0.05 m/s情况下相对波高的沿程分布；同时，为了检查本文解析模型的普适性，也考虑了ks=0，即区域II中的海床为不可渗透的特殊情况。不同渗透系数情况下复波数计算结果见表1，相对波高如图3所示。然后，计算反射系数Kr=|A0|和透射系数Kt=|D0|随ks的变化情况，结果如图4所示。

表1 不同渗透系数情况下复波数计算结果Table 1 Complex wavenumber calculated for different permeability coefficients

图3 不同渗透系数情况下相对波高分布情况Fig.3 Wave height distribution for different permeability coefficients

图4 反射系数和透射系数随渗透系数变化情况（h/l=0.375，L=6l）Fig.4 The change of reflection coefficient and transmission coefficient with permeability coefficients (h/l=0.375,L=6l)

可以看出，当ks=0时，相对波高等于1，反射系数为0，透射系数为1，即波高沿程不变；而当ks≠0时，区域I和区域II的相对波高存在波动，这是由于波浪反射造成的；波浪传过局部可渗透海床，相对波高沿程按exp[−ki(x+L/2)](−L/2≤x≤L/2)指数趋势衰减，渗透系数越大，复波数虚部ki也越大，相应地衰减程度也越 大，同时，局部可渗透海床对波浪的反射更强烈。

3.2 相对水深h/l的影响

在本算例中，渗透海床长度L=6l，渗透系数ks=0.5 m/s，入射波波长l=40 m。首先计算水深h=20 m、15 m、10 m，即h/l=0.5、0.375、0.25情况下相对波高的沿程分布。不同相对水深情况下复波数计算结果见表2，相对波高如图5所示。接着，计算反射系数Kr和透射系数Kt随h/l的变化情况，结果如图6所示。

表2 不同相对水深情况下复波数计算结果Table 2 Complex wavenumber calculated for different water depth

图5 不同相对水深情况下相对振幅分布情况Fig.5 Wave height distribution for different water depth

图6 反射系数和透射系数随相对水深变化情况（ks=0.5 m/s，L=6l）Fig.6 The change of reflection coefficient and transmission coefficient with water depth (ks=0.5 m/s,L=6l)

可以看出，随着相对水深h/l的减小，复波数虚部ki增大，波高衰减加剧，且波动加强；反射系数随相对水深的减小而增大，透射系数随相对水深的减小而减小。这说明相对水深的减小增大了渗透海床对波浪的影响。

3.3 渗透海床相对长度L/l的影响

渗透海床相对长度L/l的影响主要探讨反射系数Kr和透射系数Kt随L/l的变化规律。在本算例中，水深h=15 m，入射波波长l=40 m（相应地波浪圆频率ω=1.229 6 Hz），渗透系数ks分别取 0.10 m/s和 0.50 m/s，渗透海床长度L介于0～1 000l之间变化，则反射系数Kr和透射系数Kt的计算结果见图7a和图8a。为了进一步观察变化规律，我们对渗透海床长度0≤L≤6l范围内的计算结果进行了局部放大，结果如图7b和图8b所示。

从图7可以看出，反射系数很小，但呈现出振荡变化的规律：当渗透海床长度时，反射系数达到谷值；而当渗透海床长度2,3,4,5)时，反射系数达到峰值。随着渗透海床长度的不断增加，反射系数振荡减小，并趋于常数。

图7 反射系数随渗透海床相对长度变化情况Fig.7 The change of reflection coefficient with the length permeable seabed

从图8可以看出，透射系数随渗透海床长度的增加 而按exp(−kiL)指数减小。

图8 透射系数随渗透海床相对长度变化情况Fig.8 The change of transmission coefficient with the length permeable seabed

4 结论

局部可渗透水平海床由有限长度的可渗透海床及其两侧半无限长不可渗透海床所组成，水深为常数。本文基于势流理论，建立了线性波浪在局部可渗透海床上传播的解析模型，研究结果表明：

（1）波浪在传过局部可渗透海床的过程中，波高沿程衰减，衰减程度随渗透系数的增大、相对水深的减小和渗透海床长度的增加而变大。

（2）受“海床局部可渗透”的影响，波浪发生反射和透射。反射强度小，反射系数随渗透海床长度的变化而呈现有规律的振荡：当渗透海床长度为入射波半波长的整数倍时，反射系数达到谷值，而当渗透海床长度为入射波1/4波长的奇数倍时，反射系数达到峰值。但随着渗透海床长度的不断增加，反射系数振荡减小，并趋于常数。透射系数则随渗透海床长度的增加而指数减小，在渗透海床足够长的情况下，透射系数趋于0。