From Diracwiki

Contents 1 Occupation 1.1 .CLOSED SHELL 1.2 .OPEN SHELL 1.3 .BOSONS 1.4 .MJSELE 1.5 .AOC 1.6 .FOCC 1.7 .AUTOCC 2 Trial function 2.1 .NOBNCR 2.2 .TRIVEC 2.3 .TRIFCK 3 Convergence criteria 3.1 .MAXITR 3.2 .EVCCNV 3.3 .ERGCNV 3.4 .FCKCNV 4 Convergence acceleration 4.1 .DIISTH 4.2 .MXDIIS 4.3 .DIISMO 4.4 .DIISAO 4.5 .NODIIS 4.6 .DAMPFC 4.7 .NODAMP 5 Level shift 5.1 .LSHIFT 5.2 .OLEVEL 6 State selection 6.1 .OVLSEL 6.2 .NODYNSEL 7 Iteration speedup 7.1 .NODSCF 7.2 .FIXDIF 7.3 .CNVINT 7.4 .ITRINT 7.5 .INTFLG 8 Output control 8.1 .PRINT 8.2 .EIGPRI 9 Eliminating/freezing orbitals 9.1 .PROJECT 9.2 .PRJTHR 9.3 .OWNBAS 9.4 .FROZEN 10 Sample input 11 References

This section gives directives for Hartree-Fock calculations. Kohn-Sham calculations are performed by invoking the keyword .DFT under **HAMILTONIAN. Open-shell calculations correspond to either average-of-configurations (.AOC) or fractional occupation (.FOCC). The former is the default for Hartree-Fock calculations, whereas the latter is default for Kohn-Sham calculations. Note that average-of-configurations Kohn-Sham calculations are not well defined.

Occupation

.CLOSED SHELL

For each fermion irrep give the number of closed shell electrons.

The specification of the closed shell electrons is simple. For symmetry groups without inversion symmetry, there is only one fermion irrep, and you need only to specify the number of electrons.

Example: Minimum wave function input for water (the C_2v point group):

**WAVE FUNCTION
.DHF
*DHFCAL
.CLOSED SHELL
 10

For symmetry groups with inversion symmetry, you need to specify the distribution of the electrons in the two fermion irreps ^[1].

Example: Minimum wave function input for the Neon atom:

**WAVE FUNCTION
.DHF
*DHFCAL
.CLOSED SHELL
 4 6

.OPEN SHELL

Specification of open shell(s).

For each open shell give the number of electrons and the number of active spinors.

Example:

.OPEN SHELL
 1
 5/0,6

1 open shell with 5 electrons in 6 spinors (= 3 Kramers pairs) in irrep 2. Thus, the fractional occupation is 5/6.

The open shell module in DIRAC is based on average-of-configurations ^[2]. The simplest case is one electron in two spinors (= one Kramers pair). For this special case the average-of-configuration calculation gives the same result as the usual restricted open-shell Hartree-Fock. For all other cases the calculation gives the average energy of many states.

An example of one electron in two spinors is the Boron atom:

**WAVE FUNCTION
.DHF
*DHFCAL
.CLOSED SHELL
 4 0
.OPEN SHELL
 1
 1/0,2

The input is interpreted in the following way: We have 4 inactive electrons (in the 1s and 2s orbitals) and one open shell, where 1 electron is in two spinors (the 2p_{1 / 2} and spinors) in irrep 2. This calculation would give you the P_{1 / 2} energy of Boron. We could also do a calculation with

**WAVE FUNCTION
.DHF
*DHFCAL
.CLOSED SHELL
 4 0
.OPEN SHELL
 1 
 1/0,6

where we have the active electrons distributed in all six p-orbitals. This would give us the average energy of six states:

It is also possible to specify more than one open shell. Let us take a look at two different inputs for the open shells of the Platinum atom (6s¹5d⁹):

In the first example:

.OPEN SHELL
 2
 9/10,0
 1/2,0

we get the average energy of 20 states:

.

All 20 states are s¹d⁹ states.

In the second input:

.OPEN SHELL
 1
 10/12,0

we distribute the same ten electrons in all 12 spinors. This gives a total of states, which are a mixture of one s⁰d¹⁰ state, 20 s¹d⁹ and 45 s²d⁸ states.

Note that the order of closed and open shells are assumed to be as in the following scheme:

Virtual orbitals Not occupied
...
Open shell 2 Fractionally occupied
Open shell 1 Fractionally occupied
Closed shell Doubly occupied

that is, the lowest molecular orbitals are doubly occupied, the next ones are occupied with the electrons of open shell no. 1, etc.

Other orderings can be achieved by using .REORDER MO'S and .OVLSEL.

