光滑粒子流体力学方法:自适应分辨率算法中修正项的重要性
摘要
TREEASPH是一种新代码,这种代码可以形成重力影响流,不管他有没有无碰撞组件。在TREEASPH中,万有重力由分层树形算法来计算,流体力学性能由包含在空间分辨率h(t,r)不是常数时出现修正项的SPH方法来计算。另一个重要的特点是,时间步是用为了允许单个时间步修改的PEC项目(Predict-Enaluate-Correct)在顺序和矢量计算机上来执行的。一些作者以前注意到修正项在避免在模拟中采用非物质的熵上是很需要的。在使用了TREEASPH后,我们注意到,宇宙问题的模拟中,非物质的熵是负值。因此,当我们在模拟时忽略项后,与激振前沿有关的密度峰值会被过高估计。这种结果转而会导致恒星形成过程效率的过高估计。
关键词:方法,数值模拟,流体动力学,星系,构成
1.绪论
几乎所有尺寸的天体物理系统的演变,都是由重力和流体力学所控制的。行星和恒星形成被认为是发生在吸积盘,也表示为崩溃气体云的quasi-hydrostatic平衡。星团的构成和演变由粘性力和星际介质的相互作用所控制,还有一种公认的说法是星系和星系团是耗散的崩溃过程所形成的,在这个过程中辐射冷却和粘性发热是决定他们被观察到结构形状的必要线索。
由于他们的复杂性,此类问题的研究一般要求计算机基于拉格朗日和欧拉流体方法模拟。天体物理学运用的拉格朗日方法大多数是由Lucy、Gingold和Monaghan提出的SPH(Smooth particle hydrodynamics)技术。在SPH方法中,组成系统的流体单元是用粒子抽样表示的,且动力学方程是由拉格朗日形式的流体守恒律所获得的。SPH方法可以实现简单的自适应分辨率尺度,也因此,它可以良好适应很大密度范围的高效模拟系统。其他的物理系统,比如恒星形成、反馈、化学演变等,也都可以很好的应用SPH解决。欧拉方法代替了所谓的Godunve算法,这个Godunve算法在振动捕捉方面相对于SPH方法精确的多。他们的运用自适应网格再定义的新构想尤其适用于天体物理学。欲比较这两种不同的天体物理学的代码表示,参看文献[13,22]。
严格的SPH的可变分辨率尺度的介绍会要求在粒子运动方程中必须包括一些附加项。这些附加项要求h(r,t)的可变性,并经常被叫做“项”。多数SPH代码忽略了这样的附加项,因为在没有很大密度梯度的模拟中,他们对全局结构的影响似乎是可以忽略的。然而,一些作者注意到,如果忽略了项,粒子运动方程就不再保守,SPH方法也采用了一个不明影响的非物质的熵。原则上,这些问题也不排除出现在SPH方法中,然而他们也会出现在其他的任何有可变空间分辨率的流体算法中。
由于SPH模拟已经成为学习很多不同物理方面的必要工具,故而控制每一个可以改变模拟结果的数值的人工因素都很重要。项是否可以被忽略以节约计算时间,或者相反的,项对于获得可信的结构是否是必需的,这些问题都特别需要进行分析。本文我们将描述一种新的代码,TREEASPH,这个代码用分层树形算法(TREEcode)来计算重力,用SPH方法来计算流体性能。我们也可以用这个代码来分析项对于最终结果影响的可能。
本篇论文的计划如下,在第二部分汇总了SPH和treecode方法后,我们在第三部分描述了我们代码的基本方面。第四部分分析了SPH模拟中的项的影响。从该分析中得到的主要结论在最后的第五章。
2.SPH和treecode方法的综述
2.1光滑流体水动力学(SPH)技术
像我们在第一部分提到的,SPH方法包括用气体粒子来表示流体单元,这个气体粒子充当了确定流体特性(如密度、压力梯度等)的插值中心,在每一个粒子i的位置处。这个插值是由核函数给定每个粒子提供的资料不同的比重得到的。光滑长指定了与粒子i有关的插值量的外边界。一般来说,计算是通过以中心为粒子i,半径为,包括一个粒子或附近的固定值(典型的值为32)构成的球体得到的。
对于的密度的光滑估计是由SPH给定的:
, (1)
这里的,表示对称光滑长,
, (2)
对于避免阻碍互易原则是必要的。
则粒子i的运动是由下式决定(参考[3,21,28]):
, (3)
, (4)
这里的,是局部重力势,是比内能,()是压力,是描述了不同于冲击的非绝热过程的冷却(或加热)函数, 和是采用的由Monaghan和Gingold提出的标准格式的人工粘性:
其中 (5)
在这里,是统一命令的常值参数,是为了防止数值分歧的软化参数(通常),是当地声速,,最后.
结合上面对和的表达,我们注意到每一个粒子都有自己的h值,换言之,我们要考虑到可变的光滑长度。在核函数被计算的时候,出现了由h决定的衍生函数,这种情况必须要考虑到。导致公式(3)和(4)结果的不只有空间梯度,还有包括衍生项的附加项。这些形式的附加项叫做项,在公式(2)中的光滑长在对称的情况下由Nelson和Papaloizou给定(从[34]中看其他对称形式)。这样的项可以以紧凑的方式被写为([37]以下文章):
, (6)
这里的和是一定要分别加在公式(3)和(4)内的,其中下标表示k粒子的最远邻近粒子。我们在下面的3.2节将会看到。在TREEASPH中对于的识别已经不需要额外的计算了,因为这个任务已经在光滑长更新的时候已经完成了。
2.2分层树形算法:TREECODE算法
不像流体动力学的相互作用,一个粒子的性能实际上是由适度数量的邻近粒子所决定的,每一个粒子的重力势和加速度是由该粒子和系统中全部的剩余粒子的相互作用所决定的。因此,如果重力加速度是由所有粒子间的直接和算出的(粒子和粒子方法),顺序和矢量电脑计算时间的规模是与N2成比的。故而对于数量较多的粒子的模拟会消耗过多的计算时间。我们需要更多有效的计算重力动力学的产品。例如,基于网格的方法[11,21,31,36],或分层树形算法[8,21,39,45].
在TREEASPH中,重力相互作用由分层树形算法(Treecode)来计算,该算法是由Barnes和Hut[2]引进、由Hernquist用Fortran实现[17-19].与PP和SPH方法相似,Treecode算法是完全的拉格朗日算法并且没有使用会降低空间分辨率和/或会强加几何约束的网格。此外,用来操作粒子群的数据结构可以直接被应用到一些要求SPH计算的确定方面。Treecode算法依托将空间分层细化为规则立方体单元格的方法。单元(或是节点)是包含了系统中所有的粒子的立方体体积。这个单元又被再分成相等体积的八个立方体(三维),这些立方体组成了即时的子节点。这些单元在下一步又被再分成更小的联合体,这个过程会连续的递归直到每一个子单元中只包含了一个或没有粒子。
加在一个粒子上面的力是由附近的粒子相加的(准确的)力,再加上近似于于多重级的远距离单元的力。这些项是由从层级最高的部分(从最大的体积)开始走查决定的。在每一步,我们都要进行检查,这是为了看是否满足接下来的条件:
, (7)
这里的d是从单元重心到粒子的距离,s是单元尺寸,是单元重心到几何形心的距离,是解决公差参数(在这里,)。如果公式(7)满足,所有包括在这个单元中的粒子的影响是用多重级逼近来单元计算的。否则,单元将会被继续再分到遍历该树的下一级直到满足公差标准或达到基本单元。用这个方法,包括树结构和力的估计的所有操作都可以在与NlogN成比例的时间内很快地被执行。
3.TREEASPH
3.1 PEC单个时间步积分方案
由N个粒子组成的系统的模拟有时候在一些范围(或粒子)到另一些的时候对于计算工作量有着很大的区别。例如,密度大且屈服于强冲击波的区域的必须要用比系统剩余部分短很多的时间步来模拟。按绝大多数的文献描述的,系统中所有的粒子在每一时间步是同时前进的。粒子需要的最高时间分辨率决定其他粒子的时间步长。因此,少数的粒子就会减慢整个系统的模拟速度。
为了编写一个解决多种时间尺度问题时更有效的代码,大的计算量就一定要用在需要使用它的粒子上,要避免其他粒子的无用计算。换句话说,允许每一个粒子有不同的时间步是必要的。这样的要求在TREEASPH中已经实现,是通过以如下方式修改一个习惯性的PEC(Predict-Evalute-Correct)集成方案实现的:
1.我们输入全部N个粒子的步数n(相当于的时间)和已知的位置,速度值和加速度值。我们也知道,对于所有的气体粒子,光滑长,特定比内能和他们的衍生项。进一步的,任何一个有独立时间步的的集成方案都需要附加信息来鉴定,在每一步,粒子都需要存储两个向量和。前者包括了每一个粒子对a和再计算执行的时间。第二个向量,,相当于这些变量需要被再计算的最后时间。
2.一个清单是由j个粒子构成的,这j个粒子的值接近于当前时间。这样的粒子被标记为有效的因为它们将马上需要更完整的自身变量再计算。之后我们计算:
. (8)
3.对于每一个粒子,不管有效还是无效,我们都要预测值和
, (9)
, (10)
. (11)
4.只有对有效粒子,我们才用来计算.
5.对于所有的气体粒子,我们用和来计算他们的流体动力学变量.
6.只有有效粒子需要计算和.然后调整:
,
, (12)
,
这里需要以B=1/2得到精确的二阶速度,这时A和C的选择就相对来说任意了。选择A=1/3和C=1/2是为了维持位置和内能的二阶精度。在这些表达中,表示到最后一次估计a和所经历的时间间隔
. (13)
我们注意到,不像,的值对于每个有效粒子是不同的。
7.对于每一个气体粒子,我们用3.4部分的集成方案——辐射冷却他们特定的内能来计算变量。
8.更新全局时间:.
9.我们估计了每个个体的时间步长,又更新了每一个有效粒子的和值。
这套SPH计算概括为上述的PEC方案,还有一些没有被讨论的方面和与我们先前的PPASPH代码不同(至少在一些点)的方面:(1)光滑长的更新,(2)个体时间步长的设置(3)冷却过程的实施(4)恒星形成过程。
3.2 邻近粒子的搜索和光滑长的计算
因为我们的代码使用的是紧支撑的样条核函数,粒子i的SPH特性只由在以内距离的j粒子们决定。另一方面,个体光滑长在每一时间步都会要求在以一粒子i为圆心,为半径,并包含相邻粒子的固定值Ns的球体内更新。因此,SPH的一个基本要求就是有一个搜索相邻粒子的有效算法,无论是在围绕粒子i以(取决于相邻粒子的光滑长)还是以(与无关)为半径的球体搜索范围内。
这些任务的第一步就是在TREEASPH中按如下方式执行。树结构建立后,我们保存了包含在每一个单元内的气体粒子的最大值。然后,在SPH计算中,对于粒子i的局部特性有贡献的粒子(是那些满足的粒子)由如下方式定义。首先,我们用一个边长为的立方体体积(以后都叫搜索正方体)围住粒子i。然后用树形结构向下传递力的计算。在每一层,都用一个边长为(是单元尺寸)的体积包围用树形结构的当前节点表示的体积。然后我们来检查这个体积是否与搜索正方体重叠。如果没有,将不会继续向下传递这种独有的路径。否则,单元将被再分并且向下传递到下一层级。如果当前的单元表示一个粒子,我们就检查看它是否与粒子i的距离在之内,如果有的话就将它记录下来。这个程序将递归地进行下去直到树形结构的所有路径走完。
一个相似的算法被用来计数粒子i的邻近()粒子数Ni,但是现在,用每一个单元表示的体积是由边长为的立方体(不再是)包围的。用了这个Ni值,更新hi的过程就与Hernquist和Katz用过的过程相似了[21]。每一个SPH粒子,我们都由一个预测值开始
, (14)
该预测值作为在时间时的光滑长。我们用这个预测值来计算粒子i的邻近粒子数Ni。如果Ni比Ns多了个规定的公差,那么将会被修正以便让在附近允许范围内降低。
3.3 单个时间步长
为了维持数值积分的稳定性,每一个
剩余内容已隐藏,支付完成后下载完整资料
资料编号:[141683],资料为PDF文档或Word文档,PDF文档可免费转换为Word
以上是毕业论文外文翻译,课题毕业论文、任务书、文献综述、开题报告、程序设计、图纸设计等资料可联系客服协助查找。
您可能感兴趣的文章
- 船舶在浅水航道中航行时的岸壁效应数值研究外文翻译资料
- 基于三维面元法限制水域船体下蹲的数值研究外文翻译资料
- 关于甲板大开口船体梁极限抗扭强度的实验研究外文翻译资料
- 基于斯托克斯方程计算和系统识别 方法预估实船操纵模型参数外文翻译资料
- 水面舰艇5415在PMM演习中的基准CFD验 证数据-第二部分:平均相位的立体PIV流 场测量外文翻译资料
- 初步设计阶段船舶功率推进预测第二部分初步设计中有用的服务速度船舶功率推进数学模型外文翻译资料
- 对某高速船模湍流自由表面的数值与试验研究外文翻译资料
- 第三章水下搜救与恢复操作外文翻译资料
- 液化天然气供求关系的现状与展望:一个全球性展望外文翻译资料
- 基于CFD的高层钢结构建筑风效应数值评估外文翻译资料