espresso icon indicating copy to clipboard operation
espresso copied to clipboard

Monte Carlo displacement move is broken

Open jngrad opened this issue 2 years ago • 2 comments

There are several issues with method RE.displacement_mc_move_for_particles_of_type():

  1. it is completely untested and unused in the tutorials and samples
  2. the docstring and user guide don't mention it treats particles as distinguishable (sampling without replacement), i.e. when moving 3 particles, you can never be in a situation where you move particle 1, 2, and then 1 again
  3. the code enters an infinite loop when passing argument particle_number_to_be_changed=N when the system contains less than N particles of the chosen type
  4. when the move is rejected, the original particle position is restored, but not the original particle velocity
  5. when moving N particles, if move j is rejected, the remaining N - j - 1 moves are still attempted, even though all N moves will necessarily be reverted at the end; this leaves N particles in their original positions but with different velocities

jngrad avatar Oct 07 '22 16:10 jngrad

Hi Jean-Noel

The function implements a multi-particle MC move. It moves 1...N distinct particles (by design) where N is defined by the calling function.

Good that you catched point 4.

In point 5: in an N particle move all N particles are moved and the energy prior to the move and after the move of N particles is computed. It is not N single particle moves in a row.

In low density systems a N-particle move can be computationally cheaper than N times a 1 particle moves because the energy is less frequently computed. In high density systems the acceptance rate drops the more particles are moved in one single MC move. In the N particle MC move, we could implement am early stopping mechanism if there is particle overlap (within the exclusion radius) for optimization. But that is a separate topic.

jonaslandsgesell avatar Oct 10 '22 17:10 jonaslandsgesell

Thanks for your remarks. Regarding point 5, the way I understand it, generate_new_particle_positions() will reject a move j if the new position is inside an exclusion zone and will set the flag particle_inside_exclusion_range_touched to true. It is therefore pointless to attempt the remaining N - j - 1 moves, because that flag skips the second energy calculation and sets the new energy to $+\infty$, right?

In the N particle MC move, we could implement am early stopping mechanism if there is particle overlap (within the exclusion radius) for optimization. But that is a separate topic.

The PR implements such an early stopping mechanism, since it's a low-hanging fruit and most of code around it had to be modified anyway.

The function implements a multi-particle MC move. It moves 1...N distinct particles (by design) where N is defined by the calling function.

This is now documented.

jngrad avatar Oct 10 '22 18:10 jngrad