波浪在渗透海床上传播的数学模型及其验证

刘忠波，孙昭晨，房克照

（大连理工大学 海岸和近海工程国家重点实验室，辽宁 大连 116024）

0 引 言

波浪是近岸海域重要的水动力因素，准确掌握波况对海岸工程建设有非常重要的现实意义.研究波浪传播变形通常有现场观测、室内物理模型实验或计算机数值模拟等方法，相较而言，数值模拟省时、省费用，但其前提是存在较为合适的波浪数学模型.近年来，作为模拟非线性波浪运动的一类高效的数值模型，Boussinesq类水波方程发展很快［1－3］，但是大多模型是在假定海床不透水的情况下给出的，不能直接用来模拟波浪在由砂石、砂砾组成的可渗海床上的传播变形.Cruz等第一个给出了适合渗透海床上的改进型Boussinesq方程［4］，该方程由于仅含弱非线性项，难以描述强非线性波浪运动.Hsiao等、Chen 从欧拉方程出发，分别给出具有二阶色散和二阶非线性的高阶Boussinesq方程［5－6］，其中Hsiao等也给出了一组以沿深度积分平均速度表达的Boussinesq 方程［5］，虽然其含有二阶非线性项，但其色散性与经典Boussinesq方程一致［7］，限制了其在较深水域的应用.为了改善该方程的色散性能，刘忠波等增加了二阶色散项［8］，并从理论方面分析和讨论了不同参数值下相速度及衰减率的变化.Avgeris等直接在二阶非线性Boussinesq方程中耦合可考虑渗透介质存在时的阻力方程，建立了适用于可渗海床情况的Boussinesq类水波方程，并进行了波浪在可渗潜堤上的数值研究，数值结果表明，这一方式是有效的［9］.而这种方法对具有四阶色散精度的Boussinesq方程（适合非渗透情况）是否有效，目前尚未有相关的研究结果发表.本文针对这一问题进行探讨，在文献［1］的基础上，尝试引入Cruz等的动量方程［4］，联立组成适用于可渗海床的Boussinesq方程，从理论上分析方程的相速度以及衰减率，进而建立相关数值模型，并利用数值模拟波浪在渗透潜堤和渗透地形上传播变形，最后通过相关实验结果检验该模型的适用程度.

1 基本方程

刘忠波等推导了一组具有四阶色散精度的Boussinesq类水波方程［1］，该方程表达如下：

式中：下标t表示变量对时间的偏导数为水平二维偏微分算子，η为波面升高，U为沿水深积分的平均速度，h为静水面以下水深，g为重力加速度.

方程（1）和（2）组成了一组高阶Boussinesq水波方程，但该方程是基于海底不可渗这一假设下给出的，不能应用于可渗海床情况.为了考虑渗透海床对波浪运动产生的影响，须耦合考虑渗透海床内流体运动.具体做法如下：首先保证新建方程的连续方程是守恒的，因此在连续方程中引入n·（hsUs）项（Us为渗透介质中沿水深积分的平均速度，n为孔隙率）；其次是需引入渗透介质中流体运动方程，本文采用Cruz等的结果［4］；最后在自由水体中的动量方程中需考虑渗透耦合影响，这里也直接采用Cruz等的结果，也就是引入.最终方程具有以下形式：

式中：hb＝h＋hs，h和hs分别为上层水体的静水深和渗透介质的厚度；β为色散系数，β＝1／15；cr＝1＋cm（1－n），cm为附加质量系数；α＝a1＋分别为线性阻力系数和非线性阻力系数.方程（3）～（5）组成了一组适合渗透地形的高阶Boussinesq水波方程.

2 理论分析

为了解方程的色散性能，在一维情况下对方程的色散关系进行理论分析.不考虑水深和渗透介质厚度对空间的导数，同时忽略方程中所有非线性项，方程（3）、（4）和（5）可写为

式中：r＝hs／h，将（η，u，us）＝（η0，u0，us0）ei（kx－ωt）代入上述方程中，可以得到

式 中：A11＝－ω；A12＝kh；A13＝nrkh；A21＝gk（1＋α1（kh）2＋β1（kh）4）；A22＝－ω（1＋（α1＋1／3）（kh）2＋（β1 ＋α1／3－1／45）（kh）4）；A23＝－nrω（kh）2；A31＝gk（1＋nβr（kh）2／cr）；A32＝－ω（kh）2／2；A33＝－ω（cr＋iα／ω）（1＋r2（kh）2／3）＋nr（kh）2（1＋β＋iβα／（ωcr））.