To get the energies of the individual states present in the average-of-configurations, specify .RESOLVE (see also *RESOLVE).

To get the energies of (some) of the individual states present in the average of configurations, you can use the GOSCIP module, the DIRRCI module, or the LUCITA module.

.BOSONS

Occupation of boson irreps in spin-free calculation. For example, for the D2h symmetry eight numbers in subsequent line, for the C2v symmetry there four occupation numbers in line.

.MJSELE

Occupation of MJ-splitted fermion irreps (only for the linear symmetry). Next line contains sequence of occupation numbers, see the SCF output for how many.

.AOC

Average-of-configuration calculation (default for open-shell Hartree-Fock).

.FOCC

Fractional occupation (default for open-shell Kohn-Sham).^{[please clarify]}

.FOCC calculations are less memory-intensive than .AOC calculations. In the latter case one additional AO-Fock matrix is generated for each open shell.

.FOCC calculations are therefore an interesting option for generating start orbitals for MCSCF as well as initial convergence in open-shell Hartree-Fock.

.AUTOCC

Program is allowed to change occupation during SCF cycles. This is deactivated by default. However, the program will still try to do an automatic initial occupation if neither .CLOSED SHELL nor .OPEN SHELL is given.

Trial function

An SCF-calculation (HF or DFT) may be initiated in three different ways:

• Using MO coefficients from a previous calculation.
• Using coefficients obtained by diagonalization of the one-electron Fock matrix: the bare nucleus approach.
• Using the two-electron Fock matrix from a previous calculation; this may be thought of as starting from a converged DHF potential.

The default is to start from MO coefficients if the file DFCOEF is present. Otherwise the corrected bare nucleus approach is followed. In all three cases linear dependencies are removed in the zeroth iteration.

.NOBNCR

Switch off the bare nucleus correction. This correction is on per default and improves the screening of the nuclei by estimating two-electron repulsion via nuclear-attraction type integrals:

The coefficients a_j and the exponents in this expression are chosen according to Slater's rules to obtain an approximate atomic electronic density for the initial guess. For example, with one heavy element and without this correction (that is, with the bare nucleus Hamiltonian) all electrons will end up on that heavy element in the initial guess!

.TRIVEC

Start SCF-iterations from the vector file.

.TRIFCK

Start SCF-iterations from the two-electron Fock matrix from previous calculation (stored on file DFFCK2).

Convergence criteria

Three different criteria for convergence may be chosen:

• The norm of the DIIS error vector (in MO basis). This corresponds to the norm of the electronic gradient and is the recommended convergence criterion. When you are only interested in the energy .EVCCNV = 1.0D-5 is usually sufficient. For properties and correlated methods you should converge to .EVCCNV = 1.0D-9. Large negative energy eigenvalues lead to a loss of precision that might lead to convergence problems. Remember also that a too loose screening threshold (too many integrals neglected) will hinder convergence. You should modify .SCREEN under *TWOINT if you modify .EVCCNV or one of the other two convergence criteria.
• The difference in total energy between two consecutive iterations.
• The largest absolute difference in the total Fock matrix between two consecutive iterations.

The change in total energy is approximately the square of the largest element in the error vector or the largest change in the Fock matrix. The default is convergence on electronic gradient with 1.0D-6 as threshold. Alternatively, the iterations will stop at the maximum number of iterations.

Sometimes it may happen that the specified convergence criterion is too tight for the given basis set and/or other input parameters. In this case one needs to decide whether one should proceed with post-HF steps (like correlation calculations) or not. The program decides this by looking at a secondary convergence criterion that gives the allowed convergence. This value is by default the same as first or desired convergence criterion but can be made lower to make sure that a calculation does not abort when the convergence is slightly above threshold.

For more detailed help see SCF convergence troubleshooting or here.

.MAXITR

Maximum number of SCF iterations.

Default:

.MAXITR
 50

When restarting SCF itrations from previous molecular orbitals file (DFCMO or formatted DFPCMO), we recommend to decrease maximum number of iterations together with readjusting desired and allowed convergence thresholds. By properly set desired and allowed thresholds one can have exact number of iterations specified by .MAXITR.

.EVCCNV

Converge on error vector (electronic gradient).

2 (real) Arguments:

.EVCCNV
 desired threshold allowed threshold

.ERGCNV

Threshold for convergence on total energy.

2 (real) Arguments:

.ERGCNV
 desired threshold allowed threshold

.FCKCNV

Converge on largest absolute change in Fock matrix.

2 (real) Arguments:

.FCKCNV
 desired threshold allowed threshold

Note that the allowed threshold may be omitted. It is then made equal to the desired threshold.

Convergence acceleration

