Colinp: The COLUMBUS input facility

AO integral input section and compatibility aspects

The AO-integral input program provides a database of several frequently used basis sets. Standard segmented and generally contracted sets are available. Additional basis sets can be typed in directly. However, it is recommended to extend the data base by editing the files in the directory $COLUMBUS/source/iargos/basis. Respective changes are straightforward. The basis set file consists of a concatenated set of an arbitrary number of basis set entries:

Notation for a basis set entry:
1. Title line
2. Number_of_basis_set_blocks
3. the corresponding number of basis set blocks, with each block containing
3a. Number_of_GTOs_in_Block Quantum_Number_L+1 Number_of_contractions
3b. Exponent coefficients for contraction 1 to Number_of_contractions

Note that COLUMBUS employs symmetry in form of symmetry adapted basis functions (D_2h symmetry and subgroups) throughout the program. Thus, orbitals are usually characterized by a number labeling each irrep and an index number within each irrep.

There are three AO integral programs available:

* is not a stand-alone integral program but is integrated into SCF, MCSCF, CI and integral-transformation programs. Not available in version 5.9.

AO-direct calculations circumvent the need for full integral transformations in the MCSCF and CI calculations, thereby reducing the disk space requirements largely.

Note! The choice of the AO integral program has to be made at the very beginning of an input session in order to generate consistent input files. Switching between different integral programs is not permitted! No consistency check is made at execution time!

Hint! If you have already the geometry file in COLUMBUS format, the integral input can be performed in a faster way by using prepinp tool (See Utilities section).

SCF wave function definition

The input via colinp is internally implemented through the Fortran program makscf.x. RHF and selected ROHF calculations are supported. In general, however, you must supply the alpha and beta values for ROHF calculations.

MCSCF wave function definition

The input via colinp is implemented through the Perl script makmc.perl. The major input, which has to be made, is the definition of the distinct row table (DRT). From a practical point of input generation it is sufficient to know that the DRT characterizes the configuration space. Thus, number of electrons, spin multiplicity, spatial symmetry of the state and orbital occupation scheme has to be given. It should be noted that the orbital scheme has to be specified for each irreducible representation separately. In case of excited state calculations several it will be frequently necessary to define several DRTs, each one characterizing one class of states. For more information on the DRT click here.

Hint! We strongly recommend that you fill out the Orbital Occupation Table given in Utilities section before starting a new job. It will help to provide the needed information during the input.

The first decision to be made is whether to run a single-point calculation or a geometry optimization at MCSCF or CI level. Provided there is an old MCSCF wave function definition file (mcdrtin) of a previous run in the current directory the DRT input may be skipped. This refers to the situation to rerun a calculation with identical configuration spaces but with different number of roots etc.

No more than 8 different DRTs are supported. Subsequently the number of electrons, state symmetry, state multiplicity, (cas,ras)->aux and and ras->(cas,aux) excitation level are to be entered on a per DRT basis. (cas,ras)->aux excitation level indicates the maximum number of electrons allowed in the aux orbital space. ras->(cas,aux) excitation level indicates the maximum number of holes in the ras space. Enter zero to disable the RAS and AUX orbital space usage, respectively.

All DRTs share the same orbital scheme, i.e. different configuration spaces are solely enforced through the number of electrons, state symmetry, multiplicity and excitation level plus optionally by enforcing different group occupation restrictions. In order to assign the orbitals to the various classes, the number of orbitals per irreducible representation are entered assuming the following ordering of orbitals DOCC - RAS - CAS - AUX per irreducible representation. Note, that the orbitals are labeled through their symmetry plus an index starting with the occupied orbitals and ending with the virtual orbitals.

