How to keep a rigid body stationary

Please see the forum guidelines for suggestions about how to correctly quote input files.

That post will also explain that it is crucial to report which LAMMPS version you are using and what platform you are running on.

Since you define fix setforce after the minimize command, the minimization will move all atoms.
You have excluded the interaction between the “heavyatoms” group atoms, but the forces due to all other atoms still apply. Please also note that

You would likely want to use the “overlap” keyword here to avoid having close contacts.

This is superfluous and possibly counterproductive. You don’t need to use a rigid fix to immobilize atoms. Just just don’t include them in any time integration and they won’t move; simple as that.

These two commands would not be needed when dropping fix rigid.

Is your water model rigid (TIP3P, or SPC/E, or TIP4P, or similar)? Then there seems to be a “fix shake” command missing, otherwise, your timestep is too large for the movement of the bonds involving hydrogen atoms. Please note that with recent LAMMPS versions it is also possible to use fix shake during minimization (it will use an approximation for SHAKE with strong harmonic bonds).

I don’t see you adding or removing atoms. Why this command?

Do you want atoms to be removed from the simulation? If yes, that does not make much sense since your have periodic boundaries in all directions. Lost atoms would be a sign of some serious problems with the simulation and thus they should not be ignored.

Since you are simulating a slab geometry, you should consider using boundary p p f and kspace_modify slab 3.0 to avoid a) water molecules approaching the slab from the “wrong” side and b) decouple the long-range Coulomb interactions across the periodic boundary.