如何利用数组搜索算法自动定位delta=0对应的非初始时刻时间点？

当然可以用数组搜索算法完全自动化这个过程！你的d_values和t_values是严格等长的，我们可以借助NumPy的向量化操作快速定位目标索引，直接提取对应的时间值，再也不用手动数步数啦。

核心思路 #

1. 处理浮点数精度：由于数值计算存在微小误差，直接判断delta == 0可能会漏掉实际接近0的点，所以用np.isclose设置一个极小的容差来判断是否接近0。
2. 排除初始时刻：添加t != 0的条件，过滤掉t=0的情况。
3. 批量提取结果：一次性筛选出所有符合条件的时间点，不用逐个手动查找。

修改后的完整代码 #

我已经把自动化搜索的逻辑整合到你的代码里了，放在绘图部分之后：

import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate

masses = [1, 1, 1]
r1, v1 = [0, 0], [-2*0.513938054919243, -2*0.304736003875733]
r2, v2 = [-1, 0], [0.513938054919243, 0.304736003875733]
r3, v3 = [1, 0], [0.513938054919243, 0.304736003875733]
u0 = np.concatenate([r1, v1, r2, v2, r3, v3])

def odesys(t, u):
    def force(a):
        return a / sum(a ** 2) ** 1.5
    r1, v1, r2, v2, r3, v3 = u.reshape([-1, 2])
    m1, m2, m3 = masses
    f12, f13, f23 = force(r1 - r2), force(r1 - r3), force(r2 - r3)
    a1, a2, a3 = -m2 * f12 - m3 * f13, m1 * f12 - m3 * f23, m1 * f13 + m2 * f23
    return np.concatenate([v1, a1, v2, a2, v3, a3])

# collect data
t_values = []
u_values = []
par_1_pos = []
d_values = []

# Time start, step, and finish point
t0, tf, t_step = 0, 18, 0.0001
nsteps = int((tf - t0) / t_step)
solution = integrate.RK45(odesys, t0, u0, tf, max_step=t_step)

# The loop for running the Runge-Kutta method over some time period.
u_values.append(solution.y)
t_values.append(t0)
par_1_pos.append(((solution.y[0] - u0[0])**2 + (solution.y[1] - u0[1])**2)**0.5)
d_values.append(((solution.y[0] - u0[0])**2 + (solution.y[1] - u0[1])**2)**0.5 + ((solution.y[4] - u0[4])**2 + (solution.y[5] - u0[5])**2)**0.5 + ((solution.y[8] - u0[8])**2 + (solution.y[9] - u0[9])**2)**0.5 + ((solution.y[2] - u0[2])**2 + (solution.y[3] - u0[3])**2)**0.5 + ((solution.y[6] - u0[6])**2 + (solution.y[7] - u0[7])**2)**0.5 + ((solution.y[10] - u0[10])**2 + (solution.y[11] - u0[11])**2)**0.5)

for step in range(nsteps):
    solution.step()
    u_values.append(solution.y)
    t_values.append(solution.t)
    par_1_pos.append(((solution.y[0] - u0[0])**2 + (solution.y[1] - u0[1])**2)**0.5)
    d_values.append(((solution.y[0] - u0[0])**2 + (solution.y[1] - u0[1])**2)**0.5 + ((solution.y[4] - u0[4])**2 + (solution.y[5] - u0[5])**2)**0.5 + ((solution.y[8] - u0[8])**2 + (solution.y[9] - u0[9])**2)**0.5 + ((solution.y[2] - u0[2])**2 + (solution.y[3] - u0[3])**2)**0.5 + ((solution.y[6] - u0[6])**2 + (solution.y[7] - u0[7])**2)**0.5 + ((solution.y[10] - u0[10])**2 + (solution.y[11] - u0[11])**2)**0.5)
    # break loop after modelling is finished
    if solution.status == 'finished':
        break

# Plotting of the individual particles
u = np.asarray(u_values).T

# Plot for The trajectory of the three bodies over the time period
plt.plot(u[0], u[1], '-o', lw=1, ms=3, label="body 1")
plt.plot(u[4], u[5], '-x', lw=1, ms=3, label="body 2")
plt.plot(u[8], u[9], '-s', lw=1, ms=3, label="body 3")
plt.title('Trajectories of the three bodies')
plt.xlabel('X Position')
plt.ylabel('Y Position')
plt.legend()
plt.grid()
plt.show()
plt.close()

# Plot for d(delta_t) values
plt.plot(t_values, d_values)
plt.title('Delta number for the three bodies')
plt.xlabel('Time (s)')
plt.ylabel('Delta')
plt.grid()
plt.show()
plt.close()

# Plot of distance between P1 and IC
plt.plot(t_values, par_1_pos)
plt.title('Plot of distance between P1 and IC')
plt.xlabel('Time (s)')
plt.ylabel('Distance from origin')
plt.grid()
plt.show()
plt.close()

# -------------------------- 自动化搜索部分 --------------------------
# 转换为NumPy数组，方便操作
t_array = np.array(t_values)
d_array = np.array(d_values)

# 设置容差，处理浮点数精度问题，这里用1e-8可以根据你的计算精度调整
tolerance = 1e-8
# 筛选条件：delta接近0，且排除t=0的时刻
mask = np.isclose(d_array, 0, atol=tolerance) & (t_array != 0)

# 提取符合条件的时间值
result_times = t_array[mask]

# 输出结果
if len(result_times) > 0:
    print("找到delta=0（排除t=0）对应的时间点：")
    for idx, time in enumerate(result_times, 1):
        print(f"第{idx}个时间点：{time:.6f} s")
else:
    print("未找到符合条件的时间点（排除t=0）")

关键细节说明 #

• 浮点数容差：np.isclose的atol参数是绝对容差，如果你的delta计算精度更高，可以适当调小这个值（比如1e-10），避免误判。
• 批量结果：如果你的模拟中有多个delta=0的时刻（除了t=0），这段代码会全部输出，不用重复手动查找。
• 效率：NumPy的向量化操作比手动循环快得多，即使你的数组有上百万个元素，也能瞬间完成搜索。

内容的提问来源于stack exchange，提问作者Ross Alan Slater