On the construction of double group molecular symmetry functions

26 April 1996

CHEMICAL PHYSICS LETTERS ELSEVIER

Chemical Physics Letters 253 (1996) 20-26

On the construction of double group molecular symmetry functions Lucas Visscher Department of Chemical Physics, University of Groningen, Nijenborgh 4, 9747AG Groningen, The Netherlands Received 8 November 1995; in final form 13 February 1996

Abstract

A new procedure for constructing double group symmetry functions is presented. Using this method integrals over Hermitian operators can become real quantities, even though the integrand and the functions themselves are complex. This is especially of interest to 4-component relativistic methods that use the Dirac-Coulomb Hamiltonian directly in electronic structure calculations. I. Introduction

The treatment of relativistic effects on molecular properties is possible by several different methods that are all derived from the relativistic DiracCoulomb-Breit Hamiltonian. The last decade has seen the rapid development of 4-component relativistic methods that attempt to use this Hamiltonian directly in molecular calculations [1-4]. This is in contrast to the conventional method of first transforming to an approximate and less computationally demanding Hamiltonian. Though the increase of computing power is continuing, these 4-component calculations still demand rather formidable resources. It is therefore imperative to use the symmetry of the Hamiltonian to make the calculations feasible. This can be done using double group theory. In this work double point groups are indicated by asterisks, while the extra 'fermion' representations are named using capitals to distinguish them from the ordinary 'boson' representations of the single point groups. Double point group symmetry functions have been given for many point groups [5-13]. Wigner [14] Elsevier Science B.V. All rights reserved PII S 0 0 0 9 - 2 6 1 4 ( 9 6 ) 0 0 2 3 4 - 5

classifies the double groups in three classes according to the distribution of Kramers pairs [15] over the representations. We will be concerned with the first class in which a Kramers pair spans different rows of a multiple-dimensional irrep. It is well-known that by a judicious choice of phase factors or by applying a suitable unitary transformation, representation matrices of this class of double groups can be chosen to be real. Esser [16] noted that this could also be used to work with real arithmetic in his relativistic direct MRCI method. He did not implement the technique in order to have one code applicable to all double groups. Real CI matrices were employed in the double group CI code of Pitzer and Winter [1738] for the special cases of C~v, D~ and D~h. They consider the symmetry character of an effective spin-orbit operator f(r) 1.s in a two-component quasi-relativistic scheme. In this case the Hamiltonian matrix can be chosen to be real by inserting phase factors in the many-electron wavefunctions. An alternative approach is to use quaternion algebra [19] to block diagonalize the Hamiltonian matrix [20]. Recently, a new method was proposed by Saue

L. Visscher / Chemical Physics Letters 253 (1996) 20-26

[21] and Jensen. They used quaternion algebra to treat the so-called binary groups (D2h and subgroups). In a first application [22] they considered the DiracFock method and showed that this matrix could be made real for the double groups mentioned above. The other binary double groups could also be treated within the same framework, but the matrix then became complex or quaternion. They wrote symmetry functions as combinations of boson symmetry functions and inserted quaternion phase factors when constructing the Fock matrix, giving the matrix the desired blocked structure. This method would especially be beneficial in cases with low molecular point group symmetry (C 1 or C i*), because it fully employs Kramers symmetry. We propose a scheme similar to that of Pitzer and Winter, but introduce the phase factors already in the l-electron spinors. This guarantees that all one- and two-electron integrals that appear in the equations will be real for double groups of Wigner's first class. The advantage of such a scheme is that no special derivations are necessary for different methods. Computational savings can be obtained by simply changing the arithmetic in the time-consuming routines from complex to real arithmetic. This is easily done when employing standard libraries for such work. Many schemes express the functions in terms of atomic symmetry functions and use a predefined set of representation matrices [23] in their construction. The use of representation matrices is convenient, but lowers the flexibility in the choice of basis functions. We employ the group chain method proposed by Aerts [24]. This method produces symmetry functions that transform according to the double group under consideration and to one or more subgroups thereof. This subgroup can be chosen in different ways, giving rise to different representation matrices. For instance, in the octahedral group T* one can choose as subgroups either C 3 or C 2 . As we will see later, both these choices have their merits for this double group. The scheme proposed in this Letter only needs the character tables of the groups and can start from arbitrary non-symmetry adapted basis functions. We will show that for the double groups C2"~ and D 2 , and groups of which these are a subgroup, the symmetry functions obtained give real-valued integrals

