一种计算波浪与不规则结构物相互作用的方法技术

本发明专利技术公开了一种计算波浪与不规则结构物相互作用的方法，其特点是：该方法是基于雷诺时均N‑S方程组模拟粘性波动流体，利用湍流模型封闭N‑S方程组，作为计算方法的控制方程；通过SIMPLE算法迭代调整压力，使得内部流体单元满足连续方程，自由表面单元满足自由表面动力边界条件，同时调整速度场；同传统的计算波浪与不规则结构物相互作用的方法相比较，本发明专利技术基于考虑流体粘性的N‑S方程、采用矩形网格划分计算区域，具有数值结果更符合物理实际的优点。

A method for calculating the interaction between waves and irregular structures

The invention discloses a method for calculating the interaction between waves and irregular structures. The method is based on the Reynolds time averaged N S equations to simulate the viscous wave fluid, and to use the turbulence model to close the N S equation as the control equation, and to adjust the pressure through the SIMPLE algorithm to make the internal flow. The body element satisfies the continuous equation, the free surface element satisfies the free surface dynamic boundary condition and adjusts the velocity field simultaneously. Compared with the traditional method for calculating the interaction between the wave and the irregular structure, the invention is based on the N S equation which considers the fluid viscosity and the rectangular lattice. Conforms to the advantages of physical reality.

