Localized electric field

Dear all,
I need some help to figure out how to solve this puzzling problem.

I have to simulate the effect of an electric field on nanoconfined water. The sample is a box of SPC/E water with two parallel graphene-like sheets (perpendicular to z axis) in the center of the box, so that water fills the space between the plates and also creates a reservoir around them. The electric field must be applied only between the sheets.
I have modified the fix_efield class in order to apply the field only on a cylindrical volume, closed on the top and at the bottom by the two graphene layers. To avoid the discontinuity along the other two dimensions, I smoothed the field by using a smooth function which is continuous and also all its derivatives are continuous. Due to the presence of the impermeable graphene sheets, no discontinuity field is expected for water along z direction.
However, the total energy of the system is conserved when the field is off, while it shows a positive drift when it is switched on.

Am I missing something? Is there something on the force integration algorithm that does not allow me to do something like that?

Let me know if I missed any details…

Thanks for your help.


If you can't conserve energy with fix NVE then the force
integration algorithm is dead simple and can't be the
problem. Can you put a single water in your field
and conserve energy?


I’m sorry Steve, but I don’t understand your first statement.

I tried to put a single water molecule with zero initial velocity on the field. If it is in the central part of my field-applied volume, where the field is constant (remember that I apply the field smoothing only at the boundary between the perturbed and unperturbed volume in a shell of 4 angstrom thickness), it aligns with the field as expected. When I put on the non-uniform field shell, it starts to rotate due to the torque that comes from different forces acting on three atoms.


What I meant was that if you are not conserving energy
there are a few possibilities:

a) bad timestep
b) bad fix nve
d) your added electric field

I'm assuming it's not a. It's also not b
b/c fix nve is too simple. It could be that
(c) causes problems since the field will
tend to pull the +/- charges in each molecule
in different directions. SHAKE can correct
for it, but if the field is huge that might induce
small errors. You could test this by running with
flexible water.

That leaves (d), which is the likely culprit. How
do you know your field is not adding energy to
the system? Can you simulate a single water
moleclue in your field and conserve energy?