此网站使用Cookie。继续使用本网站即表示您同意我们使用cookies。要了解更多信息,请参阅我们的隐私和Cookie政策。关闭此通知
美国天文学会,了解更多

单击此处关闭此覆盖,或按键盘上的“退出”键。

美国天文学会(AAS)成立于1899年,总部设在华盛顿特区,是北美专业天文学家的主要组织。其成员约7000人,还包括物理学家、数学家、地质学家、工程师和其他研究和教育兴趣属于当代天文学的广泛学科的其他人。原子吸收光谱的使命是增进和分享人类对宇宙的科学认识。

https://aas.org/

物理研究所,了解更多

单击此处关闭此覆盖,或按键盘上的“退出”键。

物理研究所(IOP)是一个领先的科学协会,促进物理学的发展,并将物理学家聚集在一起,造福于所有人。它在全世界拥有约5万名会员,包括来自各个领域的物理学家,以及对物理学感兴趣的物理学家。促进物理研究、应用和教育;与政策制定者和公众接触,提高对物理学的认识和理解。其出版公司IOP publishing是专业科学传播领域的世界领先者。

https://www.iop.org

文章

Nyx:一个用于计算宇宙学的大规模并行AMR代码

,,,,和

发布日期:2013年2月13日 ©2013。美国天文学会。版权所有。
, , 引用安·S。阿尔姆格伦2013亚太及日本地区 76539

0004-637X/765/1/39

摘要

我们提出一个新的N-天体和气体动力学代码,称为Nyx,用于大规模宇宙模拟。Nyx遵循了一个离散暗物质粒子系统的时间演化过程,在膨胀的宇宙中引力耦合到一个无粘性的理想流体。该算法采用欧拉框架,采用块结构自适应网格细化;使用相同网格层次的粒子网格方案来求解自引力并推进粒子。计算结果证明了Nyx在标准宇宙学测试问题上的有效性,以及Nyx对50000核的标度行为。

导出引文和摘要 BibTeX公司 RIS公司

1.简介

了解宇宙中物质的分布在整个宇宙历史中是如何演变的是宇宙学的主要目标之一。大尺度结构的性质强烈地依赖于一个宇宙学模型,以及该模型中的初始条件。因此,原则上,我们可以通过将观测到的物质分布与最适合计算的模型宇宙相匹配来确定宇宙的宇宙学参数。说起来容易做起来容易得多。观测,除了引力透镜,不能直接接触到物质的主要形式,暗物质。相反,它们揭示了重子的信息,这些信息以某种偏颇的方式追踪暗物质。在理论方面,描述密度扰动演化的方程(如Peebles1980)只有当扰动幅值很小时,才能很容易地求解。然而,在大多数感兴趣的系统中,局部物质密度ρ可以比平均密度ρ大几个数量级0,做密度对比,δ≡ (ρ ρ0)/ρ0比统一更伟大。

由于这种非线性,以及宇宙初始条件的不确定性,单一的“英雄式”模拟不足以解决大多数感兴趣的宇宙学问题;对不同的宇宙学进行多种模拟的能力,以及对支配物理的不同近似,是必不可少的。半解析模型(例如,White&Frenk1991在星系形成的背景下)或标度关系(例如,Smith等人。2003)为昂贵的模拟提供了一个有趣的替代方案,但它们也必须通过仿真进行验证或校准,以使其可信。最后,除了作为一个无价的理论家工具外,模拟在弥合理论和观测之间的差距方面发挥着至关重要的作用,随着大型天空调查开始收集数据,这一点变得越来越重要。

为了与未来星系/星系团调查的能力相匹配,宇宙学模拟必须能够在一侧模拟约1 Gpc的体积,同时分辨率约为10 kpc甚至更小。显然,在欧拉气体动力学的背景下,这就要求能够在相关物体形成后进行多层次精细化的模拟。解决*每个星系大约有100个质量元素(粒子)在这种规模的模拟中,一个星系需要数百亿个粒子。这样的系统,跨越时间为100万年或更少,必须进化大约100亿年。另一方面,模拟的目的是预测Lyα森林观测值,其动态范围更像~104但由于森林的大部分贡献来自于接近平均密度的气体,因此必须解决大部分的计算体积。这些模拟需要更“浅”的网格细化策略来覆盖更大的区域。

除了实现足够精细化的挑战外,模拟还需要在明显的自引力之外纳入一系列不同的物理。重子的存在导致了对纯重力预测的修正;这些校正在非常大的尺度上可以忽略不计,但在较小的尺度上变得越来越重要,需要对气体进行精确处理,以便提供可直接与观测结果进行比较的实际结果。辐射冷却和加热机制与气体本身的添加对许多可观测到的影响一样大,但在这些方面仍然存在很大的不确定性,随着我们进入高度非线性的系统,我们对知识的缺乏也会增加。因此,我们可以模拟高红移下低密度Lyα吸收剂(例如,Viel等人。2005),我们可以重现宇宙星系团中最大天体的体积特性(参见,例如,Stanek等人。2010),但在更小、密度更大的尺度上,星系的形成仍然难以捉摸,是一个活跃的研究领域(例如,本森)2010).

本文介绍了一种新开发的N-身体和气体动力学代码被称为Nyx。该代码将暗物质建模为一个由拉格朗日流体元素或“粒子”组成的系统,在引力作用下与代表重子物质的无粘理想流体耦合。在欧拉框架下,采用有限体积表示法和块结构自适应网格细化(AMR)对流体进行建模。用于演化流体量的网格结构也用于通过粒子网格(PM)方法来演化粒子。为了更精确地处理高超音速运动,其中动能比气体的内能大很多个数量级,我们采用双能量公式,在每个时间步长内,内能和总能量方程都在网格上求解。

现有的有限体积,AMR代码;宇宙学界最常用的是恩佐(Bryan等人。1995),艺术(Kravtsov等人。1997),闪光(Fryxell等人。2000)和拉美西斯(Teyssier2002). Enzo、FLASH和Nyx都是相似的,它们使用块结构的AMR;相比之下,RAMSES和ART则对单个单元格进行局部优化,从而产生非常不同的细化模式和数据结构。在分块结构的AMR代码中,Enzo和FLASH在补丁之间有严格的父子关系;i、 例如,每个优化的补丁都完全包含在一个父补丁中。FLASH还强制要求所有网格补丁具有相同的大小。Nyx只要求精细斑块的联合包含在粗糙斑块的联合中,并且具有适当的嵌套。同样,这会导致不同的细化模式、数据结构、通信开销和缩放行为。

