UT theoretical chemistry code forum

Hi all

I have been running a dimer caluclation for TS state (CO to CO2 conversion on ceria 100 surface), however I have been getting consistently high torque which does not drop as the dimer rotates, as below for DROTMAX = 5 (I have tried already with DROTMAX = high values like 50):

Step Force Torque Energy Curvature Angle
1 0.28664 7.81732 -787.21676 -4.19537 6.03221
1 0.28664 8.77560 -787.21676 -4.84612 11.60515
1 0.28664 7.66106 -787.21676 -4.98581 1.08454
1 0.28664 7.11642 -787.21676 -4.86630 3.54620
1 0.28664 7.56927 -787.21676 -4.81569 4.25652
2 0.26788 12.91125 -787.21721 -3.74021 17.82935
2 0.26788 15.15547 -787.21721 -4.72224 5.21785
2 0.26788 5.95358 -787.21721 -5.44741 9.17223
2 0.26788 9.38043 -787.21721 -5.05897 3.20145
2 0.26788 6.80521 -787.21721 -5.57626 4.40569
3 0.24693 12.08755 -787.21793 -6.66529 10.03584
3 0.24693 12.79037 -787.21793 -6.52186 4.17402
3 0.24693 4.63119 -787.21793 -7.75667 1.54969
3 0.24693 7.05629 -787.21793 -8.07437 5.55998
3 0.24693 6.16010 -787.21793 -7.01064 2.20227
4 0.23944 9.12140 -787.21873 -4.76547 4.02241
4 0.23944 9.04814 -787.21873 -4.99669 6.85383
4 0.23944 11.58510 -787.21873 -4.30967 4.37353
4 0.23944 5.03026 -787.21873 -5.68529 2.47758
4 0.23944 6.96767 -787.21873 -5.71619 5.63656
5 0.23743 8.61981 -787.21962 -5.32407 6.32159
5 0.23743 11.42362 -787.21962 -5.52887 7.58141
5 0.23743 10.13277 -787.21962 -4.75671 6.26814
5 0.23743 7.27833 -787.21962 -6.07716 0.19825
5 0.23743 7.30257 -787.21962 -6.05946 5.72035
6 0.25804 14.73318 -787.22080 -5.91507 -5.37293
6 0.25804 14.36086 -787.22080 -6.38876 0.78376
6 0.25804 14.13604 -787.22080 -6.40975 0.38567
6 0.25804 13.80968 -787.22080 -6.48027 -0.32053
6 0.25804 14.14131 -787.22080 -6.40685 0.50608
7 0.51558 14.27677 -787.22144 -3.66867 -4.88039
7 0.51558 13.35692 -787.22144 -4.02141 1.08477
7 0.51558 13.34678 -787.22144 -4.04167 0.60734
7 0.51558 13.48369 -787.22144 -3.99010 0.11737
7 0.51558 13.48525 -787.22144 -3.98900 0.07190
8 0.50316 12.35017 -787.22191 -3.90686 -1.65733
8 0.50316 12.56673 -787.22191 -3.81836 1.64947
8 0.50316 12.35459 -787.22191 -3.99985 -1.32129
8 0.50316 12.57227 -787.22191 -3.82421 1.51330
8 0.50316 12.36183 -787.22191 -3.99756 -1.40912
9 0.48040 15.03833 -787.22234 -2.66833 1.27501
9 0.48040 14.97974 -787.22234 -2.77271 -0.00915
9 0.48040 14.98146 -787.22234 -2.77094 0.00190
9 0.48040 14.98142 -787.22234 -2.77089 0.00176
9 0.48040 14.98143 -787.22234 -2.77094 0.00073
10 0.45933 14.27192 -787.22279 -3.47125 1.87191
10 0.45933 14.95336 -787.22279 -5.88216 -5.37672
10 0.45933 14.48264 -787.22279 -3.02841 7.46328
10 0.45933 15.21630 -787.22279 -5.73117 -7.67854
10 0.45933 14.49447 -787.22279 -3.01070 7.53638
11 0.44147 16.55120 -787.22323 -5.20516 20.72102

The INCAR is as follows:

#Must be standard in runs that will be compared
PREC = NORMAL
LREAL = AUTO
VOSKOWN = 1


# May be changed to tune speed/convergence and still compare calcs of different parameters
ALGO = Very_Fast ##Updated 8/07 Very_Fast = IALGO 48, Normal = IALGO 38, Fast = some of each
NELM = 400
ISYM = 0
ISMEAR = 2; SIGMA = 0.2 ## broadening in eV. Should be about 0.2 for metals, 0.003 for molecules and semiconductors (oxides)

