[CPMD-list] EHAM always increases

[Axel Kohlmeyer]

On Wed, 4 Jul 2007, liu liuxiandong wrote:

dear xiandong liu,

i'm cc'ing this to the cpmd mailing list,
so that the (potential) resolution to your problem is
not lost and others that find themselves in a similar 
situation can benefit from our discussion.

XL> Thank Axel.
XL> During the 300 steps of equilibration MD, EHAM always increase. Then I
XL> performed NVT CPMD with the following settings,
XL> ...

well, 300 steps of MD does not warrant to be named equilibration.
this is 43.5 fs, i.e. nothing. equilibration means that you'll
have to wait until you _reach_ equilibrium, i.e. until 
TEMPCONTROL IONS does not rescale anymore or the number of
up and down scalings will be about the same and the equivalent
of the potential energy in a classical MD (= EKS) does not 
drift but oscillates about a constant value for a while. 

even if your initial configuration is (classically) pre-equilibrated, 
i suggest to run for at least 0.5ps to make certain that some sort 
of equilibrium is achieved (given the short timescales of FPMD
simulations, you cannot expect a full equilibration) and then switch 
to nose-hoover chains. since you wrote previously that you started 
with a geometry optimization, you have to consider that it will 
need a bit longer to (re-)heat your system until it is back in 
equilibrium with the desired 300K temperature.

the easiest way for me to rationalize this situation is to
see an MD ensemble as an ensemble of coupled (an-harmonic) 
oscillators as outlined, e.g., in 
http://www.cpmd.org/pipermail/cpmd-list/2005-November/002935.html

so for a system at a given temperature the oscillators are
on average at a given r between r_0 and r_max. for 0K (i.e.
after an exhaustive geometry optimization) all of them will
be at r=r_0, for any higher temperature (or after a less 
thorough geometry optimization) r will be larger with increasing 
(previous) temperature. since the velocity rescaling, as implemented 
by TEMPCONTROL IONS primarily modifies currently excited modes, 
it will take some time until that energy is propagates through 
the system.

to speed up equilibration, it might be helpful to run some
part of the equilibration, i.e. after EKS is not drifting
anymore, with NOSE IONS MASSIVE, i.e. to have one thermostat
chain for each degree of freedom, so that all DOFs get equally
excited. i've seen in many cases, that after switching to 
massive thermalization (note, that this does _not_ imply
strong thermalization, but only that all DOFs are coupled 
to individual heat baths) systems dropped from meta-stable
quasi-equilibrium configurations into much deeper potential minima.

of course, the massive thermalization will mess up any dynamical 
information. if you are not interested in that, it may actually 
be beneficial to do even the production run with MASSIVE, so that 
you have better ergodicity, i.e. better sampling of phase space.

as i wrote in the previous e-mail, please have a look into
stat mech textbooks to learn about all the gory details.
a good understanding of the limitations imposed by small
system size and short trajectories will help you to make
the most out of your efforts. the fact that you got worried
about the drifting of the conserved quantity is a good start.

XL>   TIMESTEP
XL> 
XL>     6.0
XL> 
XL> 
XL> EMASS
XL> 
XL>     800.0
XL> 
XL> 
XL> 
XL>    TEMPERATURE
XL>    300.0
XL>    NOSE IONS
XL>    300.0 2760.0
XL>    NOSE ELECTRONS
XL>    0.008 10000.0

how did you obtain the characteristic frequencies 
for both nose thermostats?

also, i would first switch to N-H chains for the atoms
only and run for a while to make sure that adiabacity
is preserved and only then turn on the N-H chains for
electrons. for a 20000 step trajectory, a thermostat on 
the fictitious degrees of freedom should actually
not really be needed, given that you have reasonable
parameters for time step and fictitious mass (your system
does not have any hydrogens, right? or at least you made
them deuteriums...). 

turning on the thermostat on the electrons has the big
downside, that you are essentially 'blind' against
any changes of the system, where you'd lose adiabacity,
or cases where the initial wave-function optimization
optimized to the wrong state. in other words, with the
nose thermostat on the electrons, your MD will look 
as if the run was fine, but the trajectory may be bogus.

XL> ...
XL> However, during the 20,000 steps of such simulation, EHAM continues to
XL> increase. At the end of the run, the increasement of EHAM is near 1 a.u. I
XL> think this is unphysical.
XL> Then how to improve the simulation ? Thanks for any suggestion.

well, after the remarks from above it should be obvious, that
your system was not (yet) equilibrated. i would continue from 
the final configuration, perhaps adding a QUENCH BO and then
run for a bit in NVE ensemble to see how close you are to
the desired state and then continue with either more equilibration
or (real) production as outline before.

hope this helps,
    axel.


XL> 
XL> 
XL> 
XL> 
XL> Xiandong Liu

-- 
=======================================================================
Axel Kohlmeyer   akohlmey at cmm.chem.upenn.edu   http://www.cmm.upenn.edu
   Center for Molecular Modeling   --   University of Pennsylvania
Department of Chemistry, 231 S.34th Street, Philadelphia, PA 19104-6323
tel: 1-215-898-1582,  fax: 1-215-573-6233,  office-tel: 1-215-898-5425
=======================================================================
If you make something idiot-proof, the universe creates a better idiot.