时间优化策略也各不相同。FLASH以相同的时间步长Δ在每个级别上推进溶液t,而恩佐、阿特和拉美西斯则允许在时间上进行亚循环。一个标准的次循环策略执行Δt在水平上是恒定的;艺术和恩佐在适当的时候允许在时间上比在空间上更精细。纽约证券交易所维持几种不同的亚循环策略。用户可以在运行时指定Nyx不使用子循环、标准子循环、指定的子循环模式或“最优子循环”。对于最优子循环,Nyx在运行期间以指定的粗略时间间隔选择,基于对每个选项当时的计算成本的自动评估(E。Van Andel等人,2013年,编制中)。

在计算宇宙学等已经成熟的领域,开发新的模拟代码的主要目标之一是确保现有和未来硬件的有效利用。未来大规模的宇宙模拟将需要能够有效地扩展到数十万个核心的代码,以及有效地利用多核、共享内存和单个节点的特性。其中一个只用于重力模拟的代码是HACC(Habib等人。2009). 现有的hydro代码(如Enzo和FLASH)的框架目前正在为当前和未来的架构进行大量重写。纽约,就像卡斯特罗辐射流体力学代码(Almgren等人。2010; Zhang等人。2011),以及低马赫数天体物理学的MAESTRO代码(Almgren et al。2006年a,2006年b; Gren等。2008; Zingale等人。2009; Nonaka等人。2010),是基于BoxLib(BoxLib)构建的2011)块结构自适应网格方法的软件框架,因此它利用大量的努力来实现未来仿真所需的大规模并行性能。BoxLib目前支持基于MPI和OpenMP的多核体系结构的分层编程方法。单独的例程,例如那些评估加热和冷却源项的例程,都是以模块化的方式编写的,这使得它们很容易移植到gpu上。

在下一节中,我们将介绍在Nyx中求解的方程。分段我们将对代码进行概述,并在第节中4讨论了结构化网格AMR和多级算法。剖面5详细介绍了并行实现、效率和扩展结果。剖面6包含两种不同的宇宙学测试;第一个是纯粹的N-第二个是圣巴巴拉星系团模拟,它结合了气体动力学和暗物质。在最后一节中,我们将讨论未来的算法开发和即将使用Nyx进行的模拟。

2.基本方程

2.1.宇宙膨胀

在宇宙学模拟中,重子物质是通过在与膨胀宇宙共同运动的坐标系中求解自引力气体动力学方程而演化出来的。这引入了比例因子,,转化为气体动力学的标准双曲守恒定律,其中与红移有关,z,签字人=1/(1)+z),也可以作为时间变量。进化在Nyx中由Friedmann方程描述,对于我们考虑的两分量模型,在哈勃常数下,H0,以及宇宙常数,$\Omega\Uu\Lambda公司$,具有形式

方程(1)

其中Ω0是今天宇宙的总物质含量。

2.2.气体动力学

在共动框架中,我们联系了共动重子密度ρb,至适当密度,ρ适当的,按ρb=ρ适当的,并定义美国是特有的真重子速度。然后将连续性方程写成

方程式(2)

尽可能保持方程的守恒形式,我们可以把动量方程写成

方程(3)

压力,p,与适当的压力有关,p适当的,签字人p=p适当的. 在这里g= −∇phgr公司是重力加速度矢量。

我们使用一种类似于Enzo中使用的双能量公式(Bryan et al。1995),在这个过程中,我们在每一个时间步中进化出内能和总能量。内能的演化方程,e,以及总能量,E=e+(1/2)美国2,可以写为:

方程(4)

方程式(5)

其中∧碳氢化合物表示加热和制冷的组合项。我们同步Ee在每个时间步的末尾;程序取决于e相对于E,将在下一节中介绍。

我们认为重子物质满足伽马定律状态方程(EOS),其中我们计算压力,p= (γ 1) ρe假设组成是氢和氦的混合物,其原始丰度和不同的电离状态是通过假设统计平衡计算的。EOS的细节以及它与加热和冷却条件的关系将在未来的工作中讨论。给定γ=5/3,我们可以将能量方程简化为

方程式(6)

方程式(7)

2.3.暗物质

宇宙中的物质含量(80-90%)强烈地被冷暗物质所控制,它可以被模拟为一种非相对论性的无压流体。相空间函数的演化,f暗物质的,由无碰撞玻尔兹曼(又称Vlasov)方程给出,在膨胀空间中

方程式(8)

哪里p分别是质量和动量,以及phgr公司是重力势。

与直接求解Vlasov方程不同,我们使用蒙特卡罗方法对初始时刻的相空间分布进行采样,定义粒子表示,然后将粒子演化为N-身体系统。由于粒子轨道是哈密顿量的积分,刘维尔定理保证了在以后的某个时间粒子仍然按照方程中给出的相空间分布采样(8). 因此,在我们的模拟中,我们将暗物质表示为带有粒子的离散粒子有共同的活动地点,,以及特有的固有速度,美国,并解决系统问题

方程式(9)

方程式(10)

哪里g重力加速度是在粒子所在位置计算的吗,即。,g=g(,t).

2.4.自重力

暗物质和重子都有助于引力场并被它加速。给出了伴随运动的暗物质密度,ρdm公司,我们定义phgr公司通过解决

方程式(11)

哪里G是重力常数,ρ0是(ρ)的平均值dm公司+ ρb)在域上。

3.代码概述

3.1.宇宙膨胀

在气体和暗物质向前演化之前,我们计算在下一个完整的粗略级别时间步长上。鉴于n定义为当时tn,我们计算n+1通过二阶龙格-库塔台阶。价值选择的相对差异n+1计算方法和2步数小于108. 这种精度通常是通过=8–32,以及推进的计算成本与模拟的其他部分相比微不足道。我们注意到Δt受增长率的限制就这样在最粗糙的水平上,每个时间步的变化不超过1%。一次n+1计算,值为中间时间所需的数据通过线性插值计算nn+1.

3.2.气体动力学

我们把气体的状态描述为美国= (ρb,ρb美国,2ρbE,2ρbe),然后把气体的演化写为

方程式(12)

哪里F=(1)/ρb美国, ρbUU公司,bUE公司+聚氨基甲酸酯),ρbUe公司)是通量向量,Se=(0,0,0,美联社中心点 美国)表示内能演化方程中的附加项,Sg=(0,ρbg,ρb美国 中心点 g,0)表示引力源项,以及S碳氢化合物=(0,0,Λ碳氢化合物,Λ碳氢化合物)表示热源和冷源的组合项。国家,美国,所有源项都定义在单元中心;通量在单元面上定义。

