基于光滑粒子方法的推板式波浪水槽实验研究

基于光滑粒子方法的推板式波浪水槽实验研究吴宗铎;严谨;赵勇;许斐【摘要】波浪数值水槽是波浪研究方面的一个十分重要的实验装置.目前针对造波水槽所展开的数值或实验研究,多将造波系统的设定参数作为波浪的参数,而未考虑设定参数与实际参数之间的误差.利用推板造波的实验装置进行造波,并借助视频拍摄来捕捉造波板的运动,而波浪参数则由波高仪的数据采集来得到.推板造波的过程,则用光滑粒子方法来进行数值模拟.造波过程中,造波板的运动同样采用实测的参数,并计算了不同时刻的波形曲线和水面高度随时间的变化.光滑粒子的模拟结果与实验结果比,吻合的非常好.【期刊名称】《科学技术与工程》【年(卷),期】2019(019)007【总页数】6页(P24-29)【关键词】实测波浪参数;波浪数值水槽;推板造波;光滑粒子【作者】吴宗铎;严谨;赵勇;许斐【作者单位】广东海洋大学海洋工程学院,湛江524088;广东海洋大学海洋工程学院,湛江524088;大连海事大学船舶与海洋工程学院,大连116026;上海船舶研究设计院,上海201203【正文语种】中文【中图分类】O35221世纪是世界公认的海洋的世纪，而研究波浪运动规律对于远洋运输、舰船制造、科学研究等方面，都具有相当重要的理论价值和实际价值。

造波机是波浪的实验研究中必不可少的一种实验装置，在海洋工程实验室中占有极为重要地位。

利用造波机可以结合水池或者水槽完成。

相比大型的水池，实验水槽建造成本低，且能满足大多数二维造波的实验需求。

围绕造波水槽展开的研究也比较多，滕斌[1]等利用完全非线性的波浪势流理论，建立数值波浪水槽的计算模型，实时模拟造波板的造波运动。

谷汉斌[2]等则以Navier-Stokes方程为基本控制方程，用Level-Set法追踪自由水面，模拟推板造波运动。

王修亭[3]在开源程序平台OpenFOAM上建立主动吸收推板造波的水槽模型，模拟了规则波和不规则波的波浪运动。

李晖[4]等则利用Fluent软件模拟了造波和消波的运动。

张淑华等[5]同时利用造波边界法和质量源项法实现造波的模拟，并根据工程需要考虑了水下障碍物的影响。

相比常规的网格计算体系，无网格的光滑粒子方法对网格的依赖程度较低，使用起来也比较灵活，因而在模拟波浪运动方面受到大家的喜爱，尤其是破碎[6]、飞溅[7]、消波[8]、冲刷[9]等强非线性问题中。

然而，常规的水槽造波问题，光滑粒子方法却应用有限[10,11]。

一是因为其追踪粒子的计算特性难以得到准确的波形，二是因为造波实验本身的仪器误差在计算中往往被忽略。

利用已有的水槽进行推板式的造波实验，并利用光滑粒子方法模拟造波的实验过程。

水槽实验中的造波机可以通过控制推板运动规律，来造出特定周期和特定波高的波。

但是造波机实际造出来的波浪与设定的参数之间存在一定误差。

因此，先通过造波实验的数据采集系统和人工测量，得到造波板的实际运动参数。

再利用造波板的实际运动参数去模拟推板造波的运动，并将模拟的结果与实验结果进行对比验证。

1 造波实验介绍1.1 造波板运动规律这里假设所造的波浪为微幅波，即水池中所造波的波幅相对于波长为很小。

根据的线性造波理论，生成给定波浪时造波板的运动速度可由以下传递函数公式计算：(1)式(1)中，T(ω)是造波板运动和给定波浪间的传递函数，表示造波板运动参数与波浪参数之间的传递关系；A为波幅；E为摇板运动的幅值；d代表水深，k代表波数，且有k=2π/L；L为波长。

以波浪的转播方向为正方向，x轴与水池底面重合，y轴与造波板处的垂直位置重合。

造波板的构造如图1所示。

图1 推板造波装置Fig.1 Device of the piston wavemaker根据波浪按水深的分类，波浪类型可以分为深水波、有限水深波与浅水波。

而当L/20 <d<L/2时，波动为中等水深的水波，这时候波浪存在如下色散关系：ω2-kgtanh(kd)=0(2)式(2)中，ω为圆频率，且ω=2π/T，T为周期。

当水深d超过波长一半时，双曲正切函数tanh(kd)≈1，波浪为深水波。

此时，波长L与周期、波数之间则有如下关系: L=(gT2/2π) tanh(2πd/L)。

1.2 实验水槽概况实验所用的水槽海洋工程波流水槽长33 m、宽1 m、高1.1 m。

一侧设置有造波机，基本工作参数为：推板推程为(-350 mm，350 mm)，造波的频率范围为0.2～2.0，最大波高为400 mm。

测试仪器有波高仪，水流流速仪，压强测试仪，总力测试仪器等，在水槽两侧安装有由多孔介质构成的消能设备，用于消除反射波，尽可能地减少反射波对实验数据的影响。

实验室水槽的外观见图2。

图2 波浪水槽照片Fig.2 Picture of water tank1.3 造波机控制系统造波机实验采用由北京尚水信息技术股份有限公司开发的造波软件WaveMake，造波机软件界面包括两部分，一部分为造波软件界面，一部分为网络控制界面。

初始时刻，需为驱动器使能并将造波板复位到原点。

同时设置采样参数。

采样长度设置为2的指数倍，一般设置为1 024。