It is imperative to keep the number of SCF iterations at a minimum. This may be achieved by convergence acceleration schemes:

• Damping: The simplest scheme is damping of the Fock matrix that may remove oscillations. In iteration n + 1 the Fock matrix to be diagonalized is: , where c is the damping factor.
• DIIS: Direct inversion of iterative subspaces ^{[3] [4] [5]} may be thought of as generalized damping involving Fock matrices from many iterations. Damping factors are obtained by solving a simple matrix equation involving the B-matrix constructed from error vectors (approximate gradients). Linear dependent columns in the B-matrix is removed.

In DIRAC DIIS takes precedence over damping.

.DIISTH

Change the default convergence threshold for initiation of DIIS, based on largest element of error vector.

Default:

.DIISTH
 a very large number

.MXDIIS

Maximum dimension of B-matrix in the DIIS module.

Default:

.MXDIIS
 10

.DIISMO

Activate DIIS in orthogonal basis (MO) with the error vector as described above ().

.DIISAO

Activate DALTON-like DIIS using AO-basis. The error vector is

where is given by

.NODIIS

Do not perform DIIS. The default is to activate .DIISMO for closed-shell calculations, and to activate .DIISAO for average-of-configurations calculations.

.DAMPFC

Change the default damping factor.

Default:

.DAMPFC
 0.25

.NODAMP

Do not perform damping of the Fock matrix. Damping is activated by default, but DIIS takes precedence. In case all columns in the B-matrix is removed by linear dependency, damping is activated.

Level shift

.LSHIFT

Activate level shift (for virtual orbitals). Followed by a real argument (level shift).

.OLEVEL

Activate level shift (for open-shell orbitals). Followed by a real argument (level shift), one line for each open shell {Please give example}.

State selection

Convergence can be improved by selection of vectors based on overlap with vectors from a previous iteration. This method may also be used for convergence to some excited state.

If dynamic overlap selection is used, the vector set from the previous iteration is used as the criterion. For the first iteration either restart vectors or vectors generated by the bare nucleus approach (not recommended) are used.

If .NODYNSEL is given, either the restart vectors or the bare nucleus vectors are used, i.e. the overlap selection vectors are not updated in each iteration. Please note, that overlap selection based on vectors from the bare nucleus approach is not recommended.

Overlap selection is very useful together with .REORDER MO'S. This will reshuffle the vectors within the restart coefficients. Example: First one might do a open shell calculation on Boron, this would give the P_{1 / 2} state. But if we restart on the P_{1 / 2} coefficients, interchange the p_{1 / 2} with the p_{3 / 2} orbitals, and request overlap selection, we can converge to the P_{3 / 2} state.