【技术实现步骤摘要】
一种计算波浪与不规则结构物相互作用的方法
：本专利技术涉及防波堤
，具体地讲是一种计算波浪与不规则结构物相互作用的方法；该专利技术可用于波浪与圆柱体、弧形防波墙、海洋平台、半圆型防波堤等海洋工程不规则结构物相互作用问题的研究。
技术介绍
：随着海洋经济的快速发展和国家大力推动“海上丝绸之路”经济带建设，许多大型港口、码头、海洋平台、跨海大桥、浮式飞机场和人工岛等正逐步向海况恶劣、地质条件复杂的较深海域发展。以斜坡式和沉箱重力式为主的传统防波堤结构型式存在工程造价高、施工困难、不利于掩护海域水体交换等缺点，已无法满足海洋工程建设对防波堤结构型式及其功能的要求。近年来，专家学者们陆续开展新型防波堤结构研究，并提出各种型式的透空式和浮式防波堤。为满足不同结构与功能要求，各种新型防波堤结构的几何外形呈现不规则化。关于波浪与此类不规则结构物相互作用的研究，部分学者采用基于势流理论的边界元法(Liuetal.,2009；Ningetal.,2014；Zhouetal.,2015；王双强等,2016)开展。此方法将流体概化为无粘、无旋的理想流体，与实际流体运动的有旋性、漩涡的扩散性和能量的耗散性等存在较大不同，忽略流体粘性作用的势流模型，应用受到一定限制。另有学者采用贴体网格法(Liu,etal.,2009；Li,etal.,2011；李雪艳等,2013)开展；此方法中贴体网格生成较为复杂，且需对控制方程进行坐标变换，影响模型的计算效率。
技术实现思路
：本专利技术的目的是克服上述已有计算方法的不足，提供一种计算波浪与不规则结构物相互作用的方法；主要解决现有的方法应用受限及影响模型的计算效率等问题。本专利技术的技术方案是：一种计算波浪与不规则结构物相互作用的方法，其特殊之处在于，包括以下步骤：a根据设计的数值水槽尺寸和不规则结构物所在位置，生成计算区域的网格；所述的网格由网格单元组成，包括空网格单元、内部流体网格单元和近壁区流体网格单元；b对计算区域的流场，设置初始条件；所述的初始条件包括速度和压强的初始值；初始条件指的是计算域内所有网格点上的速度均设置为零，全域的压强值按静水压强分布设置；c基于雷诺时均N-S方程组模拟粘性波动流体，利用湍流模型封闭N-S方程组，作为计算方法的控制方程；所述的雷诺时均N-S方程组包括动量方程和连续方程。d利用有限差分法对控制方程进行离散；所述的有限差分法的差分网格采用交错网格；e设置计算区域的边界条件，所述的边界条件包括造波边界条件、开边界条件、数值水槽上下边界条件；f将设定速度和压强的初始值代入动量方程，求出计算区域新时刻速度场的近似解；g通过SIMPLE算法迭代调整压力，使得内部流体网格单元满足连续方程，自由表面单元满足自由表面动力边界条件，同时调整速度场；h利用壁函数技术计算近壁区流体网格单元的紊动动能k和紊动耗散率ε；i利用上一时刻的流体体积函数F值和已知的速度场，求出新时刻的F值，据此确定流体自由表面的位置；j利用虚拟边界力法模拟波浪与不规则结构物之间的相互作用；k进行数值稳定条件和收敛条件的判断，满足条件则输出压力场和速度场的结果；l重复上述e至k步骤所述过程，直至计算时间达到程序所设定的总时间。进一步的，所述的雷诺时均N-S方程组包括增加虚拟边界力项的水平方向时均动量方程1、竖直方向时均动量方程2和连续方程3；所述的湍流模型为K-ε模型，由K方程4和ε方程5组成；其中，u为x方向的速度分量，v为y方向的速度分量，t为计算时间，gx为水平方向重力加速度，取值为零，gy为垂直方向的重力加速度，取值为9.81N/kg，p为流体压力，ρ为流体密度，ν为流体运动粘滞系数，是紊动粘性系数，k为紊动动能，ε紊动耗散率；fxvbf和fyvbf分别为虚拟边界力在x和y方向的分量，θ为部分单元体参数，即结构物在整个网格单元中所占的面积与网格单元总面积的比值，范围在0～1之间；其它参数Cu＝0.09，Cε1＝1.43，σk＝1.0，σε＝0.1643，Cε2＝1.92。进一步的，所述的交错网格是指压力与速度交错布置，压力pi,j位于网格单元中心点；和分别位于网格单元右边界和上边界中心点。进一步的，所述的利用有限差分法对控制方程进行离散，水平方向时均动量方程的差分格式详见方程6：其中，δt为时间步长，FUX代表水平方向对流项，FUY代表竖直方向对流项，VISX和TUBX分别代表运动粘性项和紊动粘性项，fxvbf为x方向的虚拟边界力项；水平方向对流项边界网格点采用一阶迎风格式和二阶中心格式线性组合的偏心差分格式，见方程7和方程8：其中α是控制迎风差分量的参数，当α取值为零时，上述差分方程式为二阶中心差分；当α取值为1时，上述差分方程退化成一阶迎风格式；式中sign是符号函数的记号：其中，水平方向对流项内部网格点采用三阶迎风差分格式；当网格右侧边界中心点的水平速度大于零时，对流项差分格式，见方程9：当网格右侧边界中心点的水平速度小于零时，对流项差分格式，见方程10：其中，当网格右侧边界中心点的竖直速度v*＞0时，具体计算式见方程11；当网格右侧边界中心点的竖直速度v*＜0时，具体计算式见方程12；网格右侧边界中心点的竖直速度v*，可由相邻网格单元上下边界中心点的竖直速度取平均值获得；具体计算式见方程13；其中，所采用各参数的具体表达式如下：运动粘性项采用二阶中心差分格式，具体表达式，见方程14：紊动粘性项，采用二阶中心差分格式，具体表达式详见方程15：进一步的，所述的利用有限差分法对控制方程进行离散，竖直方向时均动量方程的差分格式，见方程16：方程中FVX，FVY，VISY，TUBY可同理写出；K方程和ε方程采用隐式线性化处理，以保证紊动动能K和紊动耗散率ε恒为正值；隐式线性化处理后的K方程表达式，见方程17：K方程中的水平方向对流项的离散格式，见方程18：其中，同理可写出FKYi,j，运动粘性项和紊动粘性项见方程19和20；隐式线性化处理后的ε方程表达式，见方程21：其中，水平与垂直方向对流项，运动粘性项和紊动粘性见方程22和23；连续方程采用中心差分格式，具体离散形式，见方程24：其中，AR与AT为网格单元右侧边界与上边界可通过流体部分的面积系数；AC是网格单元的体积系数；上述连续方程的离散形式只能作为流场是否收敛的判定条件；计算过程中，为满足连续方程24，需同时对速度和压力进行调整，反复迭代。进一步的，所述的造波边界条件由虚拟的可吸收式造波机产生，包括行进波和抵消造波板二次反射波的附加波；造波板的运动速度Um(t)的计算表达式详见方程25：其中，η0为所需要的余弦波面，η为造波板前的实际波面，ω为行进波和附加波的角频率，W和L为传递函数，表达式见方程26：造波边界条件的差分形式，如方程27所示：进一步的，所述的数值水槽上下边界条件均设为自由可滑移边界条件；上边界条件设置见方程28，下边界条件设置，见方程29：进一步的，所述的SIMPLE算法用于内部流体网格单元，迭代过程中的压力修正，见方程30：δp＝-sβω(30)其中，s为上一时刻连续方程右边的一个不为零的源项，迭代须进行到所有网格上的s不大于0.001为止；ω是为提高计算精度而在压力修正项的右端引入的压力松弛因子，取值1.7；β是与网本文档来自技高网...