再经过凑普以后，就可以设置规则波的参数。

波浪参数的设置界面，如图3所示。

图3 规则波参数设置Fig.3 The setting of regular wave coefficients由图3可以看出，规则波是参数设置比较简单。

将波浪参数设置好以后，便可以生产后缀为.wav的文件。

在系统读取完文件以后，即可完成造波。

造波开始后，等波形稳定以后，即可开始数据采集。

采集数据来自水槽中的波高仪所测得，如图2所示。

2 SPH基本理论2.1 基本原理SPH方法是一种用运动的粒子来代替流体的无网格、纯拉格朗日粒子方法。

SPH 方法的主要特点是基于积分插值。

对于物理参数A，可以表示成关于矢量坐标r的函数，即A(r)。

同时，可以构造加权函数W来建立A(r)的数值关系：(3)式(3)中，h称为平滑长度；加权函数W(r-r′,h)称为加权函数或内核，它表示函数A(r)可以描述成附近一定范围内的离散点数值的加权。

如果将表示成离散关系式的话，则有：(4)式(4)中，a代表需要拟合函数的点，b代表a点附近影响区域内的点，mb为b点粒子质量。

将b点数值通过加权函数Wab来进行拟合，可得到a点数值。

2.2 控制方程整个SPH理论涉及如下基本方程。

2.2.1 加权函数加权函数W的形式包括以下几种：高斯型：W(r，h)=αDexp(-q2)；二次样条型：三次样条型：五次样条型：2.2.2 动量方程连续场中的动量守恒方程为P+g+Θ(5)式(5)中，Θ为扩散项，P代表压力梯度项。

2.2.3 连续性方程流体密度的变化用连续性方程表示：aWab(6)式(6)中，a代表朝着a点的梯度方向。

2.2.4 状态方程流体被视为弱可压缩的。

这有助于使用状态方程来确定流体压力：(7)2.2.5 粒子运动方程(8)式(8)中，ε=0.5，ρab为a点和b点的平均值。

2.3 边界条件采用排斥力边界条件，即水粒子在排斥力作用下无法越过固体边界。

其力学特性为：固体边界的粒子，对靠近的流体粒子产生一种类似于分子作用力的内力。

这种内力，将固体边界粒子与流体粒子隔离开，并保持一定间距r。

那么，代表水质点的粒子，将受到排斥力f的作用：f=nR(ψ)P(ζ)ε(z,u⊥)(9)式(9)中，n代表固体边界法向，R(ψ)为排斥力函数，其表达式为(10)式(10)中，q=ψ/2h，A=0.01c2，c为声速。

函数P(ξ)代表平行于墙运动时受到的持续的作用力。

函数ε(z, u⊥)则是根据水深z 和速度u⊥的修正。

速度u⊥为水质点垂直于墙的速度分量。

3 实验数据与计算结果3.1 实验数据整理对于推板造波来说，波浪的参数完全决定于推板运动的幅值与周期。

因此，对推板运动的幅值和周期均作了测量。

3.1.1 推板运动幅值的测量由于实验室条件不具备直接测量推板运动幅值的方法，这里采用人工测量的方法：先将尺子固定在推板运动的滑轨上，待推板启动造波后，开始录制推板运动的视频，视频长度约10个周期，然后记录视频从最左到最右的运动距离。

具体的测量方法如图4所示。

图4 推板幅值的测量Fig.4 The measure of amplitude of the piston wavemaker根据传递函数式(1)，推板的运动幅值可由波数k和水深d决定。

测量的结果如表1所示。

表1 推板运动幅值Table 1 The amplitude of piston plane in the motion周期T/s波数k水深d/m波高H/mE理论值/mE测量值/m11.442 30.5980.10.030 10.040 41.50.905 20.5980.10.067 00.073 020.665 30.5980.10.107 80.101 0从表1中可以看出，理论值与测量值虽然较为接近。

但是它们之间仍然有差距。

因此，在实际计算过程中，为了消除实验值与测量值之间的差距，应尽可能地采用实验值。

3.1.2 波高与周期的测量利用波高仪采集到的数据，可以读出波高与周期的测量值。

根据波浪采集的数据(水位高度随时间的变化)，时间的数据点的间隔为2 ms，可以很容易推算出周期T和波高H。

本实验中，用于测量的波高仪有两支，之间的间距S0为2.72 m，根据这一距离，以及波浪在两个波高仪之间的传播时间，可以推算出波浪传播速度，进而推算出波长。

计算的波长、波高、周期等因素如下表2所示。

表2 波浪特征参数Table 2 The parameters of wave characteristicH测量值/mH理论值/mT测量值/sT理论值/sL测量值/mL理论值/m0.92760.10.9811.7921.5370.858 80.11.471.53.0122.9860.724 10.11.9624.3174.354 从表2中可以看出，理论值与计算值的差异仍然不可忽略。

波高H与周期T为造波时设定的参数，可以看出其真实值仍然和设定的值有差异。

波长L可以通过有限水深波的色散关系[式(2)]得到。

由于有限水深波的非线性较强，可以看出波高H 与周期T的误差，导致波长L有不小的差距。

3.2 波形的SPH模拟这里用SPH方法模拟了周期1 s，水深为0.598 m的推板造波的过程。

所造波浪的所有特征参数均由推板的运动周期T和幅值E所决定，这两参数均来源于测量值。

模拟的结果用光滑散点来表示。

图5 不同时刻的波形对比Fig.5 The comparison of wave form at different instants图5为造波系统在第11 s、13 s、15 s、20 s时的波浪形态。