我们计算F使用无照明角点法和无照明耦合法1990; 萨尔茨曼1994). 见Almgren等人(2010)对于没有比例因子的算法的完整描述,; 这里包括以一种保持二阶时间精度的方式在所有方面都是合适的。流体动力积分同时支持非分段线性(Colella1990; 萨尔茨曼1994)以及非分裂分段抛物线法(PPM)格式(Miller和Colella2002)构造用于定义通量的以时间为中心的边状态。如Colella&Woodward中所讨论的,实现的分段线性版本包括包含参考状态的选项(1984)还有Colella&Glaz(1985). 我们选择在向界面传播的单元中最快特性的依赖域上重建轮廓的平均值作为参考状态(如果最快的特征从边缘传播出去,则使用边缘重建轮廓的值。)参考状态的选择将算法依赖线性化方程的程度最小化。特别是,这种选择消除了二阶Godunov算法中使用的边缘特征外推的一个组成部分。对于典型的高超声速流的宇宙学模拟,我们发现参考状态的使用提高了算法的整体鲁棒性。纽约证交所的Riemann解算器使用Colella&Glaz中描述的双激波近似迭代地解决Riemann问题(1985); 对于高超声速宇宙流,这比Almgren等人描述的更简单的线性化近似Riemann解算器更为稳健(2010). 在牢房里e小于0.01%E,我们使用e独立进化,从状态方程计算温度和压力,并重新定义E作为e+(1/2)美国2; 否则我们重新定义e作为E− (1/2页美国2).

引力源项,Sg在动量方程和能量方程中,采用预测-校正方法进行时间离散。ρ的离散化有两种选择b美国 中心点 g在总能量方程中。最明显的离散化是计算单元中心动量的乘积,(ρb美国),以及细胞中心引力矢量,g,每次。虽然这在空间和时间上是二阶精确的,但它可能会产生改变内能的非物理后果,e=E− (1/2)美国2,自更新到E在分析上,但在数值上不等同于用动量更新计算的动能。第二个备选方案定义了E作为(1/2)的更新美国2. Springel对Euler方程中包含自重力的设计选择进行了更全面的回顾(2010).

3.3.暗物质

我们用一个标准的kick-drift-kick序列来演化暗物质粒子的位置和速度,与Miniati&Colella中描述的完全相同(2007). 我们计算g粒子的位置,通过线性插值g从细胞中心到. 为了移动粒子,我们首先使粒子加速Δt/2,然后使用时间中心速度推进位置:

方程(13)

方程(14)

在计算重力之后tn+1,我们完成了粒子速度的更新,

方程式(15)

我们观察到这个格式是辛的,因此平均运动积分守恒(吉田1993). 在计算上,该方案也很有吸引力,因为在时间步长期间不需要额外的存储来保持中间位置或速度。

3.4.自重

我们定义ρdm公司在网格上使用云计算方法,我们假设粒子均匀分布在边Δ的立方体上集中在. 我们将每个粒子的质量分配给网格上的单元,与每个单元与粒子立方体相交的体积成比例,并将这些单元值除以Δ定义移动有助于方程右侧的密度(11).

我们解决phgr公司在同一个细胞中心,其中ρb和ρdm公司定义。拉普拉斯算子采用标准七点有限差分模板离散,所得线性系统采用几何多重网格技术求解,特别是具有红黑高斯-赛德尔松弛的V循环。解算器的公差通常设置为1012尽管数值研究表明这种公差可以放宽。在这里,我们给出了cosmolib-neuichrs边界条件,而cosmolic-lib和boxichrs的周期性边界条件在这里都是假设的。使用先前Poisson解的解作为新解的初始猜测,可以进一步提高效率。

鉴于phgr公司在细胞中心,我们计算g在细胞的中心差phgr公司在相邻的细胞中心之间。这种离散化,结合插值g从细胞中心到上述粒子位置,确保在均匀网格上,粒子对自身的引力相同地为零(Hockney&Eastwood)1981).

3.5.计算时间步长

时间步长是使用显式方法的标准CFL条件计算的,附加的约束来自暗物质粒子和. 使用用户指定的双曲CFL系数,0<σCFL、hyp<1,以及声速,c,根据EOS计算,我们定义

方程式(16)

哪里e是中的单位向量th方向和最大值域中所有计算网格单元的最大值。

由于粒子的时间步长约束要求粒子的移动距离不超过Δ在一个单一的时间步。因为美国n+1/2页当我们计算时间步长Δ时还不知道t当时tn,用于计算Δt,我们估计美国n+1/2页通过美国n. 对暗物质使用用户指定的CFL数,0<σ_CFL、dm<1,它与双曲CFL数无关,我们定义

方程(17)

麦克斯在哪里j是所有粒子的最大值j. σ的默认值CFL、dm=0.5。在初始粒子速度很小,但重力加速度很大的情况下,我们还约束Δtdm公司通过要求gΔt2dm公司< Δ为所有粒子防止快速加速导致粒子移动过远。

最后,我们计算Δt随着时间的推移会改变1%。代码采用的时间步长由以下三个值中的最小值设置:

方程式(18)

我们注意到非零源项也可能约束时间步长;我们将讨论的细节推迟到将来的一份文件。

3.6.单级积分算法

在单一细化级别上的算法从计算时间步长Δ开始t,前进tntn+1=tn+ Δt. 剩下的时间步骤由以下步骤组成:

  • 第一步:计算phgr公司ngn使用ρnb和ρndm公司,式中ρndm公司是根据n.我们注意到在单级算法中,我们可以使用g在上一步结束时计算,因为或ρb从那时起。
  • 第二步:预付款美国通过Δt.如果我们对加热和冷却源项使用预测-校正方法,则在显式更新中包括所有源项:
    方程式(19)
    哪里Fn+1/2页Sn+1/2页e通过预测美国n如果我们用Strang分裂来代替热源项和冷源项,我们首先提出(ρe)和ρE时间项Δ1t
    方程式(20)
    方程式(21)
    然后使用以时间为中心的通量和Se以及Sg当时tn:
    方程式(22)
    哪里Fn+1/2页Sn+1/2页e通过预测美国n,*州。加热和冷却源项如何合并的细节取决于∧表示的加热和冷却机制的具体细节碳氢化合物; 我们将讨论的细节推迟到将来的一份文件。
  • 第三步:插值gn从网格到粒子位置,然后将粒子速度提高Δt/2和粒子位置t:
    式(23)
    方程式(24)
  • 第四步:计算phgr公司n+1gn+1使用ρn+1个,*b和ρn+1dm公司,式中ρn+1dm公司是根据n+1.在这里我们可以使用phgr公司n作为最初的猜测phgr公司n+1为了减少在多重网格中花费的时间达到指定的公差。
  • 第五步:插值gn+1从网格到粒子位置,然后更新粒子速度,美国n+1
    方程式(25)
  • 第6步:对的美国以时间为中心的源项,并替换e通过E− (1/2)美国2视情况而定。如果我们对加热和冷却源项使用预测-校正方法,我们将通过有效地使所有源项居中来校正解决方案:
    方程式(26)
    如果我们用Strang分裂来计算热源项和冷源项,我们只把引力源项作为时间中心,
    方程式(27)
    然后对另一个(1/2)Δ的源项进行积分t,
    方程式(28)
    公式(29)
    我们在这里注意到,引力源项的时间离散化不同于Enzo(Bryan et al。1995),其中Sn+1/2页g通过外推法从tntn1。如果,在时间步长结束时,(ρE)n+1− (1/2)ρn+1(美国n+1)2>10个4E)n+1,我们重新定义(ρe)n+1= (ρE)n+1− (1/2)ρn+1(美国n+1)2; 否则我们使用en+1当需要计算压力或温度,并重新定义(ρE)n+1= (ρe)n+1+(1/2)ρn+1(美国n+1)2.