21

over totally symmetric hermitian operators. This leads to a 2-fold reduction in storage of the relevant quantities and to a 4-fold reduction in floating point operations relative to schemes that give complex-valued integrals.

2. Construction of the symmetry functions for D~ ^ The double group O~ consistsofthe elements E, C z, C x, Cy, R, CzR, CxR and CyR, where R is a rotation through 2 ~r. It has four boson irreps a, b l, b 2 and b 3 and one 2-dimensional fermion irrep E. Symmetry functions for the fermion irrep may be obtained with the character projection operator =_

_

(1)

Since there is only one ferrnion irrep, this projection operator only has an effect if a boson type function is present in the trial function. We define orthogonal basis functions for the two rows of the fermion irrep by considering the subgroup C~. For this purpose the choice of the 2-fold axis in C 2 is arbitrary. We take the rotation along the z-axis. This^group then consists of the elements E, Cz, ~ and CzR, and is Abelian. It has two boson irreps, a and b, and two fermion irreps, E' and E". We again use character projection operators ~E:E

~ E f f E = p^'E__ 1__~ ( E ^ -R.+

iC~^ - if~(~.)

(2)

~ ' E " E = p E ~ " E = f f ' E = ¼ ( E - - f ~ - - i e z +iRt~z)

(3)

=

to project out the symmetry functions XE'E and XiE"E from a trial function Xr After normalization each pair of functions X E:E and Xff"E can be used to define matrix representations D~(O) of the operations () of the group D~. These matrix representations will, in general, differ in the phase of the off-diagonal elements. We can conform to a standard choice [25], for instance DEyE(Cx)=i and DE,,E(Cr)= 1, by multiplying X E''E by the appropriate phase factor. This procedure has been reported before [4,24] and leads to important savings in the double-group symmetryadapted relativistic calculations. It does not assure, however, that the integrals over such functions have real values. t

^

t

^

22

L. Visscher / Chemical Physics Letters 253 (1996) 20-26