【技术保护点】
1.一种计算波浪与不规则结构物相互作用的方法，其特征在于，包括以下步骤：a根据设计的数值水槽尺寸和不规则结构物所在位置，生成计算区域的网格；所述的网格由网格单元组成，包括空网格单元、内部流体网格单元和近壁区流体网格单元；b对计算区域的流场，设置初始条件；所述的初始条件包括速度和压强的初始值；初始条件指的是计算域内所有网格点上的速度均设置为零，全域的压强值按静水压强分布设置；c基于雷诺时均N‑S方程组模拟粘性波动流体，利用湍流模型封闭N‑S方程组，作为计算方法的控制方程；所述的雷诺时均N‑S方程组包括动量方程和连续方程；d利用有限差分法对控制方程进行离散；所述的有限差分法的差分网格采用交错网格；e设置计算区域的边界条件，所述的边界条件包括造波边界条件、开边界条件、数值水槽上下边界条件；f将设定速度和压强的初始值代入动量方程，求出计算区域新时刻速度场的近似解；g通过SIMPLE算法迭代调整压力，使得内部流体网格单元满足连续方程，自由表面单元满足自由表面动力边界条件，同时调整速度场；h利用壁函数技术计算近壁区流体网格单元的紊动动能k和紊动耗散率ε；i利用上一时刻的流体体积函数F值和已知的速度场，求出新时刻的F值，据此确定流体自由表面的位置；j利用虚拟边界力法模拟波浪与不规则结构物之间的相互作用；k进行数值稳定条件和收敛条件的判断，满足条件则输出压力场和速度场的结果；l重复上述e至k步骤所述过程，直至计算时间达到程序所设定的总时间。...