这就结束了单级算法的描述。

4.AMR

块结构AMR是由Berger&Oliger提出的(1984); Berger&Colella已经证明它在气体动力学方面是非常成功的(1989)在二维和贝尔等人(1994)在三维空间中,有着悠久的成功应用于各种领域的历史。纽约的AMR方法使用矩形网格的嵌套层次结构,在空间中网格的细化是层间的两倍,而层间时间的优化则由指定的子循环算法决定。

4.1.网格层次

纽约的网格层次结构由从最粗到最粗的不同层次的细化组成(ℓ = 0)至最佳(ℓ = ℓ优美的). 允许的最大细化级别数,最大值,在计算开始(或任何重新启动)时指定。在计算中的任何给定时间,可能有少于最大值层次结构中的级别,即。,优美的只要优美的≤ ℓ最大值. 每个层都由给定分辨率的非重叠矩形网格的并集表示。每个网格在每个坐标方向上由偶数个“有效”单元组成;每个坐标方向上的单元大小相同,但栅格在每个方向上的单元数可能不同。每个网格也有“重影单元”,它们在每个网格的低端和高端的每个坐标方向上扩展相同数量的单元。这些单元格用于临时保存用于更新“有效”单元格中值的数据;当在时间上进行亚循环时,重影单元区域必须足够大,使粒子处于水平位置不能离开级别的联合在一个关卡的过程中有效的和幽灵的细胞ℓ − 1时间步长。格线是正确嵌套的,从某种意义上说,格线在标高处的并集ℓ + 1包含在标高处的轴网联合中ℓ. 此外,从严格意义上说,控制水平网格足够大,至少可以保证边界n适当的水平每层周围的细胞宽度ℓ + 1格,其中n适当的由用户指定。没有严格的亲子关系;单一级别ℓ + 1网格可以覆盖多个级别的部分网格。

我们初始化网格层次结构,并按照Bell等人概述的过程重新网格化(1994). 在标高处定义轴网的步骤ℓ + 1,我们首先在级别标记单元格满足用户指定的条件。典型的标记标准,如局部超敏,可以在运行时选择;也可以使用专门的细化准则,它可以基于任何变量或变量组合,这些变量可以从流体状态或粒子信息中构造出来。标记的单元在标高处分组为矩形网格使用Berger&Rigoutsos中给出的聚类算法(1991). 然后对这些矩形面片进行优化,以在标高处形成栅格ℓ + 1.根据用户指定的情况,将大的补丁分成更小的补丁,分发给多个处理器最大网格大小参数。

一个粒子,,被称为“活”在水平中频电平是在该级别上栅格的并集包含粒子位置的最佳级别,. 据说粒子生活在n第三层网格如果粒子处于水平面以及n第三层网格包含. 每一个粒子,连同它的位置、质量和速度,都会存储有关它所在位置的当前信息,特别是级别、栅格编号和单元(,j,k). 每当粒子移动时,将在“重新分布”过程中为每个粒子重新计算此信息。在PM方案中,粒子的大小被设置为粒子当前所在级别的网格间距。

每开始k水平时间点,在哪里k修补程序由用户在所有级别定义-在所有级别运行ℓ + 1及以上,如果ℓ < ℓ最大值. 在以前被精细网格覆盖的区域,数据只是从旧网格复制到新网格;在新细化的区域中,数据是从底层较粗的网格中插值出来的。插值过程基于最近邻在每个粗单元内构造分段线性轮廓;对这些轮廓进行限制以避免引入任何新的最大值或最小值,然后将精细网格值定义为细网格单元中心的三线轮廓值,从而确保所有插值量的守恒。所有粒子处于水平通过最大值每当网格层次结构在级别上更改时,都会重新分布ℓ + 1及以上。

4.2.子循环选项

Nyx支持四种不同的选项来关联Δt,水平上的时间步长ℓ > 0,至Δt0,最粗层次的时间步。第一种是非次循环选择;在这种情况下,所有级别都以相同的Δ提前t,即Δt= Δt0对所有人ℓ. 第二个是“标准”子循环选项,其中Δt在所有水平上都是恒定的,即Δt=rΔt0对于空间细化比,r,级别之间。这些都是多级算法的标准选项。在第三种情况下,用户在运行时指定子循环模式;例如,用户可以指定级别0和级别1之间的子循环,但不指定级别1、2和3之间的子循环。然后Δt= Δt2= Δt1=rΔt0. 最后一个选项是“最优子循环”,其中最优子循环算法在每一个粗略的时间步长或其他指定的时间间隔决定子循环模式应该是什么。例如,如果每个级别的时间步长受Δ的约束t,强制条件的限制在一个时间步内不改变超过1%,那么不进行子循环是在多级层次上推进解的最有效方法。如果二级暗物质粒子的移动速度比三级粒子快两倍tdm公司是在这些水平上选择时间步长的限制条件,用相同的Δ推进2级和3级可能是最有效的t,但相对于级别0对级别1进行子循环,而相对于级别1对级别2进行子循环。这些选择是在运行期间由最优子循环算法进行的;E中描述了计算最优子循环模式的算法。Van Andel等人(2013年,准备中)。

除上述示例外,在大规模并行计算中,可能会出现额外的考虑因素,例如在特定级别上跨级别而不仅仅是在网格上并行化的能力,这可以使更复杂的子循环模式比“标准”模式更高效。我们将对这些情况的进一步讨论推迟到未来的工作中,我们将在具体的宇宙学例子中展示计算效率的权衡。

4.3.多级算法

我们将一个完整的多级“推进”定义为将解决方案在0级推进一个0级时间步长所需的操作和在所有更精细的级别推进解决方案所需的操作的组合,0<ℓ ≤ ℓ优美的,同时。我们定义n作为在水平上采取的时间步数对于在最粗略的水平上采取的每个时间步长,即Δt0=nΔt所有级别,0<ℓ ≤ ℓ优美的. 如果用户未指定子循环,则n=1;如果用户指定标准子循环,则n=r. 如果用户定义了一个子循环模式,则该模式用于计算n在计算开始时;如果使用最优子循环,则n以指定的间隔重新计算。

