我这里有一段代码,可以发现两个粒子何时彼此近距离接触。 5个粒子呢?
I have here a block of code that finds when two particles have a close encounter with each other. What about 5 particles?
这个代码的原因是为了找到一个行星系统有2个行星(粒子)彼此近距离相遇所需的时间。但是,我想解决具有 5 个粒子的行星系统的时间问题。这将打印在最后一行 print(sim.t)
.
http://rebound.readthedocs.io/en/latest/examples.html 这里有一些关于 REBOUND 用法的文档,您可以在其中找到对我所做的一些声明的最佳解释。
sim = setupSimulation()
sim.exit_min_distance = 0.01
Noutputs = 10000
year = 2.*np.pi
times = np.linspace(0.,10E+9.*year, Noutputs)
distances = np.zeros(Noutputs)
ps = sim.particles
try:
for i,time in enumerate(times):
sim.integrate(time)
dp = ps[1] - ps[2]
distances[i] = np.sqrt(dp.x*dp.x+dp.y*dp.y+dp.z*dp.z)
except rebound.Encounter as error:
print(error)
print(sim.t)
第 2 行定义了两个粒子之间的近距离接触(达到此值时模拟将停止)。
第11行取粒子坐标差,如果<=.01,则打印近距离接触
我在考虑 if-then 语句:
if ps[1] - ps[2] <= .01:
dp = ps[1] - ps[2]
else ps[2] - ps[3] <= .01:
dp = ps[2] - ps[3]
等等...
我想在我 运行 之前确保它能正常工作,因为模拟时间是 10e+9 年,并且只能 运行 在本地超级计算机上获得结果合理的时间。
我不明白为什么 if 语句方法无法检测任何两个粒子是否彼此足够接近。编辑:编写此示例就像您的 2 个粒子代码示例正在运行一样。
根据您是否需要检查非连续粒子,您需要生成每个组合并进行检查(如示例所示)。
看来您还需要不同的距离枚举策略。如果索引不必与时间计数器匹配,您可以增加一个单独的计数器(如示例所示)。
您还可以使用 itertools
生成 if 语句以避免拼写错误并提高清晰度:
import itertools
# other setup code
counter = 0
try:
for time in times:
sim.integrate(time)
for p0, p1 in itertools.combinations(range(1, 6)): # generates combinations of the particles numbered 1 thru 5
dp = ps[p0] - ps[p1]
distances[counter] = np.sqrt(dp.x*dp.x+dp.y*dp.y+dp.z*dp.z)
counter += 1
except rebound.Encounter as error:
print(error)
我再补充一点——为了确保您的代码有效,您可能应该为它编写某种测试套件。如果它解决了与从不解决会发生什么。
您可以创建一个对象来模拟 sim
及其方法,在您 知道 将满足或不满足退出条件的位置返回主体,而无需所有开销需要一台超级计算机。您可以将它们硬编码到一个小列表中。
这个代码的原因是为了找到一个行星系统有2个行星(粒子)彼此近距离相遇所需的时间。但是,我想解决具有 5 个粒子的行星系统的时间问题。这将打印在最后一行 print(sim.t)
.
http://rebound.readthedocs.io/en/latest/examples.html 这里有一些关于 REBOUND 用法的文档,您可以在其中找到对我所做的一些声明的最佳解释。
sim = setupSimulation()
sim.exit_min_distance = 0.01
Noutputs = 10000
year = 2.*np.pi
times = np.linspace(0.,10E+9.*year, Noutputs)
distances = np.zeros(Noutputs)
ps = sim.particles
try:
for i,time in enumerate(times):
sim.integrate(time)
dp = ps[1] - ps[2]
distances[i] = np.sqrt(dp.x*dp.x+dp.y*dp.y+dp.z*dp.z)
except rebound.Encounter as error:
print(error)
print(sim.t)
第 2 行定义了两个粒子之间的近距离接触(达到此值时模拟将停止)。
第11行取粒子坐标差,如果<=.01,则打印近距离接触
我在考虑 if-then 语句:
if ps[1] - ps[2] <= .01:
dp = ps[1] - ps[2]
else ps[2] - ps[3] <= .01:
dp = ps[2] - ps[3]
等等...
我想在我 运行 之前确保它能正常工作,因为模拟时间是 10e+9 年,并且只能 运行 在本地超级计算机上获得结果合理的时间。
我不明白为什么 if 语句方法无法检测任何两个粒子是否彼此足够接近。编辑:编写此示例就像您的 2 个粒子代码示例正在运行一样。
根据您是否需要检查非连续粒子,您需要生成每个组合并进行检查(如示例所示)。
看来您还需要不同的距离枚举策略。如果索引不必与时间计数器匹配,您可以增加一个单独的计数器(如示例所示)。
您还可以使用 itertools
生成 if 语句以避免拼写错误并提高清晰度:
import itertools
# other setup code
counter = 0
try:
for time in times:
sim.integrate(time)
for p0, p1 in itertools.combinations(range(1, 6)): # generates combinations of the particles numbered 1 thru 5
dp = ps[p0] - ps[p1]
distances[counter] = np.sqrt(dp.x*dp.x+dp.y*dp.y+dp.z*dp.z)
counter += 1
except rebound.Encounter as error:
print(error)
我再补充一点——为了确保您的代码有效,您可能应该为它编写某种测试套件。如果它解决了与从不解决会发生什么。
您可以创建一个对象来模拟 sim
及其方法,在您 知道 将满足或不满足退出条件的位置返回主体,而无需所有开销需要一台超级计算机。您可以将它们硬编码到一个小列表中。