# [CP2K:4622] Re: energy drift in NPT MD simulation

Ari Paavo Seitsonen ari.p.s... at gmail.com
Fri Sep 20 05:17:05 UTC 2013

```Dear B.Satyanarayana,

My personal opinions:

1) What is determined as "converged" depends on the quantity observed and
accuracy wanted (thus on the person making the choice). For example the
pressure tensor, being a second derivative of the energy, converges much
slower than the forces (first-order derivative) that in turn converges much
slower than the energy (zero-order quantity)

2) In my opinion there are several simulations in the literature that are
not converged (and this does not apply to CP2K only). I don't make any
statement about the convergence in the article that you cited, just that it
is better to test oneself: After that one also understands the convergence,
the system, the method etc better

3) The energy on the line 'FORCE_EVAL' is the potential energy, not the
quantity that should be conserved. This shows that the energy does not
explode, not more (and the sample is very short, plus probably taken at the
beginning of the simulation when the things still vary more, due to the
rough equilibration needed

4) You see that the value of the determinant of the stress tensor diverges
_very_ fast, so there is certainly at least one problem. I was hinted that
the time constant that you have taken for the barostat is excessively
small: This would basically try to make the cell respond to the fastest
oscillations in your system (I hope that there _is_ something oscillating
at 4000 cm^-1, it would be necessary for the thermostat chain you have). I
don't know if this is the only reason, but I would switch the time constant
to something larger, for example of the order picoseconds

5) The convergence in the electronic structure indeed looks smooth on the
linear scale, but here one usually uses logarithmic scale in order to see
how the energy converges closer to the self-consistency. Another useful
check is to calculate the eigenvalues and see that the HOMO-LUMO gap (in
molecules) or between the valence and conducation band is more or less that
you expected

This after a quick glimpse on your report, maybe I missed again something
relevant though.

Good Luck, :)

apsi

2013/9/20 Bonakala Satyanarayana <satyanaray... at gmail.com>

> Dear Seitsonen,
> 1) I performed NVT simulation with same input, what I have enclosed. It
> ran well. These are the converged energy values:
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.764025944994501
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.749474147461115
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.722559495744463
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.689485930668525
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.657322690376532
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.631945457889742
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.616274832782437
>  ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
>  -7812.609299468986137
> 2) I have taken cut-off as per literature (Elucidating the Breathing of
> the Metal–Organic Framework MIL-53(Sc) with ab Initio Molecular Dynamics
> Simulations and in Situ X-ray Powder Diffraction Experiments, Duren et al.,
> JACS). Here, they have tested for  MIL MOFs and reported optimum cut-off is
> 350 Ry.
>
> Dear Hutter,
>
> 1) Yes, my system is stable in NVT simulation
>
> 2) Det( stress tensor) values are
>   Det(stress tensor)      :  -7.84528984E-04
>   Det(stress tensor)      :  -1.66559235E-02
>   Det(stress tensor)      :  -1.44160312E+00
>   Det(stress tensor)      :   2.03287049E+02
>   Det(stress tensor)      :  -4.82023324E+03
>   Det(stress tensor)      :   4.36251297E+07
>   BAROSTAT temperature is increasing:
>   BAROSTAT TEMP[K]             =          0.250067070500E+03
>   BAROSTAT TEMP[K]             =          0.407906043312E+05
>   BAROSTAT TEMP[K]             =          0.373591928441E+06
>   BAROSTAT TEMP[K]             =          0.458891514234E+08
>
> 3) SCF convergence is smooth. I enclosed the convergence vs steps
>
>
> Thank you
>
> B.Satyanarayana, Ph. D
> Molecular modeling lab
> JNCASR
> Bangalore
>
>
> On Thu, Sep 19, 2013 at 2:34 PM, Ari Paavo Seitsonen <
> ari.p.s... at gmail.com> wrote:
>
>> Dear Satyanarayana B,
>>
>>   Indeed the energy does not "drift", it explodes. What did happen in an
>> NVT simulation with otherwise the same input, did the simulation work fine?
>>
>>   To me the cut-off energy of 350 Ry looks low, even with the PBE
>> functional (smoother than BLYP); did you test that it is sufficient?
>>
>>     Greetings from Zurich,
>>
>>        apsi
>>
>> PS Does some one know, are the methods in &XC_GRID ... &END XC_GRID
>> consistent with an NPT simulation?
>>
>>
>> 2013/9/19 Matt W <MattWa... at gmail.com>
>>
>>> Hi,
>>>
>>> well your system is clearly exploding. Your cell_ref parameters look
>>> inconsistent with your actual cell - maybe the stress tensor is very bad!
>>>
>>>   &SUBSYS
>>>
>>>>     &CELL
>>>>       ABC 21.8620 16.2099 43.6020
>>>>       ALPHA_BETA_GAMMA 90.0 92.7540 90.0
>>>>      &CELL_REF
>>>>       ABC 22.8620 17.2099 44.6020
>>>>       ALPHA_BETA_GAMMA 90.0 126.0 90.0
>>>>      &END CELL_REF
>>>>
>>>>     &END CELL
>>>>
>>>>  --
>>> You received this message because you are subscribed to the Google
>>> Groups "cp2k" group.
>>> To unsubscribe from this group and stop receiving emails from it, send
>>> an email to cp2k+uns... at googlegroups.com.
>>> To post to this group, send email to cp... at googlegroups.com.
>>> Visit this group at http://groups.google.com/group/cp2k.
>>> For more options, visit https://groups.google.com/groups/opt_out.
>>>
>>
>>
>>
>> --
>>
>> -=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-
>>   Ari P Seitsonen / Ari.P.S... at iki.fi / http://www.iki.fi/~apsi/
>>   Physikalisch-Chemisches Institut der Universität Zürich
>>   Tel: +41 44 63 55 44 97  /  Mobile: +41 79 71 90 935
>>
>> --
>> You received this message because you are subscribed to the Google Groups
>> "cp2k" group.
>> To unsubscribe from this group and stop receiving emails from it, send an
>> email to cp2k+uns... at googlegroups.com.
>> To post to this group, send email to cp... at googlegroups.com.
>> Visit this group at http://groups.google.com/group/cp2k.
>> For more options, visit https://groups.google.com/groups/opt_out.
>>
>
>  --
> You received this message because you are subscribed to the Google Groups
> "cp2k" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to cp2k+uns... at googlegroups.com.
> To post to this group, send email to cp... at googlegroups.com.
> Visit this group at http://groups.google.com/group/cp2k.
> For more options, visit https://groups.google.com/groups/opt_out.
>

--
-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-
Ari P Seitsonen / Ari.P.S... at iki.fi / http://www.iki.fi/~apsi/
Physikalisch-Chemisches Institut der Universität Zürich
Tel: +41 44 63 55 44 97  /  Mobile: +41 79 71 90 935
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20130920/988b01d9/attachment.htm>
```