We can reach that goal by putting an extra constraint on the basis functions. Consider the Kramers operator I( I(= -i

02

O'y] 0 /

(4)

where I( o is the complex conjugation operator. This anti-unitary operator commutes with the Dirac Hamiltonian and transforms a function of row 'E of representation E into a function of row "E. We now require that the basis functions form Kramers pairs,

f x?

=

x,

and

f x?

=

-x?

(5)

This constraint defines the phase of all basis functions without ambiguity, once an initial choice is made for one off-diagonal element of a representation matrix of an operation that is not contained in the subgroup C~. If we now consider integrals over such functions we will find that their value will be real. Take the integral over a totally symmetric hermitian one-electron operator h and operating with (~y (~y(E, 'E, i lh IE, 'E, i)

'E,"E =

which means that the integral must be real. For operators involving more electrons the proof is similar.

3. Definition of the symmetry functions for other double groups For C2"~ the construction of symmetry functions is done by the same procedure as derived for D~, choosing again C~ as a subgroup. We now need the off-diagonal element of the matrix representation of either d'v or ~v, to define the phases of the basis functions. Higher groups that contain either D~ or C2"~ are treated using character projection operators that are a product of two character projection operators, one for the highest point group irrep, and one for the highest Abelian subgroup irrep to resolve the degeneracy of the highest irrep. The group chains used in this procedure are given in Table 1, along with the correspondence of representations. For groups with a principal rotation axis of even order, the functions obtained also define representation matrices for C2"~

Y'~ (E, F t , i / h l E , F 2 , i )

r'l,F 2 X D/~*l:E (C y)DF2, E(C y ) -- (E, "E, i lh IE, "E, i)

(6)

t~y(E, 'E, i lh IE, 'E, i) = (E, 'E, i lh IE, 'E, i)

(7) and taking the same function and operating with the Kramers operator ~ ( E , 'E, i lh IE, 'E, i) = (E, "E, i lh* IE, "E, i) = (E,"E, i lh IE,"E, i)

Table 1 Group chains used in the construction of symmetry functions. The same chains are used for groups that are a direct product of the listed groups with the inversion group, for instance (D2*h, C~h) Point group

Irreducible representation

Abelian subgroup

Irreducible representation

D~ C2"~ C~

E E El E2 E1 E2 E1 E2 El E2 U El E2 U E 'U "U E U

C~ C~ C~

'E,"E 'E,"E 'El, "El 'E2, "E2 'El, "El ' E2, "E2 'EI,"EI ' E2, "E2 'El, "El ' E2, "E2 'El, "El, 'E2,"E2 'EI,"EI ' E2, "E2 'EI,"EI,'E2,"E2 'E,"E B,'E B,"E 'E,"E 'E,"E

I)4"

(8)

D~ Td"

I((E, 'E, i lh ]E, 'E, i) = (E, 'E, i lh IE, 'E, i)* (9) In Eqs. (6) and (8) we choose to let the operator work on the labels of the integral, while in Eqs. (7) and (9) it operates on the value of the integral. Combining formulae (6) to (9) gives (E, 'E, i lh IE, 'E, i) = (E, 'E, i lh IE, 'E, i)*

O*

T*

T*

C~ S~ S~

C~

C;

C2"

23

L. Visscher / Chemical Physics Letters 253 (1996) 20-26

or D~. In this case we simply require that all these representations be identical. This again defines the phase factors and makes the integrals real. For groups with an odd principal rotation axis, like T *, the choice of the highest Abelian subgroup (C 3 ) resolves the twofold degeneracy of the 'U and "U representations. It does not give real integrals because the symmetry functions do not give a proper representation for D~. Symmetry functions for D~ can only be obtained by combining rows of the irreps 'U and "U, making a reducible representation U. This destroys some of the symmetry blocking but provides real integrals. Which choice is more advantageous depends on the particular implementation of the methods that use these symmetry functions.

4. Implementation Symmetry functions obeying the two requirements given above, i.e. giving consistent representation matrices for C2"~ or D~ and forming proper Kramers pairs, can be constructed in a simple procedure. We illustrate this procedure for the F 2 molecule and give an indication of the speed-up that is achieved by this technique. Let us first consider an F 2 molecule in a twospinor basis of only s- and p-type functions. We use

the point group D4h. The four fermion representations of this group may be labeled as E | g , E2g, Elu and E2u. We use the Abelian group C4h to produce orthogonal symmetry functions, using the character projection operators where ( F , F ' ) = (Elg,'E,g), (E2g,'E2g), (Elu, 'E,u) or (EEl, 'E2~ ) . If we take Is l,/3 ) as a trial spin orbital we obtain the function Is I + s 2, /3) as a basis function for the first row ('Elg) of irrep E2g. The Kramers partner is found, by operating with K, as - Is I + s 2, a ) . This function can be used as a basis function for the second row ('Elg) of that irrep. The representation matrix of an operation of D4~, such as CEx, that is not contained in C ~h is given by D E " ( C 2 ~ ) = ( 0-i

O i)

We continue by projecting out the functions for other irreps, using the same trial function. This gives the Kramers pair is I -$2, /3) and - I s ~ - s 2, u ) for irrep E ~ , which defines a representation matrix

for C2~ as +i

O

Table 2 Symmetry functions for the D~ point group D4~ irrep

C ~h irrep

Large comp. functions

Small comp. functions

Eig

'Elg

[s t + s 2, /3) I P I : - P2z, fl) ]Ptx - P2x, a ) - i IPly - P2y. or) - Ist + s 2, a ) -- IPlz - P2z, or) I P l x - P 2 x , / 3 ) + ilp~y--P2y, /3)