只有满足det（A）＝0，方程（9）才有非零解，进一步将其与Gu等的色散关系表达式对比［10］：

在保持n＝0.5不变的情况下，图1（a）、（b）、（c）分别给出当s＝nω／a＝0.1，r＝0.5，2.0，5.0时本文方程解（利用det（A）＝0给出的计算结果）与理论解析解（式（11）计算的结果）的对比.由图1可见：

（1）在不同渗透厚度比情况下，方程的相速度与解析解吻合程度很高，这说明方程具有较好的色散性；对比图1（a）、（b）、（c）发现，当量纲一水深h／L0＝1.0 时，随着r的增大，误差也逐渐增大，最大误差为2.5%（图1（c））.而不存在渗透介质时，由方程（1）和方程（2）组成的Boussinesq方程在h／L0＝1.0时的最大误差为2%，这反映出渗透介质厚度对相速度的影响相对较小.

图1 模型量纲一相速度、衰减率与线性波浪理论解析解的比较Fig.1 Comparisons of dimensionless phase velocity and damping rate from the model and the linear wave theory

（2）在衰减性方面，r＝2.0时，在h／L0＜0.6范围内，方程的衰减率与解析解相差不大，超过这一范围后误差迅速增大；而在其他情况下，尤其是r＝5.0时，当h／L0超过0.4以后，本文方程的结果与解析解存在较大误差.因此，理论分析表明，方程的衰减率与解析解的误差相对较大，纵观不同渗透厚度比下的情况，其可接受的最大h／L0为0.4.本文渗透方程的最高阶仅为2阶，这应该是导致衰减率对比效果不佳的主要原因.

3 一维数值模型

采用类似刘忠波等的方法［1］，对一维情况下的控制方程（3）～（5）进行数值求解.即在非交错网格上离散方程，采用高精度差分格式对方程中的时间和空间导数进行近似，采用三阶Adams－Bashforth预报格式、四阶Adams－Moulton 校正格式进行时间积分.这里也将方程（3）～（5）整理成文献［1］中的格式，具体表达如下：

式中

式中：f（x，t）是内部造波源项，具体表达形式可参见文献［11］.

进行时间积分时，预报格式为

校正格式为

求解动量方程时，采用五对角宽带解法［1］.而关于时间t求导的项、为增加程序稳定性而采用的滤波技术以及边界消波技术，均采用文献［2］给出的方法.当两次计算的值（3个变量）相对误差均小于0.000 1时，当前迭代结束，否则，重新利用式（18）～（20）进行校正计算，主要计算流程见图2.

图2 计算流程示意图Fig.2 The schematic diagram of calculation flow chart

4 数值验证

Hsiao等进行了波浪在渗透潜堤上传播变形的实验［5］.实验中，前坡坡度为1∶8.89，而后坡坡度为1∶5.93，最大水深为0.175m，三角形潜堤上最小水深为0.04m，潜堤由直径约为1.9cm的石子组成，石子间孔隙率n＝0.42；实验中的波周期为1s，波高为2.7cm.计算中，采用Hsiao给出的参数值：α1＝5.83s－1，α2＝58.85m，时间步长0.01s，空间步长为0.02m.

数值结果与实验结果的对比见图3，由图可见，当不考虑渗透作用时（文献［1］数值模拟结果），在x＝10.8 m 处及以后，出现较大的误差；随着渗透作用的进一步积累，数值结果与实验结果的相位差越来越大，而且波幅也存在较大差异.采用本文模型模拟时，在大多位置处，二者在相位和幅值上均较为吻合.而在个别位置处，如x＝11.4m 处相位开始出现一定的误差，且幅值上也存在一定差异，所计算的结果比实验结果偏小.总之，对比计算结果可知，本文模型计算结果与实验结果的吻合程度更佳，这说明多引入一个动量方程来适应渗透介质的影响的做法是有效的.

图3 数值计算的波面升高与实验结果的比较Fig.3 Comparisons of the computed surface elevation and experimental results