【技术特征摘要】
1.一种计算波浪与不规则结构物相互作用的方法，其特征在于，包括以下步骤：a根据设计的数值水槽尺寸和不规则结构物所在位置，生成计算区域的网格；所述的网格由网格单元组成，包括空网格单元、内部流体网格单元和近壁区流体网格单元；b对计算区域的流场，设置初始条件；所述的初始条件包括速度和压强的初始值；初始条件指的是计算域内所有网格点上的速度均设置为零，全域的压强值按静水压强分布设置；c基于雷诺时均N-S方程组模拟粘性波动流体，利用湍流模型封闭N-S方程组，作为计算方法的控制方程；所述的雷诺时均N-S方程组包括动量方程和连续方程；d利用有限差分法对控制方程进行离散；所述的有限差分法的差分网格采用交错网格；e设置计算区域的边界条件，所述的边界条件包括造波边界条件、开边界条件、数值水槽上下边界条件；f将设定速度和压强的初始值代入动量方程，求出计算区域新时刻速度场的近似解；g通过SIMPLE算法迭代调整压力，使得内部流体网格单元满足连续方程，自由表面单元满足自由表面动力边界条件，同时调整速度场；h利用壁函数技术计算近壁区流体网格单元的紊动动能k和紊动耗散率ε；i利用上一时刻的流体体积函数F值和已知的速度场，求出新时刻的F值，据此确定流体自由表面的位置；j利用虚拟边界力法模拟波浪与不规则结构物之间的相互作用；k进行数值稳定条件和收敛条件的判断，满足条件则输出压力场和速度场的结果；l重复上述e至k步骤所述过程，直至计算时间达到程序所设定的总时间。2.根据权利要求1所述的一种计算波浪与不规则结构物相互作用的方法，其特征在于，所述的雷诺时均N-S方程组包括增加虚拟边界力项的水平方向时均动量方程1、竖直方向时均动量方程2和连续方程3；所述的湍流模型为K-ε模型，由K方程4和ε方程5组成；其中，u为x方向的速度分量，v为y方向的速度分量，t为计算时间，gx为水平方向重力加速度，取值为零，gy为垂直方向的重力加速度，取值为9.81N/kg，p为流体压力，ρ为流体密度，ν为流体运动粘滞系数，是紊动粘性系数，k为紊动动能，ε紊动耗散率；fxvbf和fyvbf分别为虚拟边界力在x和y方向的分量，θ为部分单元体参数，即结构物在整个网格单元中所占的面积与网格单元总面积的比值，范围在0～1之间；其它参数Cu＝0.09，Cε1＝1.43，σk＝1.0，σε＝0.1643，Cε2＝1.92。3.根据权利要求1所述的一种计算波浪与不规则结构物相互作用的方法，其特征在于，所述的交错网格是指压力与速度交错布置，压力pi,j位于网格单元中心点；和分别位于网格单元右边界和上边界中心点。4.根据权利要求1所述的一种计算波浪与不规则结构物相互作用的方法，其特征在于，所述的利用有限差分法对控制方程进行离散，水平方向时均动量方程的差分格式，见方程6：其中，δt为时间步长，FUX代表水平方向对流项，FUY代表竖直方向对流项，VISX和TUBX分别代表运动粘性项和紊动粘性项，fxvbf为x方向的虚拟边界力项；水平方向对流项边界网格点采用一阶迎风格式和二阶中心格式线性组合的偏心差分格式，见方程7和方程8：其中α是控制迎风差分量的参数，当α取值为零时，上述差分方程式为二阶中心差分；当α取值为1时，上述差分方程退化成一阶迎风格式；式中sign是符号函数的记号：其中，水平方向对流项内部网格点采用三阶迎风差分格式；当网格右侧边界中心点的水平速度大于零时，对流项差分格式，见方程9：当网格右侧边界中心点的水平速度小于零时，对流项差分格式，见方程10：其中，当网格右侧边界中心点的竖直速度v*＞0时，具体计算式见方程11；当网格右侧边界中心点的竖直速度v*＜0时，具体计算式见方程12；网格右侧边界中心点的竖直速度v*，可由相邻网格单元上下边界中心点的竖直速度取平均值获得；具体计算式见方程13；其中，所采用各参数的具体表达式如下：运动粘性项采用二阶中心差分格式，具体表达式，见方程14：紊动粘性项，采用二阶中心差分格式，具体表达式，见方程15：5.根据权利要求1所述的一种计算波浪与不规则结构物相互作用的方法，其特征在于，所述的利用有限差分法对控制方程进行离散，竖直方向时均动量方程的差分格式，见方程16：方程中FVX，FVY，VISY，TUBY可同理写出；K方程和ε方程采用隐式线性化处理，以保证紊动动能K和紊动耗散率ε恒为正值；隐式线性化处...

【专利技术属性】
技术研发人员：李雪艳，王庆，朱小松，战超，张金芝，
申请(专利权)人：鲁东大学，
类型：发明
国别省市：山东,37