英语原文共 140 页,剩余内容已隐藏,支付完成后下载完整资料
DualSPHysics代码的用户指南
DualSPHysics v4.0
2016年4月
目录
1.简介
2.开发商和机构
3. SPH公式
4. CPU和GPU的实现
5.运行DualSPHysics
6. DualSPHysics开源代码
7.编译DualSPHysics
8.格式化文件
9.预处理
10.处理
11.后期处理
12.测试案例
13.如何为您的应用程序修改DualSPHysics
14.常见问题:关于DualSPHysics 的常见问题
15. DualSPHysics v4.0 新增功能
16. DualSPHysics未来
17.参考文献
18.许可证
3.SPH公式
光滑粒子流体动力学(SPH)是拉格朗日无网格方法。该技术使用一组材料点或粒子分离连续体。 当用于流体动力学模拟时,根据周围粒子的物理特性,离散化的Navier-Stokes方程在每个粒子的位置局部整合。 该组相邻粒子由基于距离的函数(圆形(二维)或球形(三维))确定,相关特征长度或平滑长度通常用h表示。 在每个时间步中,为每个粒子计算新的物理量,然后这些粒子根据更新后的值移动。
连续介质流体动力学的守恒定律从其偏微分形式转化为适合基于粒子的模拟的形式,该形式使用基于插值函数的积分方程,该方程给出了在特定点处的值的估计。 通常这个函数被称为核函数(W),可以采用不同的形式,最常用的形式是立方或五次方。 然而,在所有情况下,它都被设计成用r#39;来表示由积分近似定义的函数F(r)
平滑内核必须满足多个属性[Monaghan,1992; Liu,2003],如定义在相互作用区域内的正定性,紧密的支持,规范化和随着距离和可微性单调递减的值。 有关SPH的更完整描述,读者可以参阅[Monaghan,2005; Violeau,2013]。
在方程(1)中,函数F可以近似于一个非连续的、基于粒子集合的离散形式。在这种情况下,这个函数是在一个粒子(a)内插值的,其中对落在其紧密支持区域内的所有粒子进行求和,如由平滑长度h所定义的那样
其中下标表示单个颗粒, Vb是相邻粒子(b)的体积。 如果Delta;Vb= mb/rho;b,其中m和rho;分别是粒子b的质量和密度,则方程(2)变成
3.1 光滑粒子
SPH模型的性能在很大程度上取决于光滑粒子的选择。 粒子表示为粒子之间的无量纲距离(q)的函数,由q = r/h给出,其中r是任意两个给定粒子a和b之间的距离,参数h(平滑长度)控制粒子a周围的区域的尺寸,该尺寸将相邻的粒子考虑在内。 在DualSPHysics中,用户可以从以下内核定义中选择一个:
- 三次样条[Monaghan and Lattanzio,1985]
其中alpha;D在二维模型中等于的10/7,在三维模型中等于的1 /。
[Monaghan,2000]提出的拉伸校正方法仅在一阶导数为零且粒子距离为q的内核的情况下被有效地使用。
b)五次 [Wendland,1995]
其中alpha;D在二维模型中等于的7/4pi;,在三维模型中等于21/16pi;。
在下面的文本中,只考虑了影响域为2h(qle;2)的内核。
3.2. 动量方程
连续系统中的动量守恒方程为
其中Gamma;是指耗散项,而g是重力加速度。 DualSPHysics提供不同的选项来包含耗散的影响。
3.2.1. 人造粘度
由[Monaghan,1992]提出的人造粘度方案在使用SPH的流体模拟中是一种常用的方法,主要是因为其简单性。在SPH表示法中,方程6可以写成
其中Pk和rho;k是对应于粒子k的压力和密度(在a或b处评估)。 粘度项Pi;ab由下式给出
其中和, rk和vk分别是粒子的位置和速度。,是声音的平均速度,和alpha;是需要调整的系数,为了引入适当的耗散,alpha;= 0.01的值已被证明在验证波浪水槽的研究中取得了最好的结果,以研究波浪传播和海浪结构上施加的波浪载荷[Altomare et al., 2015a; 2015c]。
3.2.2. 层流粘度和亚颗粒尺度(SPS)湍流
动量方程中的层流粘性应力可表示为[Lo and Shao, 2002]
其中是运动粘度(对于水来说通常为)。在SPH离散符号中,这可以表示为
亚粒子尺度(SPS)的概念首先由[Gotoh et al., 2001]描述,以表示其移动粒子半隐式(MPS)模型中的湍流效应。 动量守恒方程定义为
其中层流项按照方程9,并且代表SPS应力张量。[Dalrymple and Rogers,2006]使用Favre平均将SPS引入弱可压缩SPH,可将方程11重写为
其中上标是指粒子a和b。
3.3. 连续性方程
在整个弱可压缩SPH模拟的持续时间内(如本文所示),每个颗粒的质量保持恒定,只有它们相关的密度波动。这些密度变化是通过求解SPH形式的质量守恒或连续性方程来计算的:
3.4. 状态方程
在[Monaghan,1994]的工作之后,DualSPHysics中定义的SPH形式中的流体被认为是弱可压缩的,并且状态方程被用于确定基于粒子密度的流体压力。 调节可压缩性,以便人为降低声速; 这意味着在任何一个时刻采取的时间步长(根据Courant条件确定,基于当前计算的所有粒子的声速)可以保持在合理的值。 然而,这种调整将声速限制为比最大流体速度快至少十倍,将密度变化保持在小于1%内,因此不会从不可压缩方法引入大的偏差。 在[Monaghan等,1999]和[Batchelor,1974]之后,压力和密度之间的关系遵循表达式:
其中,b=/,其中参考密度是,参考密度下的声速是
3.5. DeltaSPH
在DualSPHysics中,还可以应用delta-SPH公式,该公式引入了扩散项以减少密度波动。 状态方程描述了一个非常僵硬的密度场,并且与拉格朗日粒子的自然无序一起,发现了高频低幅度振荡以填充密度标量场[Molteni and Colagrossi,2009]。 DualSPHysics在连续性方程中使用扩散项,现在写成
这代表了[Molteni and Colagrossi, 2009]最初的delta-SPH公式,其中自由参数delta;Phi;需要归功于一个合适的值。这种修改可以解释为将密度场的拉普拉斯算子添加到连续性方程中。[Antuono et al., 2012]通过对拉普拉斯算子进行分解,观察算子的收敛以及进行线性稳定性分析来检验扩散系数的影响,对系统中该项的影响进行了仔细分析。这个公式恰好代表了域散布中的扩散项。行为变化接近开放边界,如自由曲面。由于核的截断(没有粒子在开放边界之外被采样),一阶贡献并不为零[Antuono et al., 2010],导致施加到粒子上的净力。这种效应对于非静力条件不被认为是相关的,其中该力比任何其他涉及的力量低许多个数量级。 [Antuono et al., 2010]提出了对这种影响的修正,但涉及密度梯度的重整化问题的解决方案,计算成本相当大。对于大多数应用,建议使用0.1的delta-SPH(delta;Phi;)系数。
3.6. 转换算法
各向异性颗粒间距是SPH中的重要稳定性问题,因为特别是在暴力流动中,颗粒不能保持均匀分布。 其结果是在速度和压力场中引入噪声,以及在某些情况下在水流内产生空隙。
为了对抗各向异性粒子间距, [Xu et al., 2009]提出了一种粒子移位算法来防止不稳定性。 该算法首先为不可压缩SPH创建,但可以扩展到DualSPHysics中使用的弱可压缩SPH模型[Vacondio et al,2013]。 利用移位算法,粒子被移动(“移动”)到具有更少粒子(更低粒子浓度)的区域,从而允许域保持均匀的粒子分布并消除由于噪声可能出现的任何空隙。
[Lind et al., 2012] 提出了对初始换档算法的改进,他们使用Fick的第一扩散定律来控制换档的大小和方向。菲克第一定律将扩散通量与浓度梯度联系起来:
其中J是通量,C是颗粒浓度,DF是Fickian扩散系数。
假设通量即单位时间内通过单位表面的颗粒数量与颗粒速度成比例,则可以找到颗粒移动速度和随后的颗粒移动距离。 使用粒子浓度,粒子偏移距离delta;rs由下式给出:
其中D是一个新的扩散系数,控制移动的大小和吸收比例常数。粒子浓度梯度可以通过SPH梯度算子找到:
比例系数D通过[Skillen etal., 2013]提出的形式计算。 它被设置为足够大以提供有效的粒子移动,而不会引入显着的错误或不稳定性。 这是通过执行对流扩散方程的冯诺依曼稳定性分析来实现的:
其中Delta;tmax是给定局部速度和粒子间距的CFL条件允许的最大局部时间步长。 CFL条件指出:
结合方程19和20我们可以导出一个方程来找到移动系数D:
其中A是独立于问题设置和离散的无量纲常数。 建议使用[1,6]范围内的值,其中2作为默认值。
移位算法严重依赖于完整的内核支持。 然而,自由表面处和其附近的粒子不能获得完整的内核支撑,这将在自由表面预测中引入误差,可能导致非物理不稳定性。 由于自由表面处浓度梯度大,直接应用菲克定律会导致液体颗粒从液体体积中快速扩散。
为了抵消这种影响,[Lind et al., 2012]提出了一种自由表面修正,限制扩散到表面法线,但允许在自由表面上切线。 因此,这种修正只用于自由表面附近,由粒子发散值确定,通过[Lee et al.,2008]首次提出的以下公式计算:
通过将方程(17)的移动距离乘以自由曲面校正系数,将该想法应用于DualSPHysics编码。
其中是自由表面阈值,是粒子散度的最大值。 后者取决于领域维度:
而为DualSPHysics选择自由表面阈值为:
为了确定粒子相对于自由表面的位置,使用了粒子发散与的差异。 因此,具有自由表面校正的完整移动方程(方程17)是:
更多关于变化实施的信息可以在[Mokos,2013]中找到。
3.7. 时间步
DualSPHysics包括数值积分方案的一个选择,如果动量(),密度()和位置()方程首先以下面的形式写入
这些方程使用计算简单的基于Verlet的方案或数值稳定但计算密集的二阶辛方法在时间上整合。
3.7.1. Verlet方案
这种基于普通Verlet方法[Verlet,1967]的算法被分成两部分,并且好处在于与其他一些集成技术相比提供了低计算开销,主要是因为它不需要为每一步进行多次(即预测和校正)计算。 预测步骤根据下式计算变量
其中上标n表示时间步长,和使用方程7(或方程12)和方程 13(或方程14)分别计算。
但是,每隔Ns个时间步长(建议取Ns),变量将根据以下公式计算
第二部分旨在阻止因为方程不再耦合而造成的积分值随时间发散,。 在使用Verlet方案但发现数值稳定性成问题的情况下,增加应用该方案的第二部分的频率可能是明智的,但是如果应该有必要将该频率增加超过Ns = 10,那么这可能表明该方案不能适当地捕获该情况的动态,并且应当使用辛方案。
3.7.2. 辛方案
在没有摩擦或粘性影响的情况下,辛氏积分算法可以时间可逆[Leimkuhler,1996]。 它们还可以保留几何特征,例如运动方程中存在的能量的时间反转对称性,从而改善长期解决方案行为的解决方案。 这里使用的方案是一个明确的二阶辛格式方案,其时间精度为O(Delta;),并且涉及预测器和校正器阶段。
在预测阶段期间,根据下式在时间步的中间估计加速度和密度的值
在校正阶段d期间,用于计算校正后的速度,并因此根据下式计算在时间步骤结束时颗粒的位置
最后使用和的更新值计算密度d/dt=的校正值[Monaghan,2005]。
3.7.3. 可变时间步骤
在明确的时间积分方案中,时间步长取决于Courant-Friedrichs-Lewy(CFL)条件,强迫项和粘性扩散项。 根据[Monaghan等人,1999]使用下式可计算可变时间步长Delta;t
其中Delta;tf是基于每单位质量的力(| fa |),Delta;tcv结合Courant和粘滞时间步长控制。
3.8. 边界条件
在DualSPHysics中,边界由一组被认为是与流体粒子分开设置的粒子描述。 该软件目前提供了坚实的不可渗透和定期开放边界的功能。 也存在允许边界粒子根据固定的强制功能移动的方法。
3.8.1. 动态边界条件
动态边界条件(DBC)是DualSPHysics [Crespo et al., 2007]提供的默认方法。
全文共6905字,剩余内容已隐藏,支付完成后下载完整资料
资料编号:[13579],资料为PDF文档或Word文档,PDF文档可免费转换为Word
以上是毕业论文外文翻译,课题毕业论文、任务书、文献综述、开题报告、程序设计、图纸设计等资料可联系客服协助查找。