ils I - s i lplz + i[plx + ils I - s i lplz + -ilPlx

ilst-s2,/3)

Is~ +s 2, 13)

i lplz + P2z, /3) i [ p t x + P2x, or) + IPty+ P2y, or) i ls I - s 2, a ) i l p l z + p 2 z , or) - i l P l ~ + P2x, / 3 ) + [Ply + P2y, /3) [ P l x - P2x, a ) + i l p l y - P2y, a ) IPt~ - P2x, /3) - i l p l y - P2y, fl) i lplx + P2x, or) - IPly + P2y, or) - i l P l ~ + P2x, / 3 ) - I P l y + P2r, /3)

] P l z - P2z, /3) I P l x - P2x, o r ) - i l P l y - P2r' a ) - I s I + $2. a ) - [Pl~- P2z, ct) I P l x - P2x, /3) + i [ P l y - P2y, /3) i lplx + P2x, 0~) - ]Ply + P2y, a ) - i [Pt~ + P2x, /3) - IPly 4- P2y, /3) IPlx - P2x, or) + i lPlr - P2y, a ) I P l x - P2x, / 3 ) - i l p l y - P2y, /3)

"Elg

Eju

'Et u

"Elu

E2g E2u

,~E2g Lag 'E2u "E2u

(12)

2, /3) P2z, ]3 ) P2x, or) + [Ply + P2y, or) 2, a ) P2z, or) + P2x. /3) + IPly + P2y, /3)

(13)

L. Visscher/ Chemical Physics Letters 253 (1996) 20-26

24