There also exists a keyword for reordering the converged DHF orbitals. This is useful for reordering the orbitals for the 4-index transformation and subsequent correlation calculations (CCSD, CI etc.) (see .POST DHF REORDER MO'S).

.OVLSEL

Activate dynamic overlap selection. The default is no overlap selection.

.NODYNSEL

No dynamic update of overlap selection vectors. The default is dynamic update.

Iteration speedup

The total run time may be reduced significantly by reducing the number of integrals to be processed in each iteration:

• Screening on integrals: Thresholds may be set to eliminate integrals below the threshold value ^[6]. The threshold for LL integrals is set in the basis file, but this threshold may be adjusted for SL and SS integrals by threshold factors set in the **INTEGRALS section.

• Screening on density: In direct mode further reductions are obtained by screening on the density matrix as well ^[6]. This becomes even more effective if one employs differential densities, that is . The default value for α is which corresponds to a Gram-Schmidt orthogonalization. As SCF converges, α goes towards 1, but α can also explicitly be set equal to 1 with .FIXDIF.

• Neglect of integrals: The number of integrals to be processed may be reduced even further by adding SL and SS integrals only at an advanced stage in the DHF-iterations, as determined either by the number of iterations or by energy convergence. The latter takes precedence over the former.

.NODSCF

Do not perform SCF-iterations with differential density matrix.

Default: Use differential density matrix in direct SCF.

.FIXDIF

Set α equal to 1.

.CNVINT

Set threshold for convergence before adding SL and SS integrals to SCF-iterations.

2 (real) Arguments:

.CNVINT
 CNVINT(1) CNVINT(2)

.ITRINT

Set the number of iterations before adding SL and SS integrals to SCF-iterations.

Default:

.ITRINT
 1 1

.INTFLG

Specify what two-electron integrals to include (see .INTFLG under **HAMILTONIAN).

Default: .INTFLG from **HAMILTONIAN.

Output control

.PRINT

General print level for the SCF method. For instance, value of 2 prints eignevalues during each iteration.

Default:

.PRINT
 0

.EIGPRI

Controls the print-out of positive energy and negative energy eigenvalues (1 = on; 0 = off).

Default: Only the positive energy eigenvalues are printed.

.EIGPRI
 1 0

Eliminating/freezing orbitals

In studies of electronic structure it may be of interest to eliminate or freeze certain orbitals. This option is furthermore useful for convergence, in particular to excited electronic states. A simple case is the thallium atom. The ground state ²P_{1 / 2} has the electronic configuration [Xe]. The first excited state ²P_{3 / 2} with the electronic configuration can easily be accessed by first calculating the ground state, then eliminating the from the ensuing calculation. In the final calculation the excited state is relaxed using overlap selection. The use of frozen orbitals is demonstrated in test/33.frozen: When the geometry of the water molecule is optimized with the oxygen 1s and 2s orbitals frozen, a bond angle of 96.242 degrees is found, contrary to the 90 degrees one might have expected when s-p hybridization is thus blocked ^[7].

The elimination of orbitals is achieved by projecting the selected orbitals out of the transformation matrix to orthonormal basis. The selected orbitals can be expressed either in the full molecular basis or in the basis set of the chosen fragment. In the latter case, the same set of fragment orbitals can in the case of atomic fragments be used at different molecular geometries. One may even perform a geometry optimization, but only using the numerical gradient. When freezing orbitals the selected orbitals are first eliminated from the transformation to orthonormal basis, but then reintroduced in the backtransformation step. They will appear in the output with zero orbital eigenvalues. Note that when freezing orbitals the orbitals to be eliminated must be specified as well. The frozen orbitals must be a subset of the eliminated orbitals.

Fragments are defined with respect to the list of symmetry-independent atoms appearing in the DIRAC basis file: Consider the water molecule in the full C_2v symmetry. Then there are two fragments: the oxygen atom and the H₂ moiety. However, with no symmetry there will be three fragments: the oxygen atom and the two hydrogen atoms. At the moment there are no orthonormalization of fragments on different fragments and so in practice one should only use orbitals from one fragment.

.PROJECT

Eliminate orbitals by projecting them out of the transformation matrix to orthonormal basis.

Arguments: Number of fragments (NPRJREF).

Then for each fragment read name PRJFIL of the coefficient file followed by the number of symmetry-independent nuclei in this fragment followed by an orbital string of selected orbitals for each fermion irrep.

DO J = 1,NPRJREF
  READ(LUCMD,'(A6)') PRJFIL(J)
  READ(LUCMD,'(I6)') NPRJNUC(J)
  READ(LUCMD,'(A72)') (VCPROJ(I,J),I=1,NFSYM)
ENDDO

.PRJTHR

Smallest norm accepted when eliminating orbitals.

Default:

.PRJTHR
 1.0D-10

.OWNBAS

Eliminated/frozen orbitals are given in the fragment basis. Note that the list of fragments is assumed to follow the list of symmetry independent nuclei in the DIRAC basis file.

.FROZEN

Freeze orbitals (must only be used in conjunction with .PROJECT).

For each fermion irrep, give an orbital string of orbitals to freeze.

Sample input

*DHFCAL
.PRINT
1
.NELECT
10 0
.MAXITR
30
.EVCCNV
1.0D-6
.MXDIIS
10
.CNVINT
1.0D-2 1.0D-4
.EIGPRI
1 0

References

1. ↑ T. Saue and H. J. Aa. Jensen, Quaternion symmetry of the Dirac equation, in Mathematical Methods for Ab Initio Quantum Chemistry, edited by M. Defrancheschi and C. Le Bris, Springer, Berlin, 2000.
2. ↑ J. Thyssen and H. J. Aa. Jensen, Average-of-configurations SCF manuscript, unpublished (1998).
3. ↑ P. Pulay, Convergence acceleration of iterative sequences. The case of SCF iteration, Chem. Phys. Lett. 73, 393 (1980) electronic version.
4. ↑ P. Pulay, Improved SCF acceleration, J. Comput. Chem. 3, 556 (1982) electronic version.
5. ↑ T. Hamilton and P. Pulay, Direct inversion in the iterative subspace (DIIS) optmization of open-shell, excited-state, and small multiconfiguration SCF wave functions, J. Chem. Phys. 84, 5728 (1986) electronic version.
6. ↑ ^{6.0 6.1} T. Saue, K. Fægri, T. Helgaker, and O. Gropen, Principles of direct 4-component relativistic SCF: Application to cesium auride, Mol. Phys. 91, 937 (1997) electronic version.
7. ↑ W. Kutzelnigg, Chemical Bonding in the Higher Main Group Elements, Angew. Chem. Int. Ed. Engl. 23, 272 (1984) electronic version.