多级推进从计算Δ开始t对所有人ℓ ≥ 0,首先需要计算最大可能时间步长Δt最大值 , ℓ在每个级别独立(仅包括该级别的粒子)。给定一个亚循环模式和Δt最大值 , ℓ在每个级别,时间步计算如下:

方程(30)

方程(31)

一旦我们计算出Δt对所有人先进的从时间tt+ Δt0,通过调用递归函数定义一个完整的多级“advance”,预付款(ℓ), 如下所述ℓ = 0.我们注意到,在没有子循环的情况下,整个多级推进发生在一个调用预付款(0);在标准次循环的情况下预付款(ℓ) 将为每个ℓ, 0≤ ℓ ≤ ℓ优美的. 中间子循环模式需要调用预付款(ℓ) 每个级别相对于水平的亚循环ℓ − 1.在下面的注释中,上标n指的是时间,tn,每次呼叫预付款从水平开始ℓ, 上标呢n+1表示时间tn+ Δt.

预付款(ℓ):

  • 第0步:定义作为最好的水平$\Delta t{\ell\u a}=\Delta t\ell$如果我们在两层之间进行亚循环以及ℓ + 1,那么= 这个电话预付款仅预付款级别ℓ. 然而,如果我们不进行次循环,我们可以在水平上使用更有效的多级推进通过.
  • 第一步:计算phgr公司ngn在水平面ℓ → ℓ优美的使用ρnb和ρndm公司,式中ρndm公司是根据n.该计算包括在多级网格层次上求解Poisson方程全部的更精细的层次;关于计算ρ的进一步细节dm公司下一节给出了多级层次的求解方法。如果phgr公司在步骤1中的多级求解过程中,此时已在较粗的级别上计算,则跳过此步骤。
  • 第二步:预付款美国通过Δt在水平面ℓ → ℓ.因为我们使用显式方法来推进气体动力学,所以每一级都是由Δ独立推进的t根据边界条件的需要,从较粗的层次使用Dirichlet边界条件。如果> ℓ, 然后是水平ℓ → ℓ已经提出,我们执行流体力学量的显式回流(参见,例如,Almgren等人。2010有关回流的更多详细信息)从级别向下到水平以确保保护。
  • 第三步:插值gn在水平面ℓ → ℓ然后从粒子的位置,向粒子的水平移动ℓ → ℓ通过Δt/2和粒子位置t.此步骤中的操作与单级别版本中的操作相同,不需要在网格之间或级别之间进行通信。虽然在这一步中,粒子可能从有效的单元移动到重影单元,或者从重影单元移动到另一个重影单元,但它仍然是一个级别粒子直到重新分布发生在水平完成后ℓ − 1时间步长。
  • 第四步:计算phgr公司n+1gn+1在水平面ℓ → ℓ使用ρn+1个,*b和ρn+1dm公司,式中ρn+1dm公司是根据n+1.与步骤1中的解算不同,此解算超出了级别ℓ → ℓ除非这是一个无子循环算法,否则它只是一个部分多级解;关于多级求解的进一步细节将在下一节中给出。
  • 第五步:插值gn+1在水平面ℓ → ℓ从网格到粒子位置,然后更新粒子的速度ℓ → ℓ增加一个Δt/二. 同样,这一步与单级算法相同。
  • 第6步:对的美国在水平面ℓ → ℓ以时间为中心的源项,并替换e通过E− (1/2)美国2视情况而定。
  • 第7步:如果< 优美的,预付款(+(一)$n{\ella+1}/n\ell$ 时代.
  • 第八步:执行所有最终同步

4.4.重力解算

在上述算法中,有时重力位的泊松方程,

方程(32)

在多个级别上同时求解,以及仅在单个级别上需要解决方案的时间。我们定义作为标准的七点有限差分近似2在水平面ℓ, 适当修改ℓ/(ℓ − 1) 接口ifℓ > 0;例如,见Almgren等人(1998)有关以单元格为中心的多级接口模具的详细信息。

什么时候ℓ = ℓ优美的在步骤1中,或ℓ = ℓ在步骤4中,我们解决

论平面网格的并仅限;这被称为水平解算.

如果ℓ = 0则此级别的网格覆盖整个域,所需的唯一边界条件是域边界处的周期边界条件。如果ℓ > 0则水平的Dirichlet边界条件在ℓ/(ℓ − 1) 从数据层接口ℓ − 1.价值phgr公司在水平面ℓ − 假设1个生活在水平的细胞中心ℓ − 1个细胞;这些值在ℓ/(ℓ − 1) 提供具有间距Δ的边界值的接口对于级别网格。修改后的模板消除了对垂直于界面的插值的需要。利用几何多重网格技术,特别是具有红黑高斯-赛德尔松弛的V-圈,求解得到的线性系统。我们注意到ℓ > 0,phgr公司在水平面以及ℓ − 1满足Dirichlet但不满足Neumann匹配条件。这种不匹配在步骤8中描述的椭圆同步中被纠正。

当需要多级求解时,我们定义补偿ℓ,作为复合网格近似2在水平面上通过,并定义复合解算作为解决问题的过程

在水平面上通过. 复合解的解满足

在水平面,但满足

ℓ ≤ ℓ' <只在每个级别的区域被更细的网格覆盖或靠近更细网格区域的边界。在一个级别的区域ℓ' 水平覆盖的网格ℓ' + 1个网格$\phi^{\ell^\prime}$定义为$\phi^{{\ell^\prime}+1}$; 水平ℓ' 紧靠水平联合边界的单元ℓ' + 1网格,如Almgren等人所述,是一种改进的界面运算符(1998)被使用。用带V-圈的几何多重网格法求解该线性系统;这里的水平通过作为V循环中最好的等级。在层上进行复合解算后通过,phgr公司满足水平网格间界面的Dirichlet和Neumann匹配条件通过; 然而,Neumann匹配条件仍然存在不匹配phgr公司在水平面以及ℓ − 1如果ℓ > 0.这种不匹配在步骤8中描述的椭圆同步中被纠正。

4.5.云计算方案

Nyx使用云在细胞方案来定义每个暗物质粒子的质量对ρ的贡献dm公司; 粒子对ρ的贡献dm公司只在与边Δ的立方体相交的单元上以粒子为中心。在单水平算法中,所有粒子都生活在一个水平上,并对ρ有贡献dm公司在那个水平。如果一个级别上有多个网格,并且一个粒子生活在与另一个网格边界相邻的单元中,则该粒子的一部分贡献将贡献给它所在的网格,其余的将贡献给相邻的网格。粒子对相邻网格中单元的贡献分两个阶段计算:首先将贡献添加到粒子所在网格的重影单元,然后将重影单元中的值与相邻网格的有效单元中的值相加。泊松解的右侧只看到ρ的值dm公司在每个网格的有效单元格上。

在多级算法中,云计算过程变得更加复杂。在一个水平层级求解或层级求解ℓ →,一个水平的粒子住在附近ℓ/(ℓ − 1) 界面可能会在界面上“损失”一部分质量。另外,在一个级别求解时ℓ, 或者是多层次的解决方案ℓ →这是必要的,包括来自暗物质粒子在其他水平的贡献,即使这些粒子不在粗-细边界附近。通过创建和使用事实上的粒子。

大多数粒子在水平面上的质量ℓ' < ℓ 在Dirichlet边界条件下phgr公司phgr公司ℓ − 1ℓ/(ℓ − 1) 接口。然而,水平ℓ − 1个粒子ℓ/(ℓ − 1) 在水平面上沉积部分物质的界面ℓ, 或者实际上已经从水平移动了ℓ − 1至水平在一个水平ℓ − 1时间步,未正确考虑边界条件。这些效果是通过“重影粒子”(ghost particles)包含的,这些粒子是级别上的副本水平ℓ − 1个粒子。尽管这些幽灵粒子在水平面上活动ℓ, 它们保留了尺寸Δℓ − 1粒子的拷贝。

粒子在水平面上的质量ℓ' >在一个水平水平求解或水平ℓ →多级求解,包含在级别通过创建虚拟粒子。这些都是水平的居住在水平面上的粒子的表示ℓ' >. 在纽约市,可以选择将精细粒子与水平面一一对应虚粒子;还有一个选项可以将这些精细级别的粒子聚合成更小、更大质量的虚拟粒子. 在这两种表示中,较细粒子的总质量是守恒的。这些虚拟粒子是在级别上创建的在一个水平的开始时间步进,随真实水平移动粒子;这是必要的,因为ℓ' 粒子将不会移动,直到第二个泊松解与顶层已发生。虚粒子的质量对ρ有贡献dm公司在水平面遵循与level后面相同的cloud-in-cell过程粒子,但虚拟粒子的大小为Δ+1(即使它们是水平面上粒子的副本ℓ' ℓ' >+(一)。

5.软件设计与并行性能

5.1.概述

Nyx是在BoxLib(BoxLib)中实现的2011)framework,一个混合C++/Fortran90软件系统,为开发并行块结构AMR应用程序提供支持。基本并行化策略对基于MPI和OpenMP的多核体系结构使用分层编程方法。在纯MPI实例化中,每个级别至少有一个网格分布到每个核心,每个核心只使用MPI与其他每个核心通信。在混合方法中,每个插座上都有n所有访问同一内存的内核,我们可以用更少的、更大的每个套接字的网格,与这些网格相关的工作分布在n使用OpenMP的核心。

