TURBOMOLE Users Forum
TURBOMOLE Modules => Ridft, Rdgrad, Dscf, Grad => Topic started by: Andi K on September 01, 2009, 11:19:17 AM
-
Hello All,
I've encountered some problems with the use of $fermi in TURBOMOLE V6.0
I want to optimize the structure of some hydrocarbons with DFT (using the pbe functional for instance). For this purpose I use the jobex script with relax or statpt.
Seemingly during the optimization everything works quite alright. When I check the occupation with eiger after the run has finished, I see a negative HOMO-LUMO gap. To cope with this problem I tried to use $fermi with the following parameters
$fermi tmstrt= 25.00 tmend= 25.00 tmfac=1.000 hlcrt=1.0E+00 stop=1.0E-09
.
However, as the result of an ridft run I still don't get a positive HOMO-LUMO gap, even though the occupation is calculated with $fermi until the energy converges. The problem seems to be that after convergence an additional energy calculation step takes place, where $fermi for unknown reasons does not apply any more. Here are the details of the last few steps of my ridft output file:
ITERATION ENERGY 1e-ENERGY 2e-ENERGY NORM[dD(SAO)] TOL
62 -1179.8455826780 -6122.2798253 2670.1641575 0.523D-04 0.225D-12
Exc = -161.191934466 Coul = 2831.35609200
N = 181.99989069
current damping = 1.050
Norm of current diis error: 0.39805E-04
max. resid. norm for Fia-block= 3.340D-06 for orbital 91a
max. resid. fock norm = 3.475D-06 for orbital 91a
irrep a : virtual orbitals shifted by 0.10846
-------------------- FON Calculation ---------------------
Fermi level F = -0.12539
Total number of electrons N = 182.00000
Occupation numbers calculated for T = 25.00000
Current HOMO-LUMO gap = 0.48780E-01
----------------------------------------------------------
ITERATION ENERGY 1e-ENERGY 2e-ENERGY NORM[dD(SAO)] TOL
63 -1179.8455826783 -6122.2798664 2670.1641986 0.326D-04 0.225D-12
Exc = -161.191933054 Coul = 2831.35613170
N = 181.99989069
current damping = 1.100
Norm of current diis error: 0.34798E-04
max. resid. norm for Fia-block= 3.484D-06 for orbital 91a
max. resid. fock norm = 6.557D-06 for orbital 366a
irrep a : virtual orbitals shifted by 0.06828
-------------------- FON Calculation ---------------------
Fermi level F = -0.12539
Total number of electrons N = 182.00000
Occupation numbers calculated for T = 25.00000
Current HOMO-LUMO gap = 0.15104E-01
----------------------------------------------------------
ENERGY CONVERGED !
Overall gridpoints after grid construction = 665045
ITERATION ENERGY 1e-ENERGY 2e-ENERGY NORM[dD(SAO)] TOL
64 -1179.8455516457 -6122.2798324 2670.1641956 0.326D-04 0.790D-12
Exc = -161.191902780 Coul = 2831.35609842
N = 182.00000695
current damping = 1.150
Norm of current diis error: 0.25621E-04
max. resid. norm for Fia-block= 2.367D-06 for orbital 91a
max. resid. fock norm = 2.453D-06 for orbital 91a
End of SCF iterations
convergence criteria satisfied after 64 iterations
Using the eiger script I get the following occupation:
Total energy = -1179.8455516460 H = -32105.2492441 eV
HOMO-LUMO Separation
HOMO: 92. 91 a -0.12635850 H = -3.43839 eV
LUMO: 91. 92 a -0.14829956 H = -4.03544 eV
Gap : -0.02194106 H = -0.59705 eV
!! WARNING: HOMO-LUMO Gap is negativ !!
Number of MOs= 396, Electrons= 182.00, Symmetry: c1
Nr. Orbital Occupation Energy
108. 108 a +0.085909 H = +2.338 eV
107. 107 a +0.084492 H = +2.299 eV
106. 106 a +0.074339 H = +2.023 eV
105. 105 a +0.072635 H = +1.976 eV
104. 104 a +0.065869 H = +1.792 eV
103. 103 a +0.063547 H = +1.729 eV
102. 102 a +0.055610 H = +1.513 eV
101. 101 a +0.053687 H = +1.461 eV
100. 100 a +0.039238 H = +1.068 eV
99. 99 a +0.037697 H = +1.026 eV
98. 98 a +0.001710 H = +0.047 eV
97. 97 a -0.029023 H = -0.790 eV
96. 96 a -0.050084 H = -1.363 eV
95. 95 a -0.057593 H = -1.567 eV
94. 94 a -0.069351 H = -1.887 eV
93. 93 a -0.103417 H = -2.814 eV
92. 91 a 2.000 -0.126358 H = -3.438 eV
91. 92 a -0.148300 H = -4.035 eV
90. 90 a 2.000 -0.216450 H = -5.890 eV
89. 89 a 2.000 -0.221611 H = -6.030 eV
88. 88 a 2.000 -0.223412 H = -6.079 eV
87. 87 a 2.000 -0.241772 H = -6.579 eV
86. 86 a 2.000 -0.245310 H = -6.675 eV
85. 85 a 2.000 -0.253361 H = -6.894 eV
84. 84 a 2.000 -0.268394 H = -7.303 eV
83. 83 a 2.000 -0.293210 H = -7.979 eV
82. 82 a 2.000 -0.298763 H = -8.130 eV
81. 81 a 2.000 -0.303391 H = -8.256 eV
80. 80 a 2.000 -0.310105 H = -8.438 eV
79. 79 a 2.000 -0.332052 H = -9.036 eV
78. 78 a 2.000 -0.334060 H = -9.090 eV
77. 77 a 2.000 -0.334539 H = -9.103 eV
76. 76 a 2.000 -0.349200 H = -9.502 eV
What exactly is happening here? And as a general question: Why is this last energy step calculated, after the energy is converged?
Thanks very much for your help,
Andi
-
Hi,
your start and end temperatures are identical (although quite low), the annealing factor is 1.0, and your HOMO-LUMO gap is very small. It also seems that you have switched off the automatic stop of Fermi during the SCF iterations by setting the criterias (hlcrt and stop) to numbers that will never be reached.
The last iteration is done after convergence and it is due to the gridsize settings in $dft - most likely you have been using the default 'm3' which uses the larger grid only after convergence. That should not be the reason for the wrong occupation, but perhaps setting it to 3 or 4 instead of m3 can help if you have small gaps.
Some comments:
- $fermi in geometry optimizations can be dangerous, because it does not keep the multiplicity fixed (which can result in switching the potential energy surface during optimization) - however, you have a closed shell case, so it could be acceptable in this case.
- the $fermi option was thought as a black-box method to find the lowest spin state. If you know that you have a singlet, $fermi will most likely not help a lot.
- the homo-lumo (hlcrit) and the stop criteria are reasonable because if Fermi found a state with integer occupation, you will most likely want keep it
Some ideas on how to solve your problem:
- set the gridsize in $dft to 4
- remove $fermi
- add $denconv 1d-7 to the control file for a more tight density convergence
- use Turbomole 6.1 instead of 6.0 and add the new keyword $lastdiag - it might give more accurate orbital energies in cases the HOMO-LUMO gap is very small
Hope it helps,
Uwe