An option follows to enter orbital group restrictions. The list of all orbitals belonging to the RAS, CAS and AUX classes in terms of labels composed of index and symmetry is displayed. In order to apply orbital group restrictions the listed orbitals must be regrouped by assigning them to a group. This works by entering a letter (a-z, A-Z) on the line below the corresponding orbital label followed by a question for the number of electrons in the respective groups. Up to 50 sets of orbital group restrictions, each composed of up to 20 orbital groups, are supported. From the original set of configurations fulfilling the DOCC, RAS, CAS, AUX restrictions all those are finally selected that additionally fulfill at least one set of orbital group restrictions. In this way almost arbitrary configuration spaces may be selected. This is best illustrated by an example:

             a1 b1 b2 a2
       DOCC   2  0  0  0
       RAS    1  0  0  0
       CAS    0  2  0  0
       AUX    1  1  1  1
     

       3a1 1b1 2b1 4a1 3b1 1b2 1a2 
   

 3a1 1b1 2b1 4a1 3b1 1b2 1a2 
  A   A   A   B   B   B   B
  occupation of group A : 4 group B: 0 
   

 3a1 1b1 2b1 4a1 3b1 1b2 1a2      
  A   B   C   D   D   D   D
  occupation of group A : 2 group B: 1 group C: 0 group D: 1 
   

This concludes the MCSCF configuration selection scheme. The MCSCF DRT program itself (mcdrt.x) allows for even more flexibility in configuration space selection than this input maker but is a little more tedious to use. For most calculations the input maker's capabilities are sufficient.

MCSCF energy optimization input

Next step is the input for the MCSCF energy optimization procedure (defaults are usually reasonable):