在BoxLib中,内存管理、流控制、并行通信和I/O都用程序的C++部分表示。计算中数值密集的部分,包括多重网格求解器,在Fortran90中处理。该软件支持一级数据的两种数据分配方案,以及一种动态交换方案,该方案根据一级网格的数量和处理器的数量来决定使用哪种方法。第一种方案基于Crutchfield中描述的启发式背包算法(1991)在Rendleman等人(2000). 第二种是基于Morton序空间填充曲线的使用。

5.2.颗粒和颗粒网格操作

BoxLib还支持并行粒子和PM操作。粒子根据其位置分布到处理器或节点;与粒子有关的信息在位置仅存在于拥有包含. 粒子操作是在粒子迭代器的帮助下执行的,粒子迭代器在粒子上循环;每个MPI进程都有一个迭代器,每个迭代器循环与该进程相关联的所有粒子。

因此粒子与网格数据的相互作用,如ρ的计算dm公司,或插值g从网格到,不需要处理器或节点之间的通信;它们与粒子数成线性关系,并且表现出几乎完美的弱标度。同样,使用插值重力加速度更新粒子速度和使用粒子速度更新粒子位置都不需要通信。粒子数据(位置、质量和速度)的唯一节点间通信发生在粒子穿过粗-细或细-细网格边界时;在这种情况下,重新分布过程将重新计算粒子级别、栅格编号和单元,粒子信息将发送到粒子现在所在的节点或核心。不过,我们注意到,当来自粒子的信息在网格层次结构上表示时,例如当粒子以泊松解算中使用的密度场的形式将其质量“沉积”在网格结构上时,如果粒子靠近其所在网格的边界,则基于网格的数据可能需要在节点之间传输,并将其部分质量沉积到同一水平或不同水平的不同网格中。

展望未来,我们还注意到,BoxLib中的粒子数据结构是用C++模板以足够通用的形式编写的,这些模板可以很容易地用来表示其他类型的粒子(例如那些相关的恒星形成机制),这些粒子除了位置、质量和速度之外,还携带着不同数量的属性。

5.3.并行I/O和可视化

随着模拟越来越大,高效地读写大型数据集的需求变得越来越重要。用于检查点和分析的数据由Nyx以自描述格式写入,该格式包含每个写入时间步的目录。检查点目录包含从该时间步重新启动计算所需的所有数据。Plotfile目录包含用于后处理、可视化和分析的数据,可以使用阿姆维斯,一个在LBNL开发的定制可视化软件包,用于在AMR网格上可视化数据,请访问(访问用户手册2005)或yt(Turk等人。2011). 在每个checkpoint或plotfile目录中,每个AMR级别都有一个ASCII头文件和子目录。标题描述了AMR层次结构,包括级别的数量、每个级别的网格框、问题大小、级别之间的优化比率、步骤时间等。每个级别目录中都是该级别数据的多个文件。检查点和打印文件目录按用户指定的间隔写入。

对于输出,在运行时为每个级别指定特定数量的数据文件。数据文件通常包含来自多个处理器的数据,因此每个处理器将其关联网格中的数据写入一个文件,然后另一个处理器可以将其关联网格中的数据写入该文件。指定的I/O处理器写入头文件,并协调允许哪些处理器写入哪些文件以及何时写入。处理器之间的唯一通信是在处理器可以开始写入时发信号和交换报头信息。虽然即使在一次运行期间,I/O性能也可能不稳定,但NERSC料斗机器上的并行I/O的最新计时,理论峰值为35 GB s1结果表明,在由多个均匀大小的网格组成的单一级别上运行时,Nyx的I/O性能可以达到34gbs以上1.

重新启动计算可能会给高效读取数据带来一些困难。由于文件的数量通常不等于处理器的数量,因此在重新启动期间协调输入以有效地读取数据。每个数据文件一次只能由一个处理器打开。指定的I/O处理器创建一个数据库,用于将文件映射到处理器,协调读取队列,并交叉读取自己的数据。每个处理器从当前打开的文件中读取所需的所有数据。该代码试图始终保持输入流的数量等于文件的数量。

检查点和绘图文件可移植到与写入文件的计算机具有不同字节顺序和精度的计算机。如果需要,字节顺序和精度转换在读取数据时自动完成。

5.4.并行性能

