Dear Paul, and Kunwoo,
Thank you so much for the detailed reply.
I am simulating a relatively small system (~5000 atoms)
using LJ cutoff and the cutoff is 10 Angstroms.
Calculation of long range electrostatic forces are actually
very important in my simulation, because the polymer has
ionic side chains. I am using PPPM for this calculation.
I know that the energy would fluctuate, but in this case it fluctuates
around an average value, which gradually increases. My test run was
relatively
long (100 ps), And the system was equilibrated for 200 ps before the
test run.
Probably, as you have suggested, using SHAKE is the best way to deal
with
this problem.