基于FLUENT的波浪数值仿真及其对出水物体的作用研究.docx
《基于FLUENT的波浪数值仿真及其对出水物体的作用研究.docx》由会员分享，可在线阅读，更多相关《基于FLUENT的波浪数值仿真及其对出水物体的作用研究.docx（38页珍藏版）》请在冰豆网上搜索。
基于FLUENT的波浪数值仿真及其对出水物体的作用研究
基于FLUENT的波浪数值仿真及其对出水物体的作用研究
姜涛裴金亮唐岱能源学院
指导教师：
陈浮
一、课题研究目的
物体出水运动是一个涉及气液两相问题的三维非定常过程。
在这一过程中，物体的边界条件发生剧烈变化，同时波浪的存在，对物体边界流场的压力、流线分布也起到十分重要的影响。
因此，分析波浪力对于研究水面运动体和出水物体所受应力十分关键。
目前解决该问题的研究手段主要有物理模型实验与数值模拟等。
物理模型实验主要是通过在波浪水槽中进行的实验来研究波浪，采用PIV实验对流场进行跟踪；数值模拟则是通过建立数值模型，通过GAMBIT、FLUENT等CFD软件来进行离散计算。
数值模拟可以节约人力、物力、财力和时间，而且数值模拟可重复性好，条件易于控制，比实验更灵活，此外在海洋结构物的分析和设计中，一般来说，解析解只适用与简单几何形状或线性波浪问题，因而数值解法更有普遍意义。
如果能够对高阶非线性波进行计算模拟，那么就可以用数值波浪水槽模拟各种条件下、特别是极端波况下的波浪运动特性。
所以此项目将采取以数值计算为主，微型实验为辅助的方式开展。
项目分析结果将对解决水下导弹发射等实际工程问题起到参考借鉴作用。
二、课题背景
用计算机模拟取代或部分取代海岸与海洋工程模型试验的设想近些年正逐渐成为现实.与物理模型试验相比,数值模拟不仅成本低,可以避免比尺效应,而且在工况选择以及复杂流场的分析处理等方面也具有明显的优越性.关于数值波浪水池的想法由来已久[1],其实质是构建一个数值模拟平台,在该平台上赋予通常实验室中的波浪水池所具有的功能.
基于势流理论和应用边界元方法构建数值波水池的工作已有不少尝试.目前发展了以时域高阶边界元方法求解完全非线性的势流方程,例如,Kim等和Grilli等的工作.然而,结构物附近由于粘性作用而导致的各种复杂流动状况毕竟不能用势流理论来反映.此外,边界元方法在处理复杂自由水面时难免失效.
自Harlow等提出MAC方法和Hirt等提出VOF方法以来,带自由表面粘性不可压缩流体运动的数值计算技术得到了迅速的发展.在此基础上构建数值波浪水槽的工作也受到了重视.Wang基于VOF方法建立了二维数值波浪水槽并应用所建立的数值波浪水槽开展了波浪对近海平台底部冲击过程的研究.最近,日本一研究小组推出了一个二维的CADMAS-SURF系统,其核心技术是VOF方法.较早将VOF方法推广到三维带自由表面粘性流体运动的是Torrey等.Wang和Su应用改进的VOF方法进行了圆柱容器内液体晃动问题的三维数模
在海洋工程问题中，波浪力是作用在工程结构上的最主要的外力之一。
为合理地设计处于波浪环境中的结构物，保证其安全、有效地运行，需要正确预报其波浪荷载。
另一方面，水下打捞和导弹水下发射等一系列军民项目工程也都受到波浪力的影响。
因此，建立波浪的三维数值模型分析计算便具有了很重要的意义对于波浪与物体的相互作用，人们已经进行了长期的研究。
目前，线性水动力理论已相当成熟，并以其简单、实用而为工程界采用。
然而，波浪与物体相互作用的本质是时域上的非线性问题，实际波浪和物体的运动也是非线性的。
线性理论忽略了实际波浪问题的一些重要的特征，使得很多现象无法解释。
许多海洋工程事故的发生都与非线性因素有关，结构物的破坏往往发生在连续几个大波的极端状况下。
工程实践促使人们将研究的注意力集中到波浪与物体相互作用的非线性特性上。
本课题便基于CFD技术的相关理论，建立存在海风、海流等外界因素下的波浪流体数值模型，为其他工程提供参考
三、课题研究主要内容
1、绪论：
课题的研究背景及意义、国内外在该方向的研究现状综述、课题的主要研究内容；
2、数值模拟方法：
求解器、差分格式、湍流模型、VOF模型、UDF模型、波浪流动特点等介绍；
3、不同波长、波速、水深条件下二维非线性波数值模拟及分析，对波浪的特点及不同波浪参数对波浪品质的影响进行探讨；
4、三维非线性的数值模拟及对运载器出水姿态的影响研究；
四、结论（成果介绍）
波浪的数值模拟过程包括了对造波边界的设置、VOF模型的应用以及数值计算过程中时间步控制等诸多问题。
一方面在模拟过程中要保证数值波浪的品质，从而在充分发展的波浪计算域中能够正确分析其对运载器的影响。
另一方面计算网格又不能过多，否则对于现行计算来说，计算工作量将是巨大的，时间也会很长。
本文所模拟的波浪是二阶STOKES波，主要是因为二阶STOKES波比线性波的精度更高、存在质点漂移现象，与实际海浪更为接近。
本章将通过设置速度入口的方式进行造波，并验证方法的可行性。
在二维计算模型中，将在水深分别为5米和50米的两种工况下对波长9米、波高1米、周期2.4秒的二阶STOKES波进行数值模拟。
分析水深对波浪生成的影响。
其中，主要通过5米水深方案来分析数值波浪的波形参数、压力分布、速度分布以及波浪衰减。
在二维波浪的基础之上，本文将在三维区域探讨二阶STOKES波的数值模拟方法及其可行性。
其中应用的边界条件有：
速度入口、自由出流、壁面以及对称边界。
度入口边界条件：
用于定义流动速度和流动入口的流动属性相关的标量。
这一边界条件适用于不可压缩流，如果用于可压缩流会导致非物理结果，这是因为它允许驻点条件浮动。
自由出流边界条件：
该边界条件用于模拟前无法知道出口速度或者压力的情况;出口流动符合完全发展条件，出口处，除了压力之外，其他参量梯度为0。
壁面边界条件：
对于粘性流动问题，FLUENT默认的设置是壁面无滑移条件但也可指定壁面切向速度分量，给出壁面切应力，从而模拟壁面滑移。
对称边界条件：
对称边界条件应用于计算的物理区域是对称的情况。
在对称轴或者对称平面上，没有对流通量，因此垂直于对称轴或者对称平面的速度分量为0。
在对称轴和对称平面上，没有扩散通量，即垂直方向上的速度为0。
1、二维区域5米水深波浪模拟分析
在二维区域、5米水深的条件下进行二阶STOKES波的数值模拟。
模拟过程中要考虑重力的影响，同时考虑了流体的粘性，所以数值波浪在传播过程中将会逐渐衰减。
在该方案中，将着重分析数值波形与理论波形的匹配、参数变化、流场情况以及质点运动轨迹等。
方案设计
考虑所要模拟波浪的波形参数，建立20m
10m的二维计算区域。
采用结构化网格进行区域划分。
所模拟的二维波浪参数为：
波高1米，波长9米，周期为2.4秒，水深为5米。
模型的网格划分
在网格划分的过程中，本模型采用结构化网格。
对于规则区域结构化网格的质量很高，并且比非结构化网格的数量要少。
对于本模型，结构化网格可以保证网格与流体速度方向一致，这种模型在计算过程中更易于收敛。
图1-1计算区域网格划分图
图1-2计算区域局部网格放大图
在计算区域中最小单元格面积是0.08，最大单元格面积是0.14。
总网格数量是17500。
模型的边界条件及参数设置
图1-3计算区域局边界条件设置图
将左侧边界设置为速度入口，右侧边界设置为自由出流边界。
由于计算域上部的气流运动对造波过程的影响较小，所以可以将上边界设置为壁面。
这虽然与实际情况中的无限宽广自由大气不符，但在计算域较大的情况下可以近似代替。
其中上半部分深色区域为气相，下半部分浅色区域为液相。
表1-1流场计算物性参数表
水密度
ρ
998.2
kg/m3
水动力粘度
μ
0.001003
kg/m⋅s
空气密度
ρ
1.225
kg/m3
空气动力粘度
μ
1.7894e-05
kg/m⋅s
在初始时刻设置20m
5m的区域为液体，流场中气液两相物性参数如表1-1所示。
通过二阶STOKES波的速度势函数编写相应的UDF文件，从而控制造波源处的速度变化。
同时要通过二阶STOKES波的波面方程控制速度入口的有效高度。
保证造波震荡的准确性。
造波过程中涉及到气液交接面即自由表面的处理问题，本模型选择VOF模型作为处理该问题的方法，并且考虑流体的重力。
计算过程中时间步控制在0.0001-0.0005之间，否则计算时间将会是巨大的。
数值波浪模拟结果及分析
图1-4是13.6秒时刻的波形图。
从图中可以看到计算域中包含有两个波形，这与9米波长、20米计算域这一计算条件是相一致的。
图1-413.6秒时刻波形图
图1-513.6秒时刻数值、理论波形曲线
在图1-5中，虚线所代表的数值模拟结果与实线所代表的理论波形较为吻合。
在4.3米和13.5米处处形成波峰，两者之间距离与造波时给定的9米波长一致。
理论上，在求解波面方程的过程中是忽略粘性的作用的，在本模型中考虑了实际流体的粘性，所以数值波高较理论波高有一定的衰减。
同时可看出，波峰较为尖突，波谷较为平坦，并且波峰偏离初始平衡位置的距离要大于波谷偏离平衡位置的距离。
这是符合理论上的二阶STOKES波波形的。
X=4m
X=8m
X=12m
X=16m
图1-6不同X位置处波形时变曲线图
如上图1-6所示，为水平方向不同位置处的波形随时间变化曲线。
在X=4m位置处t=1.1秒时刻开始形成波形，其时间间隔接近半个周期。
随着时间的进一步延长，波形逐渐趋于稳定，及波形充分发展。
对比X=8m、X=12m、X=16m处的波形图可以发现，距离速度入口越远其波形的发展时间越长，形成稳定波形也就越缓慢。
从初始时刻，速度入口进行周期性的变化，但由于流体并非刚性体，所以造波源的运动规律需要一定时间才能够传播到相应的位置，加之流体在运动过程中受粘性作用，其波能也在不断衰减的过程中，所以会出现当波场传播到相应位置处时，波形十分微小的现象。
造波源进经过一定时间的周期变化，当波场各位置处所获得的能量与耗散掉的能量达到平衡时，波场中便会形成稳定波形。
图1-7不同X位置处波形时变曲线比较
图1-7表示的是水平方向不同位置处的充分发展波形的比较。
在粘性耗散影响下，距速度入口较远的X=16m处，其波高最小。
为了进一步分析波高沿程衰减的情况，沿水平方向等距的取10个观测点。
在每一点跟踪其波峰和波谷的位置，并记录其周期，如表1-2所示。
在X=0m即造波源处的波高为0.984米，与设计波高相比误差仅为1.6%。
随着波浪的传播，波高逐渐衰减，在X＝18m处波高衰减达62.5%。
图1-8是波高沿程衰减曲线。
表1-2数值波形参数表
X=0
X=2
X=4
X=6
X=8
X=10
X=12
X=14
X=16
X=18
波峰
0.724
0.7
0.6
0.6
0.51
0.6
0.499
0.522
0.5
0.315
波谷
-0.26
-0.23
-0.22
-0.219
-0.126
-0.205
-0.128
-0.114
-0.098
-0.06
波高
0.984
0.93
0.82
0.819
0.636
0.805
0.627
0.636
0.598
0.375
误差
1.6%
7%
18%
18.1%
36.4%
37.5%
37.3%
36.4%
40.2%
62.5%
周期
2.4
2.4
2.2
2.3
2.5
2.4
2.5
2.4
2.33
2.3
误差
0
0
8%
4.17%
4.17%
0
4.17%
0
2.90%
4%
图1-8波高沿程衰减
同样在表1-2中显示，随着波高的沿程衰减波浪周期并没有明显的变化。
从理论上分析二阶STOKES波的周期是波长和水深的函数，与波高无关。
这一点可以从色散方程中看出。
（a）（b）
图1-9水平、垂直方向理论速度分布
（a）（b）
（c）（d）
图1-10入口水平速度分布
（a）（b）
（c）（d）
图1-11入口垂直速度分布
图1-9、图1-10分别表征的是速度入口的理论速度分布和实际模拟结果的速度分布。
在图1-9中每条曲线代表着某一时刻理论入口速度随水深的变化情况，其速度分布是周期性变化的。
通过在流场中测量入口截面在不同时刻的速度分布，可以绘制出图1-10和1-11。
在图1-10和1-11中，速度随水深逐步增大，当达到自由表面以上的气体位置处时，速度迅速减小为零。
同时从水平和垂直速度图中可以看出其周期性的变化规律。
这与理论给定的速度分布情况是一致的，从而验证了本文所应用的造波方法的合理性，以及给定的入口速度能够准确的体现在流场中。
图1-12流线分布
图1-13单一波峰流线分布
图1-12和1-13所显示的是不同位置处的流线分布。
图1-12是经过直线Y=5的所有流谱。
图1-13是单一波峰处的流谱。
在波面运动的影响之下，波峰处形成速度涡。
为了验证二阶STOKES波的质点运动轨迹为类椭圆形状，本文在（4.5,4.25）、（4.5,4.35）、（4.5,4.75）三个位置处设置初速度为零的示踪粒子，跟踪其运动轨迹如图1-14所示。
可以明显看出，各粒子经过相同时间的运动，轨迹为不封闭的类椭圆曲线。
这主要是由二阶STOKES波水平速度在进行泰勒展开后出现常数项引起的。
而不同位置注入示踪粒子，其运动轨迹也不相同。
初始位置水深越深，质点的移动速度就越小。
这一结果同样与理论分析保持了较高的一致性。
从而得出本方案所模拟的二阶STOKES波不仅在波形上符合理论，其内部质点运动所体现的漂移力也同样与理论吻合。
图1-14质点运动轨迹
图1-15压力分布图
图1-16入口压力分布曲线
图1-15、1-16显示的是全流场压力分布图和入口界面的处压力分布曲线。
在压力分布图中，压力分层分布，相对于大气压强水底最高压强层为48000pa。
符合5米水深的压强分布规律。
2、二维区域50米水深波浪模拟分析
本模型将在50米水深模拟二阶STOKES波，并分析其波形、压力分布、流线等流场特征。
方案设计
考虑所要模拟波浪的波形参数，建立40m
60m的二维计算区域。
采用结构化网格进行区域划分。
所模拟的二维波浪参数为：
波高1米，波长9米，周期为2.4秒，水深为50米。
模型的网格划分
图2-1计算区域网格划分图图2-2计算区域局部网格放大图
在计算区域中最小单元格面积是0.26，最大单元格面积是0.3。
总网格数量是30000。
模型的边界条件及参数设置
图2-3计算区域局边界条件设置图
将左侧边界设置为速度入口，右侧边界设置为自由出流边界。
由于计算域上部的气流运动对造波过程的影响较小，所以可以将上边界设置为壁面。
表2-1流场计算物性参数表
水密度
ρ
998.2
kg/m3
水动力粘度
μ
0.001003
kg/m⋅s
空气密度
ρ
1.225
kg/m3
空气动力粘度
μ
1.7894e-05
kg/m⋅s
数值波浪模拟结果及分析
图2-4是33秒时刻的相图，在此时刻，波形充分发展。
整个计算区域较大，所以该图上的波形并不明显。
但从图2-5中可以看出，数值波浪与理论波形匹配效果较好。
由于计算区域的加大，波浪在传播过程中，伴有明显的衰减。
前两个波长中，波形理想，随后波高迅速降低。
在整个计算区域中共出现了4.5个波形，这符合40米区域、9米波长理论上所应存在的波浪数目。
图2-433秒时刻波形图
图2-533秒时刻数值、理论波形曲线
在本模型中，水深大幅增加，在图2-5中的波形匹配依然较好。
证明在水深增加的情况下，本文所用的造波方法同样可以得到高品质的波浪。
从图2-6和图2-7中，可以看出，虽然水深大幅增加，但有效造波速度依然处于自由表面位置处。
在形成波形的的区域带中，速度明显增大，在水深44米以下的部分，速度基本为零。
当然，在大区域计算域中，耗散与衰减将会变得明显，但对于波浪的生成影响不大。
（a）（b）
图2-6水平方向理论、模拟速度分布
（a）（b）
图2-7垂直方向理论、模拟速度分布
（a）（b）
图2-9压力分布图和入口压力分布曲线
（a）（b）
图2-10水平和处置方向速度分布
50米水深的二阶STOKES波模拟过程与5米水深的情况相比，波形同样与理论较为吻合。
造波边界所产生的入射流体速度并没有因为水深的增加而产生数值上的较大变化。
其速度有效区域相应提高，始终保持在自由表面位置处，在波幅宽度的带状区域内作周期性变化。
在图2-10中，可以看到，波面处的速度无论方向如何，其绝对值在流场中较大，为明显的高速区。
3、三维区域5米水深波浪模拟分析
通过以上二维区域的二阶STOKES波的模拟，可以验证本文所提出的造波方法是可行的。
数值波浪无论从波形上还是内部流场的运动规律方面来看，都使与理论分析相一致的。
因此，在二维基础之上，本文将二阶STOKES波的模拟拓展到三维区域。
在三维区域内模拟二阶SOKES波具有重要意义，因为单纯的波浪模拟其应用价值较小，只有将运载器或海工建筑的受力与波浪力进行耦合计算，才对生产实践产生指导意义。
方案设计
建立50m
10m
10m的二维计算区域。
采用结构化网格进行区域划分。
所模拟的二维波浪参数为：
波高1.2米，波长12米，周期为2.8秒，水深为5米。
模型的网格划分
图3-1计算区域网格划分图
在此模型中仍然采用结构化的正方体网格。
在计算区域中最小单元格面积是0.021
，最大单元格面积是0.031
。
总网格数量是183804。
将计算域扩展到三维区域，计算网格数量大幅增加，所以对于计算机的要求也相应提高。
如果对波形品质要求较高，还需要在波面区域加密网格。
三维计算过程中时间步为0.001-0.0001秒。
模型的边界条件及参数
图3-2计算区域初始相图
三维区域的边界条件设置与二维基本相同，左侧为速度入口，右侧为自由出流边界条件，上下边界均为壁面，新增的垂直Z轴的两个面设置为对称面。
表3-1流场计算物性参数表
水密度
ρ
998.2
kg/m3
水动力粘度
μ
0.001003
kg/m⋅s
空气密度
ρ
1.225
kg/m3
空气动力粘度
μ
1.7894e-05
kg/m⋅s
数值波浪模拟结果及分析
图3-317.86秒时刻波形图
在本模型中讨论的仍然是规则的单一波，并没有考虑Z方向的速度变化，三维区域波浪实际是二维波浪在空间区域的平推。
尽管如此，波浪对流场内物体的作用分析过程仍然是无法用二维波浪模型代替的。
图3-3的波形连续平滑，与给定波浪参数一致。
图3-417.86秒时刻压力
图3-517.86秒时刻水平速度分布
图3-617.86秒时刻处置速度分布
三维区域的压力分布与二维情况一致，均为分层分布。
水底最大压力与理论值相一致。
从速度图上可以看出，自由表面位置处速度值较大，形成四个高速区域。
这四个高速区域便是波峰位置，另外流场中存在四个波形也与50米计算域和12米波长这一参数匹配。
三维计算区域的模拟为下一步进行波浪力与运载器耦合计算提供条件。
三维区域更为接近实际，可以在流场中分析运载器在Z方向的受力，这在二维模型中是无法做到的。
通过上述分析，可知三维区域波浪的品质符合要求，与理论波形匹配。
4、波浪对运载器出水过程的的影响
导弹潜射简介
潜射导弹因其兼有潜艇的隐蔽性好、机动范围大、生存能力强的特点和导弹的射程远、打击精度高、破坏性强的优点，自导弹诞生初期，潜射导弹的发展就备受关注，并逐渐成为现在核威慑和实施二次核打击的中坚力量[12]。
而配备常规战斗部的潜射导弹能实施敌防区外军事打击，对敌境内的各种军事目标实施“结构破坏”战，支援陆上或登陆作战。
导弹水下发射要经历筒内弹道、水下弹道、出水弹道、空中弹道几个阶段，其中前三个阶段是水下导弹武器发射的关键。
水下发射时，由于水的密度是空气密度的800倍，力学环境与空气中有很大不同，水的浮力、阻力、粘滞力将呈现水下独有的特点[13]。
随着水深的变化，尾喷管的背压不断发生着变化，从而也使得推力呈现出不稳定的工况。
同时导弹受各种海洋潮流、水深引起的静压载荷、动压载荷、艇速、发射方式的影响，水下弹道不易控制，甚至会直接导致导弹结构的破坏。
潜射导弹的研制工作从上世纪五十年代中期开始，但因导弹的水下发射技术难度大，目前只有俄罗斯、美国、法国等少数几个国家拥有较成熟的潜载导弹发射技术。
在出水过程中，弹体在极短的时间内经过气液交接面，此过程中，导弹受力发生剧烈变化，力学环境极其不稳定。
尤其再受到波浪力的作用，弹体将会受到很大的偏转力矩。
如果不能正确分析波浪波浪力对其影响，导弹很可能在出水过程中偏离预定轨道。
所以分析波浪对其影响显得尤为关键[14]。
方案设计
按导弹在水中所处的状态分主要有裸式导弹（无运载器）、浸湿导弹（注水运载器或称发射舱）和干式导弹（水密运载器）三种。
本文所选择的模拟方案为裸式导弹。
此外本章基于理论方法研究的目的，所采用的模型以实际弹形为基础进行等比缩小。
在这一过程中，主要分析其受力变化情况，并验证模拟方案的可行性。
通过建立20m
10m
5m的立方体区域来模拟弹体出水位置周围的水域，水深为5米，在水下1米处安置弹体。
造波参数为：
波高1米，波长7米，周期为2.1秒。
为了分析波浪力对导弹出水的影响，本模型在前9秒不考虑弹体的自由运动即将其设置为固定刚体。
在这一过程中调整参数，优化波浪品质，直致整个自由表面的波形得到充分发展。
网格划分及边界条件设置
对本模型进行网格划分应避免网格数量过多，整个计算域范围大，所以仅在弹体周围和自由表面处作带状加密。
边界条件的设置与三维区域的边界条件相同，同时将弹体设置为刚体壁面。
网格划分
图4-1计算区域网格剖面图
弹体在运动过程中涉及到动网格技术，由于在本模型中不考虑艇速，即弹体是在静止状态下由固定推力推升出水，所以四面体网格区域在水平方向适当减小。
从而达到减少网格数量的目的。
网格重构区域两侧的采用结构化网格，为保证良好波形，在中间区域进行加密。
在计算区域中最小单元格面积是
，最大单元格面积是
，总网格数量是630924。
对于这一数量级的网格，目前的计算机性能是可以接受的。
边界条件设置
图4-2计算区域初始相图
左侧边界设置为速度入口，右侧边界设置为自由出流，上下边界均为壁面，垂直于Z轴的前后两个边界设置为对称面。
流场中的弹体为壁面边界。
导弹质量为35kg,推力为400N，初速度为0。
弹体出水模拟结果及分析
出水过程全流场分析
（a）t=0.354s
（b）t=0.504s
（c）t=0.634s
（d）t=0.734s
图4-3导弹出水过程图
图4-3显示的是弹体从启动到出水的分时图，从图中可以看到，弹体在出水过程中受到波浪力作用，产生一定角度的偏转。
导弹出水时刻偏转角度与和方向与波浪的相位是密切相关的。
波浪运动是周期性的过程，所以不同周期时刻的波形对弹体的影响是有很大不同的。
在二阶STOKES波浪模型中，由于始终存在漂移速度，所以沿X轴正方向的运动趋势较为明显。
弹体压力分析
（a）t=0.289s
（b）t=0.329
（c）t=0.554
（d）t=0.784
图4-4弹体相、压力图
在图4-4中，显示了导弹在出水过程中四个不同时刻的相图与压力图，其中压力是以大气压为基准的相对压力。
图（a）中的弹体完全浸没在水中，
升级会员 