2. 大连理工大学 海洋可再生能源研究中心, 辽宁 大连 116024
2. Offshore Renewable Energy Research Center, Dalian University of Technology, Dalian 116024, China
随着海洋开发的不断深入,多个海洋工程结构的并靠作业日趋普遍。并靠结构物所形成窄缝内的水体共振现象,逐渐得到工程界和学术界的关注。该水体共振主要表现为,在某些特殊频率波浪激励下,缝内水体会产生大幅波高,危及作业安全。
针对两箱体间缝隙内水体的共振现象,Miao等[1-2]用渐近匹配法研究了带狭缝二维双箱的共振现象,给出了狭缝很小时双箱的理论共振频率;Saitoh等[3]对固定双箱窄缝内流体共振现象进行模拟研究,通过观察双箱窄缝之间的波高变化,发现窄缝内水体的共振波高可达入射波波高的5倍,窄缝流体的共振频率随着方箱吃水深度或窄缝宽度的增大而向低频方向移动。Iwata等[4]分别对不同入射波浪作用下三方箱之间窄缝的波高变化进行了试验研究,发现共振频率与箱体的吃水深度、窄缝宽度和箱体个数成一定的函数关系。宁德志等[5-6]对带窄缝双箱和三箱的共振问题开展了数值模拟,研究了非线性波浪作用下箱体数量、宽度、吃水和窄缝宽度等因素对窄缝内水体的共振频率和响应的影响。苏晓杰等[7]针对孤立波与带窄缝两结构物相互作用问题开展模拟,发现相比于规则波窄缝内不会发生明显的水体共振现象,而且箱体宽度、入射波高等因素均会导致透射波高逐渐减小。谭雷等[8-9]研究了不同半径圆底角的双箱体,在实验中发现直角方箱间的共振波高小于同等条件下圆角方箱间的波高,反映了结构尖角处粘性耗散对窄缝共振特性的影响。Lu等[10-11]通过建立粘性流模型和引入人工阻尼层的势流模型来研究双箱体和三箱体窄缝间的流体振动现象,发现这2种模型均能成功模拟出窄缝间的共振频率和共振波高。张婧文等[12]采用CIP方法对固连双浮体模型中窄缝处内水体共振特性进行了模拟研究,发现浮动情况下水面抬高均小于模型固定情况下的结果。
本文将在上述研究基础上,进一步对毗邻双箱体间水体共振的非线性特性开展研究。利用高阶边界元方法建立完全非线性时域数值波浪水槽,在窄缝内引入人工阻尼项等效水体粘性耗散,通过与实验结果对比,确定不同入射波波高对应的阻尼参数。进而通过数值实验研究波浪非线性对窄缝内共振水体运动特性的影响。
1 基本理论与数值方法考虑波浪水槽中的双箱体问题,如图 1所示。建立二维笛卡尔坐标系oxz,其中z=0位于静水面上,z轴向上为正,x轴向右为正。计算域包含自由水面SF,水底边界SD,箱体边界SB和造波源边界SQ。图中H为水槽静水深,W为箱体宽度,h为箱体吃水深度,Wg为两箱体间窄缝的宽度。取3个测试点:点1是箱体Ⅰ的迎浪点,点2是窄缝中间位置水面点,点3是箱体Ⅱ背浪点。
Download:
|
|
假定流体无粘、不可压缩,且运动无旋,则整个流域可用速度势ϕ来描述。为了能在较短计算域内模拟波浪与带窄缝双箱体相互作用问题,并避免在入射边界发生二次反射,本研究采用域内源造波技术产生入射波浪,整个流域控制方程为泊松(Poisson)方程[13]。自由水面边界条件满足完全非线性边界条件,本文采用混合欧拉-拉格朗日方法更新自由水面,并在计算域两端的自由水面分别布置人工阻尼层来吸收出流波浪。在窄缝自由面上引入人工阻尼来近似流体的粘性耗散。在箱体和水底满足固壁边界条件。由此可以得到所研究问题的控制方程和初边界条件如下:
$ {\nabla ^2}\phi = {q^*}\left( {{x_{\text{s}}}, z, t} \right) $ | (1) |
$ \left\{ \begin{gathered} \frac{{{\text{d}}X\left( {x, z} \right)}}{{{\text{d}}t}} = \nabla \phi-{\mu _1}\left( x \right)\left( {X-{X_0}} \right) + k{\mu _2}^2\phi \hfill \\ \frac{{{\text{d}}\phi }}{{{\text{d}}t}} =-g\eta + \frac{1}{2}\left( {{N_s} - 1} \right){\left| {\nabla \phi } \right|^2} - \hfill \\ \;\;\;\;\;\;\;\;{\mu _1}\left( x \right)\phi - 2{\mu _2}{\left( {gk} \right)^{0.5}}\phi \hfill \\ \end{gathered} \right. $ | (2) |
$ {\left. {\frac{{\partial \phi }}{{\partial \mathit{\boldsymbol{n}}}}} \right|_{{S_B}}} = 0 $ | (3) |
$ {\left. {\frac{{\partial \phi }}{{\partial \mathit{\boldsymbol{n}}}}} \right|_{{S_D}}} = 0 $ | (4) |
$ \phi {|_{t = 0}} = \eta {|_{t = 0}} = 0 $ | (5) |
式中:η为波面高度;阻尼系数μ1为计算域左右边界2个阻尼区的阻尼强度,其为波浪角频率ω和空间位置x的函数[6];μ2为窄缝处自由水面的阻尼系数;q*(xs, z, t)=2Vδ(x-xs)为造波源强度,δ(x-xs)狄拉克函数,造波位置x=xs(本文均取xs=0),V为流体质点水平速度,采用二阶Stokes速度解析解:
$ \begin{gathered} V =- \frac{{\partial {\phi _i}}}{{\partial x}} =- \frac{{g{A_i}k\cosh \left[{k\left( {z + h} \right)} \right]}}{{\omega \cosh \left( {kh} \right)}}\cos \left( {kx - \omega t} \right) - \hfill \\ \frac{3}{8}\left( {{N_s} - 1} \right)A_i^22k\omega \frac{{\cosh \left[{2k\left( {z + h} \right)} \right]}}{{{{\sinh }^4}\left( {kh} \right)}}\cos \left[{2\left( {kx-\omega t} \right)} \right] \hfill \\ \end{gathered} $ | (6) |
式中:ϕi为入射波速度势;Ai为入射波波幅;g是重力加速度。本文将对比非线性和线性条件下波浪爬高的不同,用参数Ns来控制,Ns=1时为线性,Ns=2时为非线性;ω是波浪角频率;k为波数,满足线性色散关系:
$ {\omega ^2} = gk\tanh \left( {kh} \right) $ | (7) |
式中:在整个流域内对速度势应用格林第二定理,可得边界积分方程[14]:
$ \alpha \phi \left( \xi \right) = \int\limits_S {\left[{\phi \left( \xi \right)\frac{{\partial G\left( {\xi, \eta } \right)}}{{\partial \mathit{\boldsymbol{n}}}}-G\frac{{\partial \phi }}{{\partial \mathit{\boldsymbol{n}}}}} \right]{\text{d}}S} + \int\limits_\varOmega {{q^*}G{\text{d}}\varOmega } $ | (8) |
式中:α为固角系数;Ω为整个计算域;S为计算域边界,S=SF+ SB,其中SF为自由水面,SB为物面。G是简单格林函数,可以表示为:
$ G\left( {\xi, \eta } \right) = \frac{1}{{2\pi }}\left( {\ln {r_1} + \ln {r_2}} \right) $ | (9) |
式中:ξ=(x0, z0)为源点;η=(x, z)为场点。r1为源点ξ和场点η的距离:
$ {r_1} = \sqrt {{{\left( {x-{x_0}} \right)}^2} + {{\left( {z-{z_0}} \right)}^2}} $ | (10) |
r2为η到ξ关于水底的镜像点的距离:
$ {r_2} = \sqrt {{{\left( {x-{x_0}} \right)}^2} + {{\left( {z + {z_0} + 2h} \right)}^2}} $ | (11) |
采用三节点高阶边界元将积分边界离散为曲线单元,单元上速度势分布用二次形状函数ki(ζ)插值得到[5]。离散后的边界积分方程可表示为:
$ \begin{gathered} \alpha \phi \left( \xi \right)-\sum\limits_{j = 1}^{{N_B}} {\int\limits_{-1}^1 {\sum\limits_{i = 1}^3 {{k_i}\left( \zeta \right){\phi _{ji}}\frac{{\partial G}}{{\partial \mathit{\boldsymbol{n}}}}\left| {J\left( \zeta \right)} \right|{\text{d}}\left( \zeta \right)} } } + \hfill \\ \;\;\;\;\;\;\sum\limits_{j = 1}^{{N_F}} {\int\limits_{-1}^1 {\sum\limits_{i = 1}^3 {{k_i}\left( \zeta \right)G\frac{{\partial {\phi _{ji}}}}{{\partial \mathit{\boldsymbol{n}}}}\left| {J\left( \zeta \right)} \right|{\text{d}}\left( \zeta \right)} } } = \hfill \\ \;\;\;\;\; - \sum\limits_{j = 1}^{{N_B}} {\int\limits_{ - 1}^1 {\sum\limits_{i = 1}^3 {{k_i}\left( \zeta \right)G\frac{{\partial {\phi _{ji}}}}{{\partial \mathit{\boldsymbol{n}}}}\left| {J\left( \zeta \right)} \right|{\text{d}}\left( \zeta \right)} } } + \hfill \\ \;\;\;\;\; - \sum\limits_{j = 1}^{{N_F}} {\int\limits_{ - 1}^1 {\sum\limits_{i = 1}^3 {{k_i}\left( \zeta \right){\phi _{ji}}\frac{{\partial G}}{{\partial \mathit{\boldsymbol{n}}}}\left| {J\left( \zeta \right)} \right|{\text{d}}\left( \zeta \right)} } } + \hfill \\ \;\;\;\;\;\;\sum\limits_{j = 1}^{{N_Q}} {\int\limits_{ - 1}^1 {\sum\limits_{i = 1}^3 {{k_i}\left( \zeta \right){q^*}\left( \zeta \right)G\left( \zeta \right)} } } \left| {J\left( \zeta \right)} \right|{\text{d}}\left( \zeta \right) \hfill \\ \end{gathered} $ | (12) |
式中:ζ为局部坐标系中自然坐标,取值范围(-1, 1);NF为自由水面单元个数;NB为物面单元个数;NQ为造波源边界上划分单元的个数;J(ζ)为雅克比(Jacobian)行列式:
$ J\left( \zeta \right) = \left( {\frac{{\partial x}}{{\partial \zeta }}, \frac{{\partial z}}{{\partial \zeta }}} \right) $ | (13) |
通过将方程(12)中未知量移到方程左端,已知量移动到方程右端,可以建立关于未知量为速度势和速度势方向导数的线性方程组,并简化为如下矩阵方程形式:
$ \mathit{\boldsymbol{X}} = {\mathit{\boldsymbol{A}}^{-1}}\mathit{\boldsymbol{H}} $ | (14) |
式中:X是未知的速度势或速度势法向导数组成的列向量;空间系数矩阵A随着计算域边界变化而改变,当考虑线性问题时,系数矩阵A各个时刻保持不变;H为已知的列向量。
在每一个时间步,物面上的速度势法向导数和自由水面上的速度势是已知的,通过求解上述方程组可以得到该时刻物面上的速度势和自由水面上的速度势法向导数。然后应用四阶Runga-Kutta法对自由水面边界条件进行时间积分,得到下一时刻的自由水面形状和自由水面上的速度势,进行时间步进[6, 14]。不同的是,对于线性问题,需要将方程(2)和(7)中非线性项忽略掉(即式中Ns=1),且积分域只是在静水面和平均物面边界上。
2 数值计算及讨论依据物理模型实验[8-9]选定参数,水槽静水深H=0.5 m,箱体宽度W=0.5 m,吃水深度h=0.252 m,入射波波高H0的范围为1.0~4.5 cm,两箱体间窄缝宽度Wg=0.05 m。在数值计算中,计算域长度取6倍波长,在水槽的左右两端各布置1.5倍波长的阻尼层,造波源位于x=0,箱体Ⅰ左侧面边界为距离造波源2.0倍波长的位置,然后依次按照Wg调整箱体Ⅱ的位置。通过数值收敛性测试,在箱体Ⅰ左侧自由水面上布置80个单元,窄缝内布置4个单元,箱体Ⅱ右侧自由水面上布置40个单元,计算域垂向边界和造波源面分别布置10个单元,箱体各侧壁边界均布置5个单元,箱底布置10个单元,时间步长设为Δt=T/60 s,每个算例模拟40个周期。
首先对数值模拟结果的稳定性进行测试。图 2分别给出了波数kh=1.55的入射波,在波高分别为H0=1.5 cm和H0=4.0 cm两种情况下t=36T和t=40T时的波形对比。图中对应两个箱体所在的位置也给出。从图中可以看出,两个时刻的波面曲线已经完全重合,包括箱体Ⅰ前的反射波、箱体Ⅱ后的透射波和两个箱体之间的窄缝内波面,并且水槽两端的波面基本趋于0,说明两端阻尼层吸收波浪的效果很好。箱体Ⅰ前的波浪已稳定,其反射回去的波浪透过造波源被前端阻尼层完全吸收,而对入射波浪没有产生影响。以上现象说明本模型的模拟结果已达到稳定。图 2中给出了窄缝内水体运动局部放大图,从图中可以看出,窄缝内水质点运动一致,即为活塞运动。除此之外,对比图 2(a)、(b)可以看出入射波高不同时,所得到的波形也有很大的区别。当入射波高较小时(H0=1.5 cm),所得的箱前的立波波峰和波谷形状规则,基本对称且很光滑。但当入射波高较大时(H0=4.0 cm),箱前的立波波峰高且窄,波谷低而宽,显示出很强的非线性。
Download:
|
|
图 3给出了入射波高H0=2.4 cm和阻尼系数μ=0.035时对应窄缝内波高非线性数值结果与实验结果[8]的对比。数值模拟中,波高Hg以测试点的稳定波面时间序列中波峰值与波谷绝对值和的平均来计算。观察结果对比入射波波高H0=2.4 cm非线性共振波高结果与实验结果吻合,说明数值模拟得到的结果是可靠的。
Download:
|
|
入射波波高的变化将导致窄缝内水体特性的变化,即共振波高和阻尼特性的变化。通过高效经济的数值模拟,可以找到不同入射波波高对应窄缝内水体相应的阻尼系数,进而求得窄缝中的共振波高。图 4给出了窄缝宽度Wg=0.05 m,对应共振条件的波数kh=1.55情况下窄缝中心位置无量纲波高Hg/H0随入射波波高H0的变化关系,及本文非线性数值结果与实验数据的比较。图 5给出了阻尼系数μ随不同入射波波高的变化规律。从图 5中可以看出,数值结果和实验结果对比吻合很好,窄缝内水体粘性随着波浪非线性增强而增大,进而导致无量纲共振波高随着波浪非线性增强而减小。
Download:
|
|
Download:
|
|
图 6给出了3个不同入射波高的共振情况对比,随入射波波高的增大, 其共振频率基本保持在同一变化范围内,而窄缝内无量纲共振波高随着入射波波高的增大而减小。由此可以看出,入射波波高的变化仅仅会使共振波高发生变化,而不会对共振频率产生影响。
Download:
|
|
图 7是在共振条件下(即kh=1.55)窄缝中心、迎浪侧和背浪侧的无量纲波高Hg/H0随不同入射波波高H0的变化关系。从图中可以看出,3个不同位置的共振波高与入射波波高均呈现出一定的函数关系,窄缝中心和背浪侧的无量纲化波高Hg/H0随入射波波高H0的增加而减小,而迎浪测的无量纲化波高Hg/H0随入射波波高H0的增加而增加。窄缝中波高线性拟合得到的函数表达式(1.0≤x0≤4.5 cm)为:
Download:
|
|
$ f\left( {{x_0}} \right) =-1.254{x_0}-0.899\;7\sin {x_0} + 8.339 $ | (15) |
背浪侧波高通过线性拟合得到的函数表达式(1.0≤x0≤4.5 cm)为:
$ f\left( {{x_0}} \right) =-0.045\;98{x_0} + 0.503\;4 $ | (16) |
而迎浪测的无量纲化波高Hg/H0随入射波波高H0的增加而增加,线性拟合的函数表达式(1.0≤x0≤4.5 cm)为:
$ f\left( {{x_0}} \right) = 0.073\;54{x_0} + 1.72 $ | (17) |
窄缝中心点的无量纲共振波高随入射波波高的增大而减小,由此反映出随入射波波高的增大,波浪的非线性逐渐增强,窄缝内的粘性耗散增大,因而共振波高呈逐渐减小的趋势。
图 8是共振条件kh=1.55时,迎浪测和背浪侧爬高分别在非线性和线性两种情况下无量纲共振波高随入射波波高的关系对比。从图中可以看出,非线性和线性两种情况下,迎浪侧的无量纲共振波高均随入射波波高的增大而增大,但当H0≥1.5 cm时,非线性的增大趋势明显大于线性结果,并随入射波波高的增大,这种增大的趋势越显著。而背浪测的无量纲共振波高随入射波波高的增大均呈减小趋势,且非线性和线性两种结果大致接近。可见随入射波高增大,双箱体迎浪侧非线性作用越明显,而背浪侧非线性并不显著。这是因为波浪高倍频部分由于波长短、透射性差,及窄缝内水体粘性耗散增强而导致穿过系统传到背浪侧的波浪能量很小。
Download:
|
|
1) 在共振频率处表征粘性耗散的阻尼项随入射波高的增大而增大;
2) 入射波波高的变化仅仅会使共振波高发生变化,而对共振频率影响有限;
3) 随着非线性的增强,窄缝处的无因次共振波高减小,双箱体迎浪侧非线性作用越明显,而背浪侧非线性并不显著,并得到了3个位置处的波高拟合曲线。
应该指出的是,本文只考虑了一种固定窄缝宽度的情况,对于不同窄缝宽度的一般性情况,未来值得进一步研究。
[1] |
MIAO Guoping, ISHIDA H, SAITOH T. Influence of gaps between multiple floating bodies on wave forces[J]. China ocean engineering, 2000, 14(4): 407-422. (0)
|
[2] |
MIAO Guoping, SAITOH T, ISHIDA H. Water wave interaction of twin large scale caissons with a small gap between[J]. Coastal engineering journal, 2001, 43(1): 39-58. DOI:10.1142/S0578563401000268 (0)
|
[3] |
SAITOH T, MIAO G P, ISHIDA H. Theoretical analysis on appearance condition of fluid resonance in a narrow gap between two modules of very large floating structure[C]//Proceedings of the 3rd Asia-Pacific Workshop on Marine Hydrodynamics. Shanghai, China, 2006: 170-175.
(0)
|
[4] |
IWATA H, SAITOH T, MIAO G P. Fluid resonance in narrow gaps of very large floating structure composed of rectangular modules[C]//Proceedings of the 4th International Conference on Asian and Pacific Coasts. Nanjing, China, 2007: 815-826.
(0)
|
[5] |
宁德志, 苏晓杰, 滕斌. 波浪与带窄缝多箱体作用共振现象的模拟研究[J]. 海洋学报, 2015, 37(3): 126-133. NING Dezhi, SU Xiaojie, TENG Bin. Numerical study of fluid resonance induced by wave action on multi-boxes with narrow gaps[J]. Haiyang Xuebao, 2015, 37(3): 126-133. DOI:10.3969/j.issn.0253-4193.2015.03.013 (0) |
[6] |
宁德志, 苏晓杰, 滕斌. 波浪与带有窄缝结构作用的非线性数值模拟[J]. 船舶力学, 2017, 21(2): 143-151. NING Dezhi, SU Xiaojie, TENG Bin. Nonlinear numerical simulation of wave action on objects with narrow gap[J]. Journal of ship mechanics, 2017, 21(2): 143-151. DOI:10.3969/j.issn.1007-7294.2017.02.003 (0) |
[7] |
苏晓杰, 宁德志, 滕斌. 孤立波与带有窄缝结构作用的数值模拟研究[J]. 哈尔滨工程大学学报, 2016, 37(1): 86-91. SU Xiaojie, NING Dezhi, TENG Bin. Numerical simulation of the interaction of solitary waves and structures with narrow slits[J]. Journal of Harbin Engineering University, 2016, 37(1): 86-91. (0) |
[8] |
谭雷.多体海洋结构间窄缝内流体共振的试验和数值研究[D].大连: 大连理工大学, 2014. TAN Lei. Experimental and numerical investigation of fluid resonance in narrow gap between multiple maritime structures[D]. Dalian: Dalian University of Technology, 2014. http://cdmd.cnki.com.cn/Article/CDMD-10141-1015571557.htm (0) |
[9] |
TAN Lei, LU Lin, TENG Bin. Experimental and numerical modelling dissipative effects of resonant waves in confined space formed by floating box in front of vertical wall[C]//Proceedings of 12th International Offshore and Polar Engineering Conference. Shanghai, China, 2014.
(0)
|
[10] |
LU Lin, CHENG Liang, TENG Bin, et al. Numerical investigation of fluid resonance in two narrow gaps of three identical rectangular structures[J]. Applied ocean research, 2010, 32(2): 177-190. DOI:10.1016/j.apor.2009.10.003 (0)
|
[11] |
LU Lin, TENG Bin, CHENG Liang, et al. Modelling of multi-bodies in close proximity under water waves-fluid resonance in narrow gaps[J]. Science China physics, mechanics and astronomy, 2011, 54(1): 16-25. DOI:10.1007/s11433-010-4194-8 (0)
|
[12] |
张婧文, 何广华, 梅蕾, 等.固连双浮体窄缝间水体共振研究[C]//第十四届全国水动力学学术会议暨第二十八届全国水动力学研讨会文集(上册).长春, 2017: 196-201. ZHANG Jingwen, HE Guanghua, MEI Lei, et al. Fluid resonance from narrow gap betweenhorizontal-motion-free twin floating bodies[C]//Proceedings of the 14th National Congress on Hydrodynamies & the 28th National Conference on Hydrodynamies. Changchun, 2017: 196-201. (0) |
[13] |
BRORSEN M, LARSEN J. Source generation of nonlinear gravity waves with the boundary integral equation method[J]. Coastal engineering, 1987, 11(2): 93-113. (0)
|
[14] |
NING D Z, TENG B, EATOCK TAYLOR R, et al. Numerical simulation of non-linear regular and focused waves in an infinite water-depth[J]. Ocean engineering, 2008, 35(8/9): 887-899. (0)
|