Table 3 Comparison of CPU times (in seconds on a Cray-J916) for DF calculations in different point groups. The test calculation concerns the F2 molecule in an uncontracted (L: 9s 5 p ldlS: 5s 10p 5d 1O scalar basis Point group

Scalar integrals

Fock matrix construction

Fock matrix diagonalizatton

ql'otal CPU time per DF iteration

C~* C~* C~ C2"v C~a D2*a C 2h D,~

13 804863 (100%) 6 905 571 (50%) 13 804863 (100%) 13 804 863 (100%) 6 905 571 (50%) 6 905 571 (50%) 3 645 286 (26%) 3 645 286 (26%)

58 (100%) 29 (51%) 59 (103%) 56 (98%) 28 (48%) 27 (47%) 15 (26%) 15 (26%)

95 (100%) 22 (23%) 10 (11%) 7.5 (7.9%) 4.3 (4.5%) 3.0 (3.2%) 2.1 (2.2%) 1.6 (1.6%)

156 (100%) 55 (35%) 72 (46%) 66 (42%) 34 (22%) 31 (20%) 18 (l 1%) 18 (11%)

This choice, however, would lead to an inconsistency in the matrix representations of the D 2 subgroup as both Elg and Elu reduce to E in this subgroup. We therefore multiply the first E~u symmetry function with i and the second with - i . Since these phase factors are each other's complex conjugate the functions remain a proper Kramers pair, (i ISl-S 2, /3), i l s l - s 2 , or)). The matrix representation is, however, changed by a factor of - 1 for the off-diagonal elements and equal to the one for Elg that is given in Eq. (12). We continue this procedure for the p-functions, multiplying the Kramers pairs obtained by projection by a phase factor to adhere to the previously defined matrix representations. The resulting functions are given in the third column of Table 2. If we work with 4-component spinors, the large component basis functions are simply the functions given in Table 1, where a stands for the first component of the 4-spinor,/3 for the second component and the two lower (small) components are zero.

The small component functions have opposite parity, which leads to a factor of i relative to the large component functions when we make their representation matrices identical. The small component s- and p-type functions are given in the last column of Table 2, with ot and /3 now the third and fourth component of the 4-spinor and the upper (large) components zero. The simple procedure described above does not always work when the E representation of C~v or D~ is spanned more than once by a group of symmetryrelated trial functions. In this case, projection of a trial function to the ' E row of the irrep, followed by operation with I( to get the partner in row "E, will not necessarily give a basis for that representation. If, instead o f apiflying K, we generate the partner by operating with C2x (in the case of D~ ) or 6"v (in the case of C2"~), we get a proper representation 'E,"E Oxi E'rl = E DE r~,r, ( 6 ) XiE'r: (14) 1"2

Table 4 Comparison of CPU times (in seconds on a Cray-J916) for CCSD(T) calculations in different point groups. The test calculation concerns the F2 molecule in the molecular 4-spinor basis (116 functions). The ls spinors were frozen Point 4-spinor integrals CPU time (s) CPU time (s) group per CCSD it. (T) corrections C~" 31 560 165 (100%) 408 (100%) 2830 (100%) C~* 15 781 429 (50%) 118 (29%) 812 (29%) C; 15 781 479 (50%) 108 (27%) 765 (27%) C~v 15 781 479 (50%) 37 (9%) 307 (11%) C~h 7 892 111 (25%) 37 (9%) 298 (11%) D~h 7 892 111 (25%) 14 (3%) 130 (5%) C~a 4013 833 (13%) 11 (3%) 114 (4%) D,~ 4013 833 (13%) 6 (1%) 56 (2%)

L. Visscher/ Chemical Physics Letters 253 (1996) 20-26

but the Kramers pairs may be mixed. We find, in general, the following structure

J

J

(151 with U a unitary matrix. When U is not equal to the unit matrix we have to transform our initial symmetry functions XiE'E and XiE:E by U +t/2. This transformation does not alter the representation matrices, but provides a new set of functions that fulfils Eq. (5) and gives real integrals. After this choice of basis function no special procedures have to be followed to work with real algebra in the above-mentioned point groups. All hermitian matrices that are expressed in this basis will become real symmetric. Molecular spinors will be linear combinations of large and small component basis functions with real coefficients. They therefore automatically fulfil the same conditions as were imposed on the basis functions.

5. Sample calculation The scheme described above has been implemented in the MOLFDIR [4] suite of codes for 4-component Dirac-Fock (DF), restricted active space configuration interaction (RASCI) and coupled cluster (CCSD(T)) calculations. We have chosen the fluorine molecule as an illustrative example of the speed-ups that may be achieved due to the use of real arithmetic. At the DF level (Table 3) we see that the number of integrals is the same for the point groups Ci*, C2h and D2*h. The speed-up going from Ci* to C2h comes from the reduction of the diagonalization time. The Fock matrix in Ci* symmetry consists of two blocks of dimension 208, while it consists of four blocks of dimension 104 in C2"h and D~h symmetry. In Ci* these blocks are independent and need to be diagonalized using complex arithmetic, in C2"h and D2, only half of the blocks need to be diagonalized, using complex and real arithmetic respectively. This reduces the diagonalization time per iteration from

25

21.9 to 4.7 and 3.0 seconds, respectively. For larger runs the diagonalization time is usually unimportant because the n 4 process of constructing the Fock matrix dominates the time. This is already seen in the reduction in CPU time going from D2h to D4*h that is almost linear in the reduction of the number of integrals. Results for electron correlation, CCSD(T) type, calculations are presented in Table 4. In the CCSD iterations the speed-up going from C l* to C i* and from D~h to D4*h again comes from a reduction of the dimensions of the matrices. The speed-up going from C2'~ to C 2 , from Cah to D2h or from C4h to D4h, on the other hand, is entirely due to the use of real instead of complex arithmetic. In this case we can divide our operations into 3 types: (a) manipulation with integer pointers and other overheads, (b) integral and amplitude sorting routines and (c) matrix-matrix and matrix-vector multiplications. The first class of operations will not benefit from the change from complex to real arithmetic, while classes (b) and (c) will give rise to a speed-up by a factor of 2 and 4, respectively. Class (a) type operations do not depend on the number of basis functions (N), but only on the symmetry. Class (b) type operations scale at most as N 4. Class (c) type operations scale as N 6 for the CCSD step and as N 7 for the triples correction. For large runs class (c) should therefore dominate, and the speed-up should become close to 4. In this example we find factors of 3.0 (C 2 --* C~v), 2.6 (C2h ~ D2h) and 1.8 (C4h "-~ D4h ).

6. Conclusions A new method of constructing symmetry functions for 2- or 4-compent spinors is presented. The method is applicable to double groups that contain either C2"~ or D~ as a subgroup and gives real integrals and Fock matrices. As an example, calculations on the fluorine molecule are presented, using different point groups in the calculation. Use of real instead of complex arithmetic is found especially profitable in CCSD(T) calculations where a speed-up between 1.8 and 3.0 is observed.

26

L. Visscher / Chemical Physics Letters 253 (1996) 20-26

Acknowledgements I would like to thank Trond Saue and Ken Dyall for valuable discussions, and the former for providing me with a preprint of his paper on quatemion symmetry [18]. Thanks are due to an unknown referee for his constructive remarks. A Cray Research Grant is gratefully acknowledged.

References [1] [K.G. Dyall, K. Faegri Jr. and P.R. Taylor, in: The effects of relativity on atoms, molecules, and the solid state, eds. I.P. Grant, B. Gyorffi and S. Wilson (Plenum, New York, 1991) p. 167. [2] T. Saue, Ph.D. Thesis (University of Oslo, Oslo, 1996). [3] A.K. Mohanty and E. Clementi, Int. J. Quantum Chem. 39 (1991) 487. [4] L. Visscher, O. Visser, P.J.C. Aerts, H. Merenga and W.C. Nieuwpoort, Comp. Phys. Commun. 81 (1994) 120. [5] Y. Onodera and M. Okazi, J. Phys. Soc. Jpn. 21 (1966) 2400. [6] H.T. Toivonen and P. PyykkO, Int. J. Quantum Chem. 11 (1977) 695. [7] J. Oreg and G. Malli, J. Chem. Phys. 61 (1974) 4349. [8] J. Oreg and G. Malli, J. Chem. Phys. 65 (1976) 1755. [9] J. Oreg and G. Malli, J. Chem. Phys. 65 (1976) 1746. [10] J. Meyer, Int. J. Quantum Chem. 33 (1988) 445.

[11] J. Meyer, Int. J. Quantum Chem. 52 (1994) 1369. [12] J, Meyer, W.-D. Sepp, B. Fricke and A. Rosen, Comp. Phys. Commun. 21 (1989) 2559. [13] J.G. Snijders, Ph.D. Thesis (Free University of Amsterdam, Amsterdam, 1979). [14] E.P. Wigner, Group theory and its application to the quantum mechanics of atomic spectra (Academic, New York, 1959). [15] H.A. Kramers, Proc. Acad. Amsterdam 33 (1930) 959. [16] M. Esser, Int. J. Quantum Chem. 26 (1984) 313. [17] R.M. Pitzer and N.W. Winter, J. Phys. Chem. 92 (1988) 3061. [18] The algorithm described in reference 17 was further optimized by Morokuma et al. See: K. Morokuma, K. Yamashita and S. Yabushita, in Supercomputer algorithms for reactivity, dynamics and kinetics of small molecules, ed. A. Lagan~ (Kluwer, Dordrecht, 1989), p. 37. [19] S.L. Altmann, Rotations, quatemions and double groups (Clarendon, Oxford, 1986). [20] N. RSsch, Chem. Phys. 80 (1983) 1. [21] T. Sane, results presented at the ESF conference Relativistic quantum theory of many-electron systems, II Ciocco, Italy, 1995. [22] T. Sane and H.-J./~. Jensen, personal communication (1995). [23] Explicit representations of 32 point groups and of the linear and fivefold symmetry groups can be found in: P. Pyykk~ and H. Toivonen, Acta Acad. Aboensis Se. B 43 (1983) 1. [24] P.J.C. Aerts, Ph.D. Thesis (Groningen, The Netherlands, 1986). [25] M. Hamermesh, Group theory and its application to physical problems (Dover, New York, 1962).