python - 计算 1d 内单个粒子的均方位移

Question

我有一维粒子轨迹，j=[]以及 time=np.arange(0, 10 + dt, dt)哪里dt是时间步长。我已经根据this article计算了MSD .

我在谷歌和这里搜索了Python中的1d MSD，但没有找到任何合适的，因为我的Python知识非常初级。我编写了一段代码，它工作正常，没有任何错误，但我不确定它是否根据给定的文章代表相同的事情。这是我的代码，

j_i = np.array(j)
MSD=[]
diff_i=[]
tau_i=[]
for l in range(0,len(time)):
    tau=l*dt
    tau_i.append(tau)
    for i in range(0,(len(time)-l)):
        diff=(j_i[l+i]-j_i[i])**2

        diff_i.append(diff)

    MSD_j=np.sum(diff_i)/np.max(time)
    MSD.append(MSD_j) 

任何人都可以检查验证代码，如果有误请提出建议。

Answer 1

代码大部分是正确的，这是一个修改版本，其中:

• 我简化了一些表达式(例如range)
• 我直接使用 np.mean 修正了平均值，因为 MSD 是平方位移 [L^2]，而不是比率 [L^2]/[T]。

最终代码:

j_i    = np.array(j)
MSD    = []
diff_i = []
tau_i  = []

for l in range(len(time)):
    tau = l*dt
    tau_i.append(tau)

    for i in range(len(time)-l):
        diff = (j_i[l+i]-j_i[i])**2
        diff_i.append(diff)

    MSD_j = np.mean(diff_i)
    MSD.append(MSD_j)

编辑:我意识到我忘记提及它，因为我专注于代码，但是论文中用 <.> 表示的整体平均值应该，顾名思义，在几个上执行粒子，优先将每个粒子的初始位置与其经过一段时间 tau 后的新位置进行比较，而不是像使用时间运行平均值那样进行比较

编辑2:这里的代码展示了如何进行适当的整体平均以准确实现本文中的公式

js     = # an array of shape (N, M), with N the number of particles and
         # M the number of time points
MSD_i  = np.zeros((N, M))
taus   = []

for l in range(len(time)):
    taus.append(l*dt)  # store the values of tau

    # compute all squared displacements at current tau
    MSD_i[:, l] = np.square(js[:, 0] - js[:, l])

# then, compute the ensemble average for each tau (over the N particles)
MSD = np.mean(MSD_i, axis=0)

现在您可以绘制 MSD 与 taus 的关系图，而 Bob 是您的叔叔