• Number of MCSCF iterations (#iter)
• Number of micro iterations per MCSCF iteration solving the PSCI equations (#miter)
• Number of subspace iterations per MCSCF iteration solving for the the eigenvalues of the HMC matrix (#ciitr)
• Orbital correction vector convergence threshold (knorm)
• State correction vector convergence threshold (wnorm)
• Energy lowering convergence threshold (DE)
• HMC matrix diagonalization: valid choices are
  
  • build explicitly, diagonalize direct or iteratively (Davidson method). For modest sizes of configuration space only.
    
  • do not build explicitly, diagonalize iteratively
• Include orbital-state coupling (quadratic cnvg.) starting from iteration number #iter. It is sometimes useful to suppress orbital-state coupling for the first 10 iteration, in particular for bad starting orbitals and large configuration spaces.
• Subspace resolution for each subspace CAS, RAS and AUX. The choices are no resolution (none), natural orbital resolution (NO), FAA or QAA Fock matrices (FAA, QAA). The possibility to resolve redundant orbital rotations (see also the discussion of fciorb in the full documentation) is important if well-defined orbitals are required in follow-up steps. This resolution is recommended if one wants to perform a CI calculation based on MCSCF orbitals and it is required if an MR-CI gradient calculation is to be performed (see the discussion on the analytic gradient in the methods section). Note that applying orbital group restrictions will in general partly or completely remove redundancies from the MCSCF wave function. Hence the correct orbital resolutions have to be inserted in this case by hand into the mcscfin file. Incorrect orbital resolution might result in convergence problems. Refer to the full MCSCF documentation.
• State-averaging: the number of states to be calculated per DRT. Only equal weights are allowed for states within one DRT.

MR-CI Wave Function Definition

There is a choice between three different CI programs: the standard CI program (ciudg), recommended for most applications, the AO-direct version thereof (ciudg_d) (not available in version 5.9) and the parallel CI program (pciudg). For an overview of available features see the CI documentation. For running the parallel CI program see the instructions given in the section on the parallel CI.

Hint! We strongly recommend that you fill out the Orbital Occupation Table given in Utilities section before starting a new job. It will help to provide the needed information during the input.

The purpose of this step is to compute the MR-CI DRT based on a list of internal orbitals and a set of reference configurations. In contrast to earlier implementations, a modified orbital group selection scheme (see below) also allows to add individual reference configurations without resorting to the internal step vector representation of configurations.

MOs containing frozen orbital indices are not written to the MO integral file.

In principle, the CI calculation is completely independent from a preceding SCF or MCSCF calculation as concerns the number of electrons, the spin multiplicity and the symmetry of the state. In practice, however, most of the calculations will be closely coupled to the structure of a SCF or MCSCF calculation.

For the calculation of transition moments in general it is necessary to calculate several CI vectors, which may belong to different irreducible representations. Hence several DRTs must be defined subject to the following restrictions:

• same number of electrons
• same multiplicity
• same internal orbital space
• same set of frozen orbitals

Up to 8 different DRTs are supported. Subsequently, the number of electrons, state symmetry, state multiplicity, the maximum number of electrons in the AUX space are to be entered on a per DRT basis. Enter zero to disable the AUX orbital space usage.

All DRTs share the same set of internal orbitals. Different configuration spaces are solely enforced through the state symmetry and reference space restrictions (through the number of reference doubly occupied orbitals plus optionally by enforcing different group occupation restrictions). Orbitals are assigned to the various classes by entering the number of orbitals per irreducible representation while assuming the counting scheme FCORE - DOCC - ACTIVE - AUX - EXT - FVIRT per irreducible representation. Note, that the orbitals are labeled by their symmetry plus an index starting with the occupied orbitals and ending with the virtual orbitals.

Invoking the CI input option you may choose among six major modes of operation:

• Taking over data from a preceding MCSCF input
• Taking over data from a preceding MCSCF input and preparing for geometry optimization
• Independent CI input in the one-DRT mode
• Independent CI input in the one-DRT mode preparing for geometry optimization
• Independent CI input in the multi-DRT mode

MCSCF related options are available only if the preceding MCSCF input did not make use of the state-averaging feature including more than one DRT

Subsequently, you need to decide how to construct the DRT. The DRT is separated into the external and internal orbital part. Frozen core or frozen virtual orbitals are removed at the AO-MO integral transformation step and do not show up in the DRT. The DRT program writes the configuration space information for the internal part only, whereas the external part is constructed implicitly in the diagonalization step. Each DRT level is associated with an molecular orbital, schematically:

       level          orbitals           BOTTOM mode      TOP mode

         n          internal orbital n   auxiliary orb.   REFDOCC orb.
         ..           ..                 active orb.      active  orb.
         1          internal orbital 1   REFDOCC orb.     auxiliary orb.
         0          external orbital     ext. orb.        ext. orb.
    

Usually it is expected that DRTs built in TOP mode form a more compact DRT representation of the configuration space. However, generalized space restrictions in TOP mode are not so efficient. As a rule of the thumb, choose TOP mode if you have a large number of REFDOCC orbitals and relatively few active/auxiliary orbitals. Otherwise prefer the BOTTOM mode (default).

MCSCF data based inputs

MCSCF data based inputs have many options in common. Hence, they are treated together. The state symmetry, number of electrons and spin multiplicity along with the basis set information and the orbital classification is immediately read from the MCSCF DRT definition file. Only the orbital classification with respect to doubly occupied and active orbitals (joined CAS and AUX spaces) can be extracted. Now you can assign and optionally split the orbital spaces to the orbital classification used in the CI DRT definition. Only orbitals located at the edges of each subspace can be moved in blocks. The shift of individual orbitals is not possible in that way.

• shift the lowest orbitals from the list of MCSCF doubly-occupied orbitals to the frozen core orbitals
• shift the highest unoccupied orbitals to the frozen virtual orbitals
• initially all remaining MCSCF doubly occupied orbitals are shifted to the ref. DOCC orbital group, all MCSCF active orbitals are put into the CI active orbital group and all remaining not frozen orbitals belong to the external orbital group.
• you may now shift additionally some of the lowest MCSCF active orbitals to the ref. DOCC orbital group
• you may shift additionally some of the highest MCSCF active orbitals to the CI auxiliary orbital group and enter the maximum electron occupation for the CI AUX orbitals.
• alternatively some of the highest MCSCF active orbitals can be shifted to the CI external or virtual orbitals.

Enter the excitation level, i.e. reference configuration CI, single-excitations included or all single and double excitation included. Generalized interacting space restrictions are available with MR-SDCI, only. This leaves you with the following summary information on screen:

       MCDRT input data:
       spatial symmetry: 1  mult:1 #el:16
       count order (bottom to top): fc-docc-active-aux-extern-fv
       irreps         ag    b3u   b2u   b1g   b1u   b2g   b3g   au
       MCSCF docc     3     1     0     0     2     0     0     0
       MCSCF active   0     2     0     0     0     0     2     0
       # basis fcts   11    7     4     2     11    7     4     2

       frozen core    1     0     0     0     0     0     0     0
       frozen virt    0     0     0     0     0     0     0     0
       external       10    5     4     2     11    7     2     2
       ref docc       2     2     0     0     2     0     0     0
       ci active      0     1     0     0     0     0     1     0
       ci auxiliary   0     0     0     0     0     0     1     0
                                        
       exc.level:2 gen.space:y #aux.elec:1
       Enter the allowed reference symmetries  1

The configuration space selection scheme allows you to include not only reference configurations matching the state symmetry but also those of different symmetry. Hence, you generate a much larger set of configurations out of which of course only those with the proper symmetry are selected. This feature comes in handy if you would like to match energy surfaces calculated with different point group symmetry. Also this is necessary to suppress symmetry breaking when doing calculations on systems belonging to point groups with multi-dimensional irreducible representations, e.g. to suppress the lifting of degeneracy of equivalent configurations in linear molecules.

Now you can optionally enter orbital group restrictions. You will be prompted with the list of all orbitals belonging to the CI active and CI auxiliary classes in terms of labels composed of counting index and symmetry. In order to apply orbital group restrictions the listed orbitals must regrouped by giving them the label of a group. This works by entering a letter (a-z, A-Z) on the line below the corresponding orbital label. Then you must enter the number of electrons for each orbital group. You may enter up to 50 sets of orbital group restrictions each composed of up to 20 orbital groups. From the original set of reference configurations that fulfill the CI active, and CI auxiliary space (and AUX occupation) restrictions all those are finally selected that additionally fulfill at least on set of orbital group restrictions. In this way almost arbitrary configuration spaces may be enforced. This is best illustrated with an example:

Assume the following orbital classes to be defined with AUX electron occupation set to 1, and 4 electrons in total:

                 a1 b1 b2 a2
       ci active  1  2  0  0
       ci aux.    1  1  1  1
     

This produces the label line in the group restrictions submenu:

       1a1 1b1 2b1 2a1 3b1 1b2 1a2
   

The orbital classes now e.g. include all single excitations from 1a₁ to the AUX orbital space. For some reason we want to retain only all active excitations plus those active-AUX excitations of type (1a₁)²(1b₁)¹(AUX)¹ in our reference configuration space. This may be achieved by the following two sets of orbital group restrictions:

orbital group restriction set 1 retains all active excitations

       1a1 1b1 2b1 4a1 3b1 1b2 1a2
        A   A   A   B   B   B   B
       occupation of group A : 4 group B: 0
   

orbital group restriction set 2 retains the specific active-aux excitations

       1a1 1b1 2b1 4a1 3b1 1b2 1a2
        A   B   C   D   D   D   D
       occupation of group A : 2 group B: 1 group C: 0 group D: 1
   

Based on these informations, one DRT is set up. Note that this DRT is not identical to the DRT created in the MCSCF step. The CI DRT is defined only for the internal orbitals since the external extensions to the CI matrix coupling coefficients (GUGA loop values) are created directly during the Davidson iteration process. Since the CI is restricted to singles and doubles there are four vertices at the bottom of the DRT connecting to the external space, which, for historical reasons, are denoted Z (no external excitations), Y (one-external), X (two-external, the two external electrons are coupled to triplet) and W (two-external, the two external electrons are coupled to singlet).

Independent CI input (one-DRT mode)

This mode allows the wavefunction definition fully independent of prior MCSCF calculations. Hence, one has to supply the full information manually

• number of electrons
• spin multiplicity
• spatial symmetry
• FCORE orbitals per irred. repr.
• FVIRT orbitals per irred. repr.
• internal orbital space, number of reference doubly occupied orbitals
• optionally, definition of AUX orbital space and maximum number of electrons in AUX space
• optionally, orbital group restrictions as explained above.
• optionally, multiple allowed reference configuration symmetries

generating the following summary output.

Multiplicity:1 #electrons:10 Molec. symm.:a1


count order (bottom to top): fc-docc-active-aux-extern-fv
irreps         a1    b1    b2    a2
# basis fcts   11    7     4     2

frozen core    1     0     0     0
frozen virt    1     0     0     0
internal       3     2     2     0
ref. docc.     1     1     0     0
ci active      2     1     2     0
ci auxiliary   0     0     0     0
external       6     5     2     2

exc.level:2 gen.space:y allowed ref. syms:1

Independent CI input (multi-DRT mode)

This mode allows the wavefunction definition fully independent of prior MCSCF calculations and allows for transition moment calculations. Supply the following data

• number of electrons
• spin multiplicity
• FCORE and FVIRT orbital space
• internal orbital space
• optionally, the auxiliary orbital space plus maximum number of electrons in here

count order (bottom to top): fcore-docc-active-aux-extern-fvirt
irreps         a1    b1    b2    a2
# basis fcts   11    7     4     2
frozen core    1     0     0     0
frozen virt    1     0     0     0
internal       4     3     2     0
ci active      4     3     2     0
ci auxiliary   0     0     0     0  
external       5     4     2     2


Input shared by ALL DRTs:
#DRTs:2 mult:1 #el:10 genspace:y

• spatial symmetry
• reference doubly occupied orbitals
• optionally, multiple allowed reference configuration symmetries
• optionally, orbital group restrictions

DRT specific input

irreps         a1    b1    b2    a2                            exc   sym   refsym
# basis fcts   11    7     4     2
ci active      4     3     2     0
ref docc DRT#1 1     0     0     0                             2     a1    1
ref docc DRT#2 1     0     0     0                             2     b1    2
ref docc DRT#3 1     0     0     0                             2     b2    3

In the standard case the symmetry of the reference CSFs is set to the symmetry of the state to be calculated. This choice is sufficient, when calculations (e.g. geometry optimizations) within a selected symmetry are performed. However, slightly different results will be obtained if for the same nuclear geometry calculations using different symmetry (e.g. C_2v and subgroups) are performed. The question of consistency and continuity between different symmetry choices is important when geometry changes from a higher to a lower symmetry are considered (as is the case for antisymmetric nuclear displacements). For this purpose it is necessary to include reference symmetries other than the state symmetry. Even though those references have the wrong symmetry and will not contribute to the wave function, excitations therefrom will have the correct symmetry. For example, if all reference symmetries are allowed, results will be identical to a calculations in C₁ symmetry.

Furthermore, the interacting space restriction (Bunge (1970), McLean and Liu (1973)) can be applied, which restricts the CSFs to those having a non-vanishing matrix element with one of the reference configurations. Interacting space restrictions work only in case of one reference symmetry.

The selection of all reference symmetries and dismissing the use of the interacting space restriction is important for obtaining smooth PES. This fact has to be considered, for example, in case of numerical calculation of force constants from analytic gradients.

DRTs are limited to at most 1027 rows and 2**20-1 internal walks.

CI diagonalization input for the standard CI program

The following screen-shot corresponds to the selection of the standard CI program ciudg.x:

   Type of calculation:               CI  [Y]  AQCC [N]  AQCC-LRT [N]
   LRT shift:                         LRTSHIFT [0                ]
   State(s) to be optimized           NROOT [1 ] ROOT TO FOLLOW [0]
   Reference space diagonalization    INCORE[Y] NITER [   ]
                          RTOL  [                                             ]
   Bk-procedure:                      NITER [1  ] MINSUB [1 ] MAXSUB [11]
                          RTOL  [1e-3                                         ]
   CI/AQCC procedure:                 NITER [20 ] MINSUB [1 ] MAXSUB [11]
                          RTOL  [1e-3                                         ]

                             FINISHED [ ]

• choice of the computational method (CI or AQCC or AQCC-LRT). Transition moment calculations are available with CI and AQCC-LRT, only. (translates to NTYPE and GSET entries in ciudgin)
• LRT shift applies to AQCC-LRT, only. Usually this is the energy difference between AQCC and reference configuration energy for the ground state. (translates to LRTSHIFT entry in ciudgin)
• Select the state to be optimized. There are two optimization schemes. (i) the first NROOT states are calculated (ii) the root with the maximum overlap with a reference state is calculated. (translates to NROOT and FROOT entries in ciudgin)
• Choose either incore or iterative reference space diagonalization. In the latter case choose a sufficient large number of iterations to converge the diagonalization to sufficient precision in terms of RTOL (defaults are mostly allright) (translates to IVMODE entry in ciudgin)
• In the Bk-procedure the Hamiltonian matrix elements not containing a reference configuration are discarded reducing the computational effort substantially. The first CI iteration is identical to the first Bk iteration if the trial vector contains only 0-external configurations (e.g. if you start from a reference diagonalization). Thus, carrying out exactly one Bk iteration saves time (translates to NBKITR, NBKMX, NBKMN, RTOLBK entries in ciudgin).
• Enter the number of CI iterations, minimum and maxiumum subspace dimensions and the required accuracy of the results in terms of residuum vector length (RTOL). (translates to NCIITR, NCIMX, NCIMN, RTOL entries in ciudgin)