为进一步验证模型的适用性，以Sawaragi等的实验数据验证本文模型［12］.实验中的水槽长30 m，宽0.7 m，高0.9 m.水槽中铺设一段长3.5 m，厚0.15m 的平底粗颗粒海床，泥沙颗粒粒径采用了1.8cm 和3.07cm 两种，海床上的水深为0.15m.本文选取入射波高Hi为3.58cm，周期T＝1.0s，孔 隙 率n为0.4 进 行 数 值 模 拟.Karunarathna等曾采用该实验数据对Navier－Stokes数学模型进行验证［13］，实验数据取自Karunarathna和Lin的文献.数值模拟结果与实验数据的比较见图4，由图可见，数值结果（波高）沿波浪传播方向存在一定的震荡，但整体来看，与实验数据的吻合程度良好，这反映出本文模型也能较好地模拟波浪在渗透海床上的传播变形.

图4 渗透海床上波浪衰减的数值结果与实验结果的比较（T＝1.0s）Fig.4 Comparison of numerical results against experimental data for the wave damping over porous seabed（T ＝1.0s）

5 结 论

在一组具有四阶色散精度的Boussinesq 方程基础上，直接耦合Cruz等导出的渗透介质中流体的动量方程，构建了一组新的方程.理论分析了方程的相速度和衰减率，并与解析解进行了对比，结果表明该方程可期望用于模拟波浪在渗透地形上的传播变形；进一步在非交错网格下对方程的一维形式进行离散，建立了基于有限差分法的数值模型，利用数值模型模拟了波浪在渗透地形上的传播变形，并将计算结果与实验结果进行了对比，二者的吻合程度较高.这也从侧面反映出本文的改进方法是有效的，为非渗透地形下的其他Boussinesq类方程扩展为适用于渗透海床情况提供了重要参考.

［1］ 刘忠波，邹志利，孙昭晨.加强的适合复杂地形的水波方程及其一维数值模型验证［J］.海洋学报，2008，30（3）：117－125.LIU Zhong－bo，ZOU Zhi－li，SUN Zhao－chen.Enhanced Boussinesq equations for rapidly varying topographies and their one－dimensional numerical validation ［J］.Acta Oceanologica Sinica，2008，30（3）：117－125.（in Chinese）

［2］ Kirby J T.Boussinesq models and applications to nearshore wave propagation，surf zone processes and wave－induced currents ［M］／／ Advances in Coastal Engineering.New York：Elsevier Science，2002：1－41.

［3］ Madsen P A，Fuhrman D R.High－order Boussinesq－type modeling of nonlinear wave phenomena in deep and shallow water［M］／／Advances in Numerical Simulation of Nonlinear Water Waves.Singapore：World Scientific Publishing Co.Pte.Ltd.，2010：245－285.

［4］ Cruz E C，Isobe M，Watanabe A.Boussinesq equations for wave transformation on porous beds［J］.Coastal Engineering，1997，30：125－156.

［5］ Hsiao S，Liu P L F，Chen Y，et al.Nonlinear water waves propagating over a permeable bed［J］.Proceedings of Royal Society，London A，2002，458：1291－1322.

［6］ Chen Q.Fully nonlinear Boussinesq－type equations for waves and currents over porous beds ［J］.Journal of Engineering Mechanics，2006，132（2）：220－230.

［7］ Peregrine D H.Long waves on a beach［J］.Journal of Fluid Mechanics，1967，27（4）：815－827.

［8］ 刘忠波，孙昭晨.波浪在渗透海床上传播的数学模型［J］.中国科技论文在线，2011，6（5）：374－379.LIU Zhong－bo，SUN Zhao－chen.Wave propagating model over a porous seabed ［J］.China Science Paper，2011，6（5）：374－379.（in Chinese）

［9］ Avgeris I，Karambas T V，Prinos P.Boussinesq modeling of wave interaction with porous submerged breakwaters ［C］ ／／ Proceedings of 29th International Conference on Coastal Engineering.Lisbon：ASCE，2004：604－616.

［10］ Gu Z，Wang H.Gravity waves over porous bottoms［J］.Coastal Engineering，1991，15：497－524.

［11］ Gobbi M F，Kirby J T.Wave evolution over submerged sills：Tests of a high－order Boussinesq model［J］.Coastal Engineering，1999，37：57－96.

［12］ Sawaragi T，Deguchi I.Waves on permeable layers［C］／／Proceedings of 23rd Conference on Coastal Engineering.Venice：ASCE，1992：1531－1544.

［13］ Karunarathna S，Lin P.Numerical simulation of wave damping over porous seabeds［J］.Coastal Engineering，2006，53：845－855.