# System/calculation specific and extremely important to confirm
ISPIN = 2 #spin polarized calculation? 1= no, restricted, 2=yes, unrestricted
ENCUT = 500


# Add system specific or run specific parameters that stand out here (U terms, NELECT, Dipole corrections)
LVTOT = .FALSE.


#SYSTEM = Typical Dimer Run
IBRION=3
POTIM=0.0
ISYM=0
NSW=800
EDIFF=1E-7
EDIFFG=-0.05
LWAVE=.FALSE.
LCHARG=.FALSE.


# DIMER PARAMETERS
ICHAIN=2
DdR=0.005
DRotMax=5
DFNMin=0.01
DFNMax=1.0

#OPTIMIZER PARAMETERS
IOPT=3

#U correction for Ceria
LDAU = .TRUE.
LDAUTYPE = 2
LDAUL = -1 3 -1 -1
LDAUU = 0 5 0 0
LDAUJ = 0 0 0 0
AMIX = 0.2
BMIX = 0.0001
AMIX_MAG = 0.8
BMIX_MAG = 0.0001
LMAXMIX = 6

I have already tried with IOPT =2 as well, but since the torques doesnt drop, the forces started increasing from ~0.15 eV to 0.35 eV when dimer started the next iteration.
Let me know whats is wrong here. Have already tried low enough EDIFF of 10^-7, not sure going even lower would help.

Thanks

Something is not right. I would start by checking to see if the electronic structure is converging to your 1e-7 EDIFF value within your specified 400 NELM steps. If that is not happening, there will be errrors in the force and the torque values will be meaningless.

If your forces are converged, then I might need to look at the calculation. Just a .tar.gz file of the input parameters would be sufficient.

Dr. Henkelman

Thanks for your reply. Yes the SCFs are converging well here. Please find attached the current run below with the input parameters. Just to tell you about the system that I am running it on :it is a PdO on ceria100 surface for CO to CO2 conversion where oxygen from one side are taken to the other side for dipole correction (suggested in the literature).

Thanks in advance.

I did a calculation using what I would consider pretty standard settings for a TS calculation. I cut your slab down a bit, froze the bottom layers, reduced the kpoints, went to soft O and C potentials and lowered the cutoff to 300 eV.

In the attached you can see reasonable convergence. I was not going for fast, I just wanted to make sure that the system would converge systematically. Once you have a good approximation to the saddle, you can scale up your system size or parameters to increase precision.

Thanks for taking time to converge the dimer. However, as I see - I need to redo a geometric optimisation with the original kpoints and original slab width to get the correct energy.

I still do not understand the origin of high torque as I described earlier, though lower kpoints and lower ENCUTs and softer potentials as done above partially removed that (torques as high as ~2 eV are still there in the attached DIMCAR). were these the origin of high torques ? I have been getting high torques in other systems as well with Fe, Ni or Co as monolayer over Pd with respective U corrections and I do not know if doing as above should be the right approach to deal with each of theses cases.

Also found out the following warning errors through out the log files and OUTCAR files.

WARNING in EDDRMM: call to ZHEGV failed, returncode = 6 3 4 with different numbers each time during each SCF cycles. The SCF cycles still converge though even after these warnings.

Such errors have been suggested to go away with change in ALGO tags in this forum. Should that be the main cause of high torque ?

I'm not sure what was going wrong with your original settings. I used the cheap settings, not to solve the problem, but to debug the problem. When I didn't see the high forces with my normal settings, I left it as a solved problem. For some reasons, your original settings appear to give noise in the forces so that the curvature calculations are unreliable. It could be related to the messages about the linear algebra routines, but I also doubt it.

If this becomes an ongoing problem, I can look into it further. But I strongly recommend finding saddles as quickly as possible with as cheap settings as you can get away with, and then reconverging them for higher precision. This will ultimately be a faster approach and you will get a measure of uncertainty having calculations at two accuracy settings.

Hey Graeme

I agree with the NEB/dimer at cheap settings, however sometimes the reconverging the TS state back to the original settings and higher kpoints take time (I generally use POTIM 0.05/0.1 for reconverging).
As for the high torque, I have been getting them on other systems as well. I will see and share soon so that you might share some insights on the peculiar behaviour.

Thanks