Transition dipoles calculated by the determinantal CIS code are incorrect.

[serguei-patchkovskii]

Summary:

The transition dipole moments (and therefore the oscillator strengths) calculated by the determinantal CIS code are incorrect. For the ground-state to an excited-state transition, the dipoles are too small by exactly the factor SQRT(2.) For the transition dipoles between excited states, the dipole are too small by exactly the factor of 2.

Details:

As a very simple illustration of the problem, this input file:
he_cis.txt
produces the 1s² -> 1s2p transition matrix element:

 TRANSITION FROM THE GROUND STATE TO EXCITED STATE  6

 STATE MULTIPLICITIES =   1  1
 STATE ENERGIES =         -2.8611535740       -0.9710282657
 EXCITATION ENERGY =  1.2436E+16 [1/SEC] =   414834.55 [1/CM] =     51.43 [EV]
                            X           Y           Z           NORM
 TRANSITION DIPOLE =    0.229181    0.093584    0.510214    0.567098 E*BOHR
 TRANSITION DIPOLE =    0.582524    0.237869    1.296843    1.441430 DEBYE
 

It is a pretty weird thing to do, but one can also calculate this matrix element using the general GUGA CI code, using this input:
he_guga.txt
The result is (the state indices and the quantization axes for the 1s2p are obviously different, but the energies match exactly, so this is the same transition):

 CI STATE NUMBER=  1  3 STATE MULTIPLICITY=  1  1
 NUMBER OF CSF-S=        15        15
 STATE ENERGIES           -2.8611535740       -0.9710282657
 TRANSITION ENERGY=  1.2436E+16 [1/SEC] =   414825.06 [1/CM] =       51.43 [EV]
                          X [Z]       Y [C]       Z [C]       NORM
 CENTER OF MASS    =    0.000000    0.000000    0.000000             BOHR
 TRANSITION DIPOLE =   -0.000000    0.124976   -0.792200    0.801997 E*BOHR
 TRANSITION DIPOLE =   -0.000000    0.317660   -2.013587    2.038489 DEBYE

Now, 0.801997/0.567098=1.41421, which is the square root of 2.

A similar comparison between the DETS-CIS and the general GUGA CI for the 1s2s -> 1s2p transition gives, for the determinantal CIS code (the input files are the same):

 TRANSITION BETWEEN EXCITED STATES   2 AND  6

 STATE MULTIPLICITIES=  1  1
 STATE ENERGIES =         -1.6922551879       -0.9710282657
 TRANSITION ENERGY =  4.7454E+15 [1/SEC] =   158291.01 [1/CM] =       19.63 [EV]
                            X           Y           Z           NORM
 TRANSITION DIPOLE =    0.056429    0.023042    0.125625    0.139631 E*BOHR
 TRANSITION DIPOLE =    0.143430    0.058568    0.319310    0.354910 DEBYE
 

and for the general GUGA CI code:

 CI STATE NUMBER=  2  3 STATE MULTIPLICITY=  1  1
 NUMBER OF CSF-S=        15        15
 STATE ENERGIES           -1.6922551877       -0.9710282657
 TRANSITION ENERGY=  4.7454E+15 [1/SEC] =   158287.39 [1/CM] =       19.63 [EV]
                          X [Z]       Y [C]       Z [C]       NORM
 CENTER OF MASS    =    0.000000    0.000000    0.000000             BOHR
 TRANSITION DIPOLE =   -0.000000   -0.043518    0.275851    0.279262 E*BOHR
 TRANSITION DIPOLE =   -0.000000   -0.110612    0.701149    0.709820 DEBYE

This time, 0.279262/0.139631=2..

This example is sufficiently simple, so that one can also perform the calculation by hand, using the wavefunctions extracted from the ORMAS CI code (I do not believe the details and inputs are relevant, but are available if anybody cares). The results match the GUGA CI answer exactly - so I am pretty certain that the GUGA CI answer is correct, and the determinantal CIS code is wrong.

Exactly the same discrepancy is found in larger molecular calculations - except that these are becoming increasingly hard to run through the GUGA CI (again, inputs and outputs available on request).

As far as I can tell, the bug is ancient, and probably goes back to the original CIS implementation. I get exactly the same answers from both the current public version [30 SEP 2020 (R2)], going back at least to 2014, and probably before.