I have been running some MD simulations using fix efield, and noticed that
"something strange" was going on with the pressures that I was getting
out... so I did some tests.
Consider a single, rigid, diatomic molecule, aligned with the z-axis, in a
non-periodic simulation cell, at zero temperature. (I.e. all of the atoms
are at rest.)
The molecule carries no net charge, but the individual atoms have charges
of +1e and -1e respectively. There are no forces defined between the atoms,
and the geometry of the molecule is kept constant using fix shake.
What do I do?
I run the simulation, from this starting point, using the NVE time
integrator, and monitor the pressure.
What do I expect?
That the pressure should be zero, because the atoms are at rest, and there
are no forces between them. Even under an applied (constant, uniform)
field, a rigid dipole could experience a torque, but no net force.
What do I observe?
That the pressure returned by the "thermo" keyword scales linearly with
the applied field.
Why might this be?
From reading the LAMMPS manual and previous emails on this list, I have
gathered that (i) External forces are not included in the virial of the
pressure. (ii) Forces due to fixes that apply constraints are. Whilst I am
not questioning the reasons for which these choices where originally made,
I believe that in the current scenario this results in unphysical
behaviour. (Pressure in a box with one rigid molecule at rest.)
Is pressure a well-defined quantity in simulations under an applied
Maybe not in general, however... in a system that satisfies the following
criteria: (i) no "free" charge carriers (ii) homogeneous (iii) electric
field is constant (iv) electric field is uniform; I cannot see why it
What would be a solution?
It might be that including the forces due to the electric field in the
virial of the pressure would solve the problem. Alternatively, a more
advanced treatment of the shake forces in which the response to
particle-particle interactions can be distinguished from the response to
external forces could probably achieve the same result.
Am I thinking along the right lines here, and do you think that a solution
would be feasible? Any help would be greatly appreciated.
I can post the input and output files of my test job if they are of
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
lammps-users mailing list