This tutorial explains how to optimize the geometry by using different algorithms. In this tutorial, a benzene molecule is placed in a cubic box with the cell edges of 10 Angstrom. pot_C.pbe1 and pot_H.lda1 are used with the cutoff energies of 36 Ry (GMAX=6) and 400 Ry (GMAXP=20) for the wave functions and the argmentation charge, respectively.
0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY) 6.00 20.00 2 12 12 : GMAX GMAXP NTYP NATM NATM2 25 0 : NUM_SPACE_GROUP TYPE_BRAVIS_LATTICE 18.89726878 18.89726878 18.89726878 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA 1 1 1 1 1 1 : K_MESH 1 0 : NCORD NINV 0.0000000000 2.6366036317 0.0000000000 1 1 1 -2.2833569864 1.3182923672 0.0000000000 1 1 1 -2.2833569864 -1.3182923672 0.0000000000 1 1 1 0.0000000000 -2.6366036317 0.0000000000 1 1 1 2.2833569864 -1.3182923672 0.0000000000 1 1 1 2.2833569864 1.3182923672 0.0000000000 1 1 1 0.0000000000 4.6909824123 0.0000000000 1 1 2 -4.0624970473 2.3454912062 0.0000000000 1 1 2 -4.0624970473 -2.3454912062 0.0000000000 1 1 2 0.0000000000 -4.6909824123 0.0000000000 1 1 2 4.0624970473 -2.3454912062 0.0000000000 1 1 2 4.0624970473 2.3454912062 0.0000000000 1 1 2 6.00 0.15 28.29 3 1 0.00 : TYPE IATOMN ALFA AMION ILOC IVAN ZETA1 1.00 0.15 28.29 3 1 0.00 : TYPE IATOMN ALFA AMION ILOC IVAN ZETA1 0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC 0 1 : IPRE IPRI 300 300 0 28000.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP 6 1 : WAY_MIX MIMX_WHAT 0 8 0.80 : ITER_START KBXMIX MIX_ALPHA 0.60 0.50 0.60 0.70 1.00 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM 500.00 4 1 0.10D-08 : DTIO IMDALG IEXPL EDELTA 0.0010 0.10D-03 0 : WIDTH FORCCR ISTRESS ggapbe 1 : XCTYPE KSPIN 1.00 : DESTM 101 : NBZTYP 4 4 4 : (DUMMY) 4 4 4 : (DUMMY) 18 : KEG 1 : NEXTST 0 : (DUMMY) 2 : IMSD 0 : EVALUATE_EKO_DIFF 0 : NPDOSAO 0 0.00 : SM_N DOPPINGIn this example, the GDIIS method is used (IMDALG=4) with the time step of 500 atomic unit. The GDIIS method works fine, when the atomic positions are close to the equilibirium. A rule of thumb is to switch from quenched molecular dynamics (IMDALG=2) to GDIIS when the maximum force becomes smaller than 1.e-2. Note that in the current implementation, GDIIS method does not work when the optimization step exceeds the 3 * N, where N is the number of atoms.
Monitoring the maximum force by using grep as (supposing the output file is "nfout_1")
grep -A1 f_max nfout_1
we obtain
NIT TotalEnergy f_max f_rms edel vdel fdel 1 -40.14215406 0.002809 0.001947 0.10D-08 0.12D-06 0.10D-08 -- NIT TotalEnergy f_max f_rms edel vdel fdel 2 -40.14215395 0.005734 0.004041 0.22D-08 0.17D-06 0.22D-08 -- NIT TotalEnergy f_max f_rms edel vdel fdel 3 -40.14224266 0.001880 0.001393 0.47D-08 0.13D-06 0.47D-08 -- NIT TotalEnergy f_max f_rms edel vdel fdel 4 -40.14225368 0.000871 0.000552 0.11D-08 0.15D-06 0.11D-08 -- NIT TotalEnergy f_max f_rms edel vdel fdel 5 -40.14225523 0.000046 0.000030 0.51D-09 0.85D-07 0.51D-09
and we can see that the geometry optimization converges within 5 steps.
Then, let us use the quenched molecular dynamics (QMD) for comparison. In this case we use IMDALG=2 and the time step (DTIO) of 200 atomic unit as
200.00 2 1 0.10D-08 : DTIO IMDALG IEXPL EDELTA
and the resulting maximum forces are
NIT TotalEnergy f_max f_rms edel vdel fdel 1 -40.14215406 0.002809 0.001947 0.10D-08 0.12D-06 0.10D-08 -- NIT TotalEnergy f_max f_rms edel vdel fdel 2 -40.14218382 0.001971 0.001448 0.22D-08 0.10D-06 0.22D-08 -- ... -- NIT TotalEnergy f_max f_rms edel vdel fdel 14 -40.14225511 0.000195 0.000135 0.15D-08 0.66D-07 0.15D-08 -- NIT TotalEnergy f_max f_rms edel vdel fdel 15 -40.14225523 0.000071 0.000033 0.21D-09 0.23D-07 0.21D-09
The convergence is achieved with 15 steps and we can see that the GDIIS method is much more efficient than QMD.