How can I implement the computation of the diffusion coefficient $D$ using periodic boundary conditions (PBC)?
I use molecular dynamics of a set of $nboby$ particles with positions $pos(3,nbody)$ in a box of length $length$. The implemetion of the PBC is
do k = 1, 3
pos(k,i) = modulo(pos(k,i),length)
end do
At now I'm using for $D$ the following code
do it=1,nstep
diff=0
do i = 1,nbody
pos2(:)=pos2(:)+pos(:,i)
diff=diff+dot_product(pos(:,i),pos(:,i))
end do
diff=diff/nbody-dot_product(pos2(:),pos2(:))/nbody**2
end do
diff=diff/nstep/6
which I think it corresponds to
$D=\lim_{t\to\infty}\dfrac{<x^2>-<x>^2}{6t}$
but I'm not very sure that the PBC are taking into account in the right way.
Can someone help me?
Thanks Matteo