在数字上1我们在NERSC的Hopper机器上使用MPI和OpenMP展示了Nyx代码的缩放行为。我们进行了一项弱标度研究,因此每一次试验都有一个128每个NUMA节点的网格,其中每个NUMA节点有6个核心。4

图1。

图1。Santa-cranyx-Hopper(Santa-crany6)弱标度问题的模拟。

标准图像高分辨率图像

该问题是使用圣巴巴拉问题的复制初始化的,将在下一节中描述。首先,原来的256粒子被读入,然后这些粒子根据域的大小在每个坐标方向上被复制。例如,运行1024次×1024×2048个细胞,域在-以及是的-方向和八次z-方向。一个128网格每NUMA节点,本次运行将使用1024个处理器。物理域用索引空间缩放,因此分辨率(Δ)在弱标度下,问题的性质没有改变,例如,当问题变大时,可能影响多网格V循环数的问题的特征也没有改变。

图中的时间安排1是用于均匀网格计算的每个时间步在算法的不同部分所花费的时间。“水力”计时表示推进水动力状态所花费的所有时间,不包括重力计算。重力的计算分为两部分;多重网格解算器的设置和初始化,以及实际的多重网格V循环本身。设置阶段包括使用“单元内云”方案将质量从粒子分配到网格。在本次运行中,每个多重网格解算需要7个V循环才能达到10的公差12. 我们注意到,对于这个平均每个单元有一个粒子的问题,移动和重新分布粒子对总时间的贡献可以忽略不计。

我们在这里看到,Nyx代码可以从48个处理器扩展到49152个处理器,在处理器增加1000倍的情况下,总时间增加了不到50%。CASTRO的流体动力核心与Nyx基本相同,在OLCF(BoxLib)的Jaguarpf机器上,它的处理器几乎可以达到211K2011),并且使用AMR进行子循环的开销很小,从8个处理器的大约5%到64K处理器的19%(Almgren等人。2010). MAESTRO是另一个基于BoxLib的代码,它使用了与Nyx相同的线性解算器框架,它在Jaguarpf(BoxLib)上展示了出色的扩展到96K核的能力2011).

6.验证

纽约的流体动力学积分器是通过将卡斯特罗的积分器扩展到宇宙膨胀方程而建立的;如Almgren等人所述,Nyx流体力学已使用相同的试验和相同的结果进行了验证(2010). 为了分离和研究更新粒子位置和速度的算法组件的行为,我们对一个简单的两粒子轨道进行了分辨率研究。在一个边上16个单位的无量纲立方域中,我们初始化了两个质量相等的粒子,每个粒子的距离为1个单位-从域中心的方向。为了将数值结果与精确的圆轨道解进行比较,将初始粒子速度指定为与理想圆轨道相对应的速度,并在区域边界处施加由两个单极子展开(每个粒子一个)之和构造的Dirichlet边界条件。进行了六次数值试验:三次试验的基础网格为64,以及0、1或2级优化;一个基本网格为128再加上0或1级的精细化,最后是一个统一的256案例。接下来是三个关键的观察结果。首先,在观察到的10个轨道上,粒子的轨道是稳定的,随着时间的推移,测量到的性质几乎没有变化。第二,每个粒子的轨道半径和动能随网格间距收敛;对于最粗糙的情况,轨道半径和粒子动能与其初始值的最大偏差分别不超过1.1%和2.4%,对于最好的情况分别不超过0.17%和0.36%。三是同一分辨率下的单级和多级模拟的轨道半径偏差小于半径的0.0002%;动能偏差小于初始动能的0.0004%。

我们还使用已知的分析解进行了几个简单的测试,例如MacLaurin球体和Zel'dovich煎饼,以进一步验证Nyx重力解算器的粒子动力学和预期二阶精度和收敛特性。最后,在本节中,我们使用宇宙数据ArXiv(Heitmann et al。2005)以及圣巴巴拉集群比较项目(Frenk等人。1999).

6.1.仅暗物质模拟

为了在真实的宇宙学环境中验证Nyx,我们首先将Nyx纯重力模拟结果与Heitmann等人描述的比较数据集进行比较(2005). 域是256 Mpch1在一边,宇宙学是Ω=0.314,ΩΛ=0.686,h=0.71,σ8=0.84。初始条件在z=50,应用Klypin和Holtzman(1997)传递函数和使用Zel'dovich(1970)近似值。模拟进化到256暗物质粒子,质量分辨率为1.227×1011. 两份比较文件(Heitmann等人。2005海特曼等人。2008)在10个不同的,广泛使用的宇宙学代码中,几个导出的量有很好的一致性。

在这里,我们将我们的结果与Heitmann等人论文中提出的Gadget-2模拟结果进行比较;我们让感兴趣的读者参考这些论文,看看Gadget-2模拟与其他代码的比较。我们给出了三种不同的Nyx模拟的结果:在256和1024,以及一个256基本网格和两个细化级别,每一级都是2倍。在数字上2,我们显示了相关函数(上面板)和Gadget-2结果的比率(下面板),以便进行更仔细的检查。我们看到Gadget-2和Nyx结果之间有很强的匹配,并且观察到随着分辨率的增加Nyx代码向Gadget-2结果的收敛。特别是,我们注意到有效的1024个使用AMR获得的结果与使用uniform1024获得的结果非常接近网格。

图2。

图2。两个暗物质的256点关联函数h1框(上面板),以及相对于Gadget-2结果的比率(下面板)。我们展示了两个均匀网格模拟和一个AMR模拟的Nyx结果基本网格和两个细化级别,每一级都是2倍。

标准图像高分辨率图像

接下来,我们将研究光环的质量和空间分布,这两个重要的统计指标在整个宇宙学中以不同的方式使用。为了生成光环目录,我们对所有跑步都使用相同的光环查找器,它查找朋友的朋友光环(FOF;Davis等人。1985)使用链接长度b=0.2。为了关注代码之间的差异,我们考虑了只有10个粒子的光晕,但在实际应用中,我们只考虑至少有数百个粒子的晕,以避免由于FOF晕的有限采样而导致的大误差。在数字上,我们展示了质量函数的结果,并且我们确认了高分辨率Nyx运行和Gadget-2之间的良好一致性(2005),Heitmann等人(2005),以及Lukić等人(2007),常用的细化策略在低质量端抑制晕-质量函数。这是因为小晕形成的时间很早,并且贯穿整个模拟领域;捕捉它们需要如此早期和如此广泛的细化,基于块的细化(如果不是一般的AMR)几乎没有比固定网格模拟有任何优势。此处显示的Nyx AMR结果与Lukić等人给出的标准完全一致(2007),和类似于其他AMR代码(Heitmann等人。2008).

图3。

图3。光环质量函数也表现出良好的收敛性,随着力的分辨率增加,Gadget-2的结果,运行在更高的分辨率。实线是(Sheth&Tormen1999)适合,这里显示的纯粹是为了引导眼睛;与其他图中一样,对Gadget取比率。

标准图像高分辨率图像

在数字上4我们给出了晕的相关函数。我们不是把所有光环放在一起看,而是把它们分成三个大箱子。对于最小的光环,我们在低分辨率下看到了偏移,但即使这样也会很快收敛,尽管1024个光环run的小质量光环仍然比使用力分辨率小几倍的小工具run少(图). 最后,在图中5我们研究晕的内部结构。我们观察到AMR的运行和统一的1024非常一致运行,确认AMR模拟中已解析晕的结构与固定网格运行中的结构匹配的预期。总的来说,Nyx和Gadget-2的结果符合预期,并且在这个分辨率上与其他AMR宇宙学代码相当,如Heitmann等人的论文所示。

图4。

图4。暗物质晕的两点相关函数,分为三个质量群。我们看到这两个代码之间有很好的一致性,而且即使是对有几十个粒子采样的晕,也能快速收敛。我们强调,在实际应用中,10-100个粒子范围内的光晕不会被认真考虑,但我们在这里展示它们只是为了检查代码之间的差异。

标准图像高分辨率图像
图5。

图5。Gadget-2中三个最大质量物体的密度分布图。这些线条是最适合Gadget-2运行的NFW配置文件。我们可以看到,来自AMR模拟的轮廓与uniform1024的轮廓非常吻合模拟。

标准图像高分辨率图像

6.2.圣巴巴拉集群模拟

圣巴巴拉星团已经成为现实宇宙环境中星系团形成的标准代码比较测试问题。为了进行有意义的比较,所有的代码都从同一组初始条件开始,并将气体和暗物质演化为红移z=0假设理想气体状态方程没有加热或冷却(“绝热流体”)。初始条件受到约束,使得一个簇(3σ稀有)在盒的中心部分形成。模拟域一侧为64mpc,SCDM宇宙学:Ω=Ω=1,Ωb=0.1,h=0.5,σ8=0.9,以及ns=1.暗物质密度和气体温度的3D快照z=0如图所示6. 对于Nyx模拟,我们使用Heitmann等人的精确初始条件(2005),并在开始模拟z=63与256粒子。圣巴巴拉测试问题的代码比较研究的详细结果发表在Frenk等人的文章中(1999); 我们在这里介绍的不是全部,而是该论文中最重要的诊断,用于将Nyx与其他可用代码进行比较。我们让感兴趣的读者参考那篇文章,了解不同代码的仿真结果的全部细节。

图6。

图6。圣巴巴拉集群模拟。我们展示了三个正交的切取框中心,显示了蓝白颜色的暗物质密度。红色叠加的是10的等温面7K气体在整个模拟领域显示出丰富的星系际气体结构。

标准图像高分辨率图像

我们首先检查全局簇的属性,在内部计算r200关于红移的临界密度z=0。我们发现团簇的总质量为1.15×1015; Frenk等人(1999)报告平均值和1σ偏差为(1.1±0.05)×1015所有代码。从纽约运行的星团半径为2.7 Mpc;Frenk等人(1999)报告2.7±0.04。NFW浓度为7.1;原始文件给出了7.5作为浓度应该是多少的粗略指导。

在数字上7我们展示了几个量的径向分布的比较:暗物质密度,气体密度,压力和熵。由于比较项目的原始数据不公开,因此我们仅与平均值以及根据Kravtsov等人的数据估算的Enzo和ART数据进行比较(2002). 在初始条件下,我们也没有,和Frenk等的初始条件完全相同(1999)或Kravtsov等人(2002)我们的红移和恩佐的红移不同(z=30),所以我们不期望完全一致。

图7。

图7。左图:暗物质(上)和气体(下)的密度;右图:圣巴巴拉星系团的压力(上)和熵(下)分布图。除了纽约时报,我们还展示了原始Frenk等人的平均值。1999论文,以及两个AMR代码的结果,Enzo和ART。

标准图像高分辨率图像

在原始代码比较中,簇中心的定义由每个模拟器自行决定,这里我们采用引力势最小来定义簇中心。如Frenk等人(1999),我们在15个球壳中径向存储量,在对数半径上均匀分布,覆盖从10 kpc到10 Mpc的三个数量级。报告的压力为p= ρ气体T,熵的定义为s=单位 (T2/3页气体). 我们发现Nyx和其他AMR编码有很好的一致性,包括在簇中心的熵平坦。与网格方法形成鲜明对比的是,平滑粒子流体力学(SPH)代码发现中心熵继续向更小的半径上升。Mitchell等人(2009)对这一问题进行了详细的研究,发现两者之间的差异与网格代码中的网格尺寸或SPH中的粒子数无关。相反,这种差异在本质上是基本的,并且被证明是由于抑制了SPH中的涡流和流体不稳定性(另见亚伯2011总的来说,我们发现Nyx与现有的AMR规范有很好的一致性。

7.结论和今后的工作

我们提出了一个新的N-用于大规模宇宙模拟的天体和气体动力学代码。Nyx被设计成高效利用数万个处理器;代码的计时(多达50000个处理器)显示出优秀的弱伸缩行为。用绝热流体力学方法验证了Nyx在纯暗物质运行和暗物质中的有效性。未来的论文将更详细地介绍源项在Nyx中的实现,并将根据不同应用中提高保真度的需要,提供模拟结果,其中包括气体的不同加热和冷却机制。纽约大学正在进行的科学研究包括对Lyα森林和星系团的研究。此外,我们计划扩展Nyx以允许模拟替代宇宙模型到∧CDM,最有趣的是动态暗能量和爱因斯坦引力的修正。现有的网格结构和多网格Poisson解算器可以直接将现有的能力扩展到求解非线性椭圆型方程的迭代方法。

这项工作得到了伯克利实验室的实验室指导研发资金(PI:Peter Nugent),由美国能源部科学办公室主任根据合同号DE-AC02-05CH11231提供。通过SciDAC FASTMath Institute为BoxLib的改进提供了额外的支持,以支持Nyx代码和其他代码,由美国能源部高级科学计算研究办公室(和基础能源科学/生物和环境研究/高能物理/聚变能源科学/核物理办公室)资助的科学发现高级计算(SciDAC)计划资助。本文中的计算使用了国家能源研究科学计算中心(NERSC)的资源,该中心由美国能源部科学办公室(Office of Science of Science of Energy)提供支持,合同编号DE-AC02-05CH11231。

我们感谢彼得·纽金特和马丁·怀特进行了许多有益的讨论。

脚注

  • 我们注意到,对于小到中等数量的内核,纯MPI运行比使用混合MPI和OpenMP的运行要快得多。然而,为了这个标度研究的目的,我们只报告混合结果。

请稍候…引用正在加载。
10.1088/0004-637X/765/1/39