swiftsim(swift.dur.ac.uk)python的I/O例程。
swiftsimio的Python项目详细描述
雨燕
使用swift天体物理模拟代码(http://swift.dur.ac.uk)
广泛的。从SWIFT中读取数据的方式有很多种。
HDF5文件。从直接使用h5py
读取到使用复合
像yt
这样的系统;但是它们都不令人满意(例如
读取hdf5中的单元信息),或者对于大多数用例来说过于复杂。这个
(瘦)包装器提供了一个面向对象的api来读取(动态)数据
来自斯威夫特。
要求
这需要python3.6.0
或更高版本。我们将不遗余力地支持
下面是python版本。请更新您的系统。
python包
h5py
unyt
用法
示例用法如下所示,它绘制了密度-温度阶段 图表,密度和温度以cgs单位表示:
importswiftsimioassw# This loads all metadata but explicitly does _not_ read any particle data yetdata=sw.load("/path/to/swift/output")importmatplotlib.pyplotaspltdata.gas.density.convert_to_cgs()data.gas.temperature.convert_to_cgs()plt.loglog()plt.scatter(data.gas.density,data.gas.temperature,s=1)plt.xlabel(fr"Gas density $\left[{data.gas.density.units.latex_repr}\right]$")plt.ylabel(fr"Gas temperature $\left[{data.gas.temperature.units.latex_repr}\right]$")plt.tight_layout()plt.savefig("test_plot.png",dpi=300)
在上面,必须注意以下几点:
- 调用
load
函数时读取所有元数据。 - 只有密度和温度(对应于
PartType0/Temperature
)数据集被读入。 - 只有调用
convert_to_cgs
方法时,才会读入该数据。 convert_to_cgs
就地转换数据;即返回None
。- 当调用
plt.scatter
时,数据被缓存而不是重新读入。
写入数据集
为宇宙学代码编写有效的数据集可以是
困难的,特别是在考虑如何最好地使用单位时。斯威夫特使用
不同的内部单位集(在参数文件中指定)可以
不一定要和初始条件相同
在中指定。尽管如此,重要的是确保
初始条件彼此都是一致的。为了促进这一点,
我们使用unyt
数组。下面的示例生成随机放置的气体
密度均匀的粒子。
fromswiftsimioimportWriterfromswiftsimio.unitsimportcosmo_unitsimportunytimportnumpyasnp# Box is 100 Mpcboxsize=100*unyt.Mpc# Generate object. cosmo_units corresponds to default Gadget-oid units# of 10^10 Msun, Mpc, and km/sx=Writer(cosmo_units,boxsize)# 32^3 particles.n_p=32**3# Randomly spaced coordinates from 0, 100 Mpc in each directionx.gas.coordinates=np.random.rand(n_p,3)*(100*unyt.Mpc)# Random velocities from 0 to 1 km/sx.gas.velocities=np.random.rand(n_p,3)*(unyt.km/unyt.s)# Generate uniform masses as 10^6 solar masses for each particlex.gas.masses=np.ones(n_p,dtype=float)*(1e6*unyt.msun)# Generate internal energy corresponding to 10^4 Kx.gas.internal_energy=np.ones(n_p,dtype=float)*(1e4*unyt.kb*unyt.K)/(1e6*unyt.msun)# Generate initial guess for smoothing lengths based on MIPSx.gas.generate_smoothing_lengths(boxsize=boxsize,dimension=3)# If IDs are not present, this automatically generates x.write("test.hdf5")
然后,对生成的test.hdf5
运行h5glance
会产生:
test.hdf5
├Header
│ └5 attributes:
│ ├BoxSize: 100.0
│ ├Dimension: array [int64: 1]
│ ├Flag_Entropy_ICs: 0
│ ├NumPart_Total: array [int64: 6]
│ └NumPart_Total_HighWord: array [int64: 6]
├PartType0
│ ├Coordinates [float64: 32768 × 3]
│ ├InternalEnergy [float64: 32768]
│ ├Masses [float64: 32768]
│ ├ParticleIDs [float64: 32768]
│ ├SmoothingLength [float64: 32768]
│ └Velocities [float64: 32768 × 3]
└Units
└5 attributes:
├Unit current in cgs (U_I): array [float64: 1]
├Unit length in cgs (U_L): array [float64: 1]
├Unit mass in cgs (U_M): array [float64: 1]
├Unit temperature in cgs (U_T): array [float64: 1]
└Unit time in cgs (U_t): array [float64: 1]
注意您确实需要注意您选择的单元系统 not允许值超过2^31,即您需要确保 当将写入文件时,提供的值(带单位)是安全的 解释为(单精度)浮点数。唯一的例外 这是以双精度存储的坐标。
提前掩蔽
swift快照包含允许我们在空间上屏蔽
数据提前。swiftsimio
提供了许多对象来帮助
这个。请参阅下面的示例。
importswiftsimioassw# This creates and sets up the masking object.mask=sw.mask("/path/to/swift/snapshot")# This ahead-of-time creates a spatial mask based on the cell metadata.mask.constrain_spatial([[0.2*mask.metadata.boxsize[0],0.7*mask.metadata.boxsize[0]],None,None])# Now, just for fun, we also constrain the density between 0.4 g/cm^3 and 0.8. This reads in# the relevant data in the region, and tests it element-by-element.density_units=mask.units.mass/mask.units.length**3mask.constrain_mask("gas","density",0.4*density_units,0.8*density_units)# Now we can grab the actual data object. This includes the mask as a parameter.data=sw.load("/Users/josh/Documents/swiftsim-add-anarchy/examples/SodShock_3D/sodShock_0001.hdf5",mask=mask)
当这个数据对象的属性被访问时,透明的只有那些 属于遮罩区域(在密度和空间上)被读取。也就是说,如果我要求 粒子温度,它将接收包含粒子温度的数组 它位于[0.2,0.7]区域,密度在0.4到0.8g/cm^3之间。
用户定义的粒子类型
现在可以添加尚未添加的用户定义粒子类型
存在于swiftsimio
元数据中。你只需要指定三个
名称(见下文)以及在中提供的粒子数据集
SWIFT将被自动读取。
importswiftsimioasswimportswiftsimio.metadata.particleasswpfromswiftsimio.objectsimportcosmo_factor,aswp.particle_name_underscores[6]="extratype"swp.particle_name_class[6]="Extratype"swp.particle_name_text[6]="Extratype"data=sw.load("extra_test.hdf5",)
以前,有一个生成器函数api执行此任务;现在 已删除(最新版本为v0.8.0)。
图像创建
swiftsimio
提供了一些非常基本的可视化软件
用X和Y创建整个盒子沿Z轴的投影
代表最后的正方形。注意,它目前只支持square
盒。这些是用numba
加速的,我们不建议使用这些
除非你安装了程序。最后,您可以使用并行版本
通过传递参数
parallel=True
。
你可以做以下内容:
fromswiftsimioimportloadfromswiftsimio.visualisationimportproject_gasdata=load("my_snapshot_0000.hdf5")# This creates a grid that has units K / Mpc^2, and can be transformed like# any other unyt quantitytemperature_map=project_gas(data,resolution=1024,project="temperature")frommatplotlib.pyplotimportimsavefrommatplotlib.colorsimportLogNorm# There is huge dynamic range usually in temperature, so log normalize it before# saving.imsave("temperature_map.png",LogNorm()(temperature_map.value),cmap="twilight")
也可以以相当简单的方式创建视频:
fromswiftsimioimportloadfromswiftsimio.visualisationimportproject_gas_pixel_grid# project_gas_pixel_grid does not perform the unit calculation.fromp_tqdmimportp_mapdefload_and_make_image(number):filename=f"snapshot_{number:04d}.hdf5"data=load(filename)image=project_gas_pixel_grid(data,1024,"temperature")returnimage# Make frames in parallel (reading also parallel!)frames=p_map(load_and_make_image,range(0,120))#... You can do your funcAnimation stuff here now you have the frames.
也可以创建slice图,而不是投影。
fromswiftsimioimportloadfromswiftsimio.visualisationimportslice_gasdata=load("my_snapshot_0000.hdf5")# This creates a grid that has units K / Mpc^3, and can be transformed like# any other unyt quantity. The slice variable gives where we want to slice through# as a function of the box-size, so here we slice through the centre of the box.temperature_map=slice_gas(data,slice=0.5,resolution=1024,project="temperature")frommatplotlib.pyplotimportimsavefrommatplotlib.colorsimportLogNorm# There is huge dynamic range usually in temperature, so log normalize it before# saving.imsave("temperature_map.png",LogNorm()(temperature_map.value),cmap="twilight")