When both using fix evaporate and compute msd in the same group

[notsweetie]

Hi! Lammps experts!

i am trying to track an Ar group diffusion charateristic in an open system.

when i am trying to apply this code

 ...
fix deletion ar evaporate 100 100000 dele 12345  #apply deletion 
...
compute msd_ar ar msd
...

The simulation runs without error. However, I’m trying to better understand how LAMMPS handles atoms that are deleted after the compute msd command is applied.

As I understand, the reference positions for the MSD are set at the time compute msd_ar is defined. Then, as the simulation progresses, some atoms are deleted due to fix evaporate. My question is:

How are the deleted atoms treated in the MSD calculation? Which of the following interpretations is more accurate?

1. (x, y, z) → (x, y, z):
   The last position of a deleted atom is retained, and since it is no longer in the compute group, its contribution to MSD becomes zero.

2. (x, y, z) → (0, 0, 0):
   Deleted atoms are treated as having position (0, 0, 0), so their MSD becomes artificially large due to dx² + dy² + dz².

The LAMMPS documentation for compute msd states:

“The value of the displacement will be 0.0 for atoms not in the specified compute group.”

But I’m a bit confused — does this mean the MSD of deleted atoms is excluded (i.e., zeroed out), or could there still be residual contributions?

Any clarification on how compute msd handles group membership changes due to deletion would be greatly appreciated!

Thanks in advance!

[akohlmey]

@notsweetie

The simulation runs without error.

It did not error out, but it should have. Using compute msd requires that the number of atoms in the compute group does not change. For example, it refuses to run with a dynamic group.

@notsweetie

How are the deleted atoms treated in the MSD calculation? Which of the following interpretations is more accurate?

1. (x, y, z) → (x, y, z):
   The last position of a deleted atom is retained, and since it is no longer in the compute group, its contribution to MSD becomes zero.
2. (x, y, z) → (0, 0, 0):
   Deleted atoms are treated as having position (0, 0, 0), so their MSD becomes artificially large due to dx² + dy² + dz².

Neither of these two interpretations are correct. What happens is that deleted atoms are no longer in the group and thus their contribution to the msd is skipped, but since the result is divided by the total number of atoms in the group at the beginning of the run, the result is incorrectly divided by a too large number.

We will correct the code to double check the size of the group and error out if it changes. Changing the number of monitored atoms is ill defined behavior for the MSD algorithm.

[Germain]

Hi @notsweetie,

To complete on @akohlmey, this is typically the type of computation that is better treated in post-processing than on-the-fly for two reasons.

1. The number of atoms in the evaporate group is unpredictable during the simulation. This is why Axel says the code should throw an error. compute MSD takes the initial position of a set of atoms and takes the average distance between their current position and their initial position at each invocation. Atoms not in the group do not contribute to the compute value.
2. The reference in time is limited to the initial initialization of the compute. If you post-process the data yourself, you can compute individual diffusion for each atom coming in the simulation until they leave, and even use FFT methods that computes running correlation along the trajectory.

[Germain]

I think the cited phrasing of the documentation is a mis-phrasing copied from the compute displace/atom command doc page that does not apply here. I shall push something in your PR of small corrections.