E_pe_tot = 0 do i = 1, n_atoms-1 do j = i+1, n_atoms if (mol(i) /= mol(j)) then dx = xu(i) - xu(j) dy = yu(i) - yu(j) do while (dx > L/2.d0) dx = dx - L end do do while (dx < -L/2.d0) dx = dx + L end do do while (dy > L/2.d0) dy = dy - L end do do while (dy < -L/2.d0) dy = dy + L end do dx_min = dx dy_min = dy dis = sqrt(dx_min**2+dy_min**2) E_pe = 0 delta_d = 0 if (dis < radius(i)+radius(j)) then delta_d = (radius(i)+radius(j)-dis) r_eff = sqrt(radius(i)*radius(j)/(radius(i)+radius(j))) E_pe = 0.4d0*r_eff*Kn*(delta_d ** 2.5) end if E_pe_tot = E_pe_tot + E_pe end if end do end do

这是一个Fortran代码片段，可以用于计算分子系统中原子之间的相互作用能量。该代码使用了一个双重循环来遍历所有原子对，并计算它们之间的相互作用。其中，mol(i)代表第i个原子所属的分子，xu(i)和yu(i)代表第i个原子的x和y坐标，L代表系统的尺寸，radius(i)代表第i个原子的半径，Kn代表弹性系数。在计算过程中，将两个原子之间的距离限制在半径之内，使用了一个6-12势函数来计算相互作用能量。总能量E_pe_tot累加了所有原子对之间的相互作用能量E_pe。

需要注意的是，该代码并没有使用空间分割法，而是采用了暴力计算的方法。这种方法的缺点是计算量较大，计算时间较长，但是可以得到较为精确的结果。如果要使用空间分割法来优化计算效率，则需要对系统进行空间划分，并只计算处于同一或相邻立方体中的原子对之间的相互作用。