MD (Nose Thermostat) decaying constant of motion and another error regarding timestep

[dodo]

I am working on molecular dynamics of InGaZnO4, and I am facing two critical issues that I hope you can help me resolve. This is my INCAR file below:

SYSTEM = a-IGZO
NWRITE = 2
ISTART = 0
ENCUT = 205
NELM = 100
NELMIN = 4
NELMDL = -9
LREAL = Auto
ROPT = 0.01 0.01 0.01 0.01
PREC = Low
EDIFF = 1E-4
MAXMIX = 50
NSW = 2500
ALGO = NORMAL
NBLOCK = 1
KBLOCK = 10
ISYM = 0
IBRION = 0
ISMEAR = 0
SIGMA = 0.1
POTIM = 1.0
ISIF = 0
TEBEG = 3000
TEEND = 3000
SMASS = 1.02
NPACO = 104
APACO = 10.4
NGX = 46
NGY = 54
NGZ = 42
IWAVPR = 12

The lines related to SMASS/Nose frequency in OUTCAR are:

SMASS = 1.02 Nose mass-parameter (am)
estimated Nose-frequenzy (Omega) = 0.11E+15 period in steps = 57.41 mass= 0.781E-27a.u.

Where, the phonon mode frequency is chosen for ZnO thin film(J. Appl. Phys. 93, 126 (2003); http://jap.aip.org/resource/1/japiau/v93/i1/p126_s1)

Phonon mode frequency = 589 /cm = 0.111E+15 radian-Hz

Problem 1:

I cannot increase the timestep above 1 fs. Doing so gives error (POTIM = 2.0)-

POSCAR found : 4 types and 84 ions
LDA part: xc-table for Ceperly-Alder, standard interpolation

-----------------------------------------------------------------------------
| |
| WARNING
| |
| VASP found 249 degrees of freedom |
| the temperature will equal 2*E(kin)/ (degrees of freedom) |
| this differs from previous releases, where T was 2*E(kin)/(3 NIONS). |
| The new definition is more consistent |
| |
-----------------------------------------------------------------------------

POSCAR, INCAR and KPOINTS ok, starting setup
FFT: planning ...( 1 )
WAVECAR not read
prediction of wavefunctions initialized - no I/O
entering main loop
N E dE d eps ncg rms rms(c)
DAV: 1 0.546521894990E+04 0.54652E+04 -0.21486E+05 1024 0.116E+03
...
DAV: 18 0.111207663315E+03 0.50018E-05 -0.19260E-04 1152 0.352E-02
1 T= 24584. E= 0.37851039E+03 F= 0.11120766E+03 E0= 0.11120766E+03 EK= 0.26375E+03 SP= 0.00E+00 SK= 0.36E+01
bond charge predicted
N E dE d eps ncg rms rms(c)
DAV: 1 -0.313031662721E+03 -0.42424E+03 -0.14536E+04 1152 0.188E+02 0.442E+01
...
DAV: 15 -0.390078595642E+03 -0.33335E-04 -0.56904E-04 1280 0.819E-02
2 T= 72277. E= 0.49649121E+03 F= -.39007860E+03 E0= -.38999131E+03 EK= 0.77544E+03 SP= 0.37E+01 SK= 0.11E+03
Information: wavefunction orthogonal band 331 0.8934
Information: wavefunction orthogonal band 362 0.8987
Information: wavefunction orthogonal band 489 0.8493
...
Information: wavefunction orthogonal band 505 0.6247
Information: wavefunction orthogonal band 512 0.4685
bond charge predicted
prediction of wavefunctions
...

DAV: 25 0.241990413362E+02 0.39615E-04 -0.30969E-04 1280 0.308E-02
22 T= ****** E= Infinity F= 0.24199041E+02 E0= 0.24308749E+02 EK= 0.62872+163 SP= 0.14E+04 SK= Infinity
Information: wavefunction orthogonal band 7 0.8504
Information: wavefunction orthogonal band 8 0.8451
Information: wavefunction orthogonal band 10 0.8792
Information: wavefunction orthogonal band 22 0.8431
...
Information: wavefunction orthogonal band 508 0.6961
Information: wavefunction orthogonal band 509 0.5490
Information: wavefunction orthogonal band 510 0.5804


Is there some way I can increase the timestep and hence reduce the runtime ?

Problem 2:

The system (tell me if I am wrong) is an NVT ensemble. Hence, there must exist a constant of motion that is:

E_Ham = E_electrons + KE_ions + PE_NHC + KE_NHC

Here E_Ham should be constant.

On the basis of following lines in the OUTCAR file:

ENERGIE OF THE ELECTRON-ION-THERMOSTAT SYSTEM (eV)
---------------------------------------------------
% ion-electron TOTEN = 111.207673 see above
kinetic Energy EKIN = 92.874206 (temperature 8656.66 K)
nose potential ES = 0.000000
nose kinetic EPS = 0.176551
---------------------------------------------------
total energy ETOTAL = 204.258430 eV

I concluded that the ETOTAL should be this constant of motion.

However, in my OUTCAR it fluctuates in the beginnning and then consistently decays with time:

total energy ETOTAL = 204.258430 eV
total energy ETOTAL = 164.989970 eV
total energy ETOTAL = 172.552853 eV
total energy ETOTAL = 178.791949 eV
total energy ETOTAL = 182.925255 eV
total energy ETOTAL = 187.386634 eV
total energy ETOTAL = 191.054123 eV
?
total energy ETOTAL = 242.684919 eV
total energy ETOTAL = 244.698359 eV
total energy ETOTAL = 246.483614 eV
total energy ETOTAL = 248.056996 eV
...
total energy ETOTAL = 255.737343 eV
total energy ETOTAL = 255.904217 eV
total energy ETOTAL = 256.027848 eV
total energy ETOTAL = 256.115376 eV
total energy ETOTAL = 256.173518 eV
total energy ETOTAL = 256.208188 eV
...
total energy ETOTAL = 256.164627 eV
total energy ETOTAL = 256.167937 eV
total energy ETOTAL = 256.169711 eV
?
total energy ETOTAL = 251.107695 eV
total energy ETOTAL = 251.105754 eV
?
total energy ETOTAL = 231.711548 eV
total energy ETOTAL = 231.697482 eV
?
total energy ETOTAL = 244.481393 eV
total energy ETOTAL = 244.475985 eV
?
total energy ETOTAL = 226.332257 eV
total energy ETOTAL = 226.334046 eV
?
total energy ETOTAL = 220.245661 eV
total energy ETOTAL = 220.242658 eV
total energy ETOTAL = 220.241179 eV517 eV

Am I correct when I consider ETOTAL to be constant of motion for a Nose thermostat? If not, what is the actual constant of motion? If yes, why is my constant of motion decaying and how can I correct the error?

This is a rather long post though I have tried to keep it short. I hope however that will help to resolve these issues.

Thanks for reading.

[admin]

In the first step you are getting positive energy (and probably huge forces). This could be caused by:
1. Not correct order of atoms in POTCAR
2. Starting geometry is far away from realistic structure.
Is your configuration of atoms pre-relaxed? If not try it. With realistic structure MD sohuld work well also with increased timestep.

[dodo]

Increasing potim to 2.0 worked after I used a relaxed structure. Still some there were some occasional occurrences of:

Information: wavefunction orthogonal band 362 0.8987

in my OUTCAR. Do I need to be worried about these sporadic messages ?

[dodo]

Also, can somebody help me out with the constant of motion problem: Problem 2 mentioned above ?

[rabi]

I'm doing MD in VASP using vasp 5.3.2 and I'm getting same
Information: wavefunction orthogonal band 431 0.8849
Information: wavefunction orthogonal band 432 0.8855
Information: wavefunction orthogonal band 432 0.8821

I change IWAVPR = 11 as it was suggest on one post but this message still coming,

Is there any way to resolved this problem or just we can ignore it?

Thank you.