PHYSICAL REVIEW A 70, 012709 (2004)
Resonances for Coulombic potentials by complex scaling and free-reflection complex-absorbing potentials Shachar Klaiman, Ido Gilary, and Nimrod Moiseyev Department of Chemistry and Minerva Center for Nonlinear Physics of Complex Systems, Technion–Israel Institute of Technology, Haifa 32000, Israel (Received 16 March 2004; revised manuscript received 4 May 2004; published 21 July 2004) Analytical expressions for the resonances of the long-range potential (LRP), V共r兲 = a / r − b / r2, as a function of the Hamiltonian parameters were derived by Doolen a long time ago [Int. J. Quant. Chem. 14, 523 (1979)]. Here we show that converged numerical results are obtained by applying the shifted complex scaling and the smooth-exterior scaling (SES) methods rather than the usual complex coordinate method (i.e., complex scaling). The narrow and broad shape-type resonances are shown to be localized inside or over the potential barrier and not inside the potential well. Therefore, the resonances for Doolen LRP’s are not associated with the tunneling through the potential barrier as one might expect. The fact that the SES provides a universal reflection-free absorbing potential is, in particular, important in view of future applications. In particular, it is most convenient to calculate the molecular autoionizing resonances by adding one-electron complex absorbing potentials into the codes of the available quantum molecular electronic packages. DOI: 10.1103/PhysRevA.70.012709
PACS number(s): 03.65.Nk, 02.70.⫺c, 31.15.⫺p
I. INTRODUCTION
The method of complex scaling has long been implemented for calculation of resonance positions and lifetimes [1]. The calculation can be done directly (i.e., by solving for the eigenvalues of the complex scaled Hamiltonian) or indirectly by scattering theory (i.e., calculating the poles of the scattering matrix). The direct method entails using an analytical continuation of the Hamiltonian. Some of these continuations are complex scaling (CS), exterior scaling, and smooth exterior scaling (SES). Applications of these methods to Coulombic potentials are based on the fundamental mathematical grounds given by Balslev and Combes [2] and Simon [3]. Another way to calculate the resonances directly is by adding a complex absorbing potential (CAP) to the Hamiltonian [4,5]. The connection between the CAP and the SES methods was shown in Ref. [6] by looking for a reflection free FR-CAP. The approximated FR-CAP derived by Riss and Meyer [6] is based on a transformation of the kinetic energy operator (the so-called T-CAP). It is a point of interest that by taking the reverse approach, i.e., deriving FR-CAP by applying the SES method, a FR-CAP has been obtained [7]. Only when the CAP is introduced in a region where the potential energy has vanished is the T-CAP derived by Riss and Meyer [6] equal to the FR-CAP that has been derived analytically without any approximations by Moiseyev [7]. Note that in such a case although the T-CAP and the universal FR-CAP look slightly different due to the chain rule, identical differential operators may have a different appearance. It is quite easy to fulfill this condition (i.e., the product of the energy potential and the CAP operators is equal to zero) for short-range potentials. However, it is hard and often impossible to satisfy this condition for long-range potentials. The analytically derived FR-CAP has been shown to be useful not only in the calculations of resonance positions and widths but also in the motion of wave packet calculations [8]. Zavin and co-workers [8] have shown that by 1050-2947/2004/70(1)/012709(6)/$22.50
using the FR-CAP derived in Ref. [7] the reflections from the artificial “walls” imposed by the use of a finite-grid or basisset approximation are avoided in numerical calculations of the motion of wave packets. The resonances of Coulombic LRP (i.e., complex poles of the scattering matrix [9]) have been calculated by the complex coordinate method for atomic and molecular (neutral) Feshbach autoionizing resonances (see Ref. [1], and references therein) and to shape-type resonances of molecular anions [10,11]. Shape-type resonances have also been calculated by using different type of CAP’s. See, for example, the calculations of resonances of metastable dianions [12–14] by using CAP’s or of a model LRP as was done in Ref. [15]. Moreover, to the best of our knowledge, numerical exact calculations of resonance positions and widths by scattering theory (i.e., calculating the poles of the scattering matrix) for LRP were carried out only for LRP’s that decay more rapidly than 1 / R [16,17]. The direct approach usually requires the use of a large grid. This need causes severe numerical difficulties in the study of many particle problems and, in particular, when the electronic correlation effects are taken into consideration in the calculations. The main challenge then is to find a way to calculate these resonances directly without the need to scale the potential. The problem is that usually CAP’s are implemented in the region where the atomic or molecular potentials do not vanish numerically. It is almost impossible to satisfactorily introduce the CAP’s in regions where the potential vanishes in the study of autoionization of real neutral molecules where the asymptote of the ionized electron state is an hydrogenlike wave function (i.e., an electron in a Coulombic electric field). The results of our study (for the one electron analytically soluble LRP model Hamiltonian) presented here leads us to the conclusion that it might be possible to overcome this numerical difficulty by adding a new term to the FR-CAP [which has been derived from the smooth exterior scaling (SES) transformation] which takes into consideration the long-range Coulombic interaction effect.
70 012709-1
©2004 The American Physical Society
KLAIMAN, GILARY, AND MOISEYEV
PHYSICAL REVIEW A 70, 012709 (2004)
The results presented here, in addition to directing us to the type of CAP that should be developed in order to enable one to insert CAP’s at regions where the Coulombic interactions are still “on,” disclose the nature of the shape-type resonances for a model Hamiltonian that, to the best of our knowledge, is the only known analytically soluble problem. Moreover, based on our results a general conclusion is obtained concerning the derivation of analytical expressions for resonance positions and widths for any type of interaction (not necessarily for LRP’s): imposing of the squareintegrable boundary condition on the eigenfunctions of the Hamiltonian does not specify the contour of integration in the complex coordinate space. The only requirement from the chosen contour is that ⌿res → 0 as x → ⬁. It can be associated with CS, exterior scaling, or SES transformations. As shown for the Doolen model Hamiltonian, all of these transformations are not always applicable.
FIG. 1. CS and SCS.
II. COMPLEX SCALING, SMOOTH-EXTERIOR SCALING TRANSFORMATIONS, AND FREE-REFLECTION CAP’s
For the sake of clarity, a brief review of the methods of complex scaling transformations which enable the calculations of resonance positions and lifetimes will be given in this section. Let us consider the time-independent Schrödinger equation (TISE) although, in principle, the same arguments hold for the time-dependent Hamiltonian when the 共t , t⬘兲 method is implemented [18]: ˆ 共r兲res共r兲 = E res共r兲, H n n n
En = ⑀n − 共ı/2兲⌫n .
共1兲
Solving the TISE with an outgoing boundary condition provides resonance eigenfunctions (associated with the complex poles of the scattering matrix) which diverge exponentially res n 共r → ⬁兲 → ⬁. By using a complex scaled internal coordinate r → rei, the resonance eigenfunctions become square integrable provided the scaling parameter (usually defined as the rotational angle of the coordinates into the complex plane) exceeds a critical well-defined value. This type of scaling will be referred to as complex scaling (CS). It is, however, possible to use different scaling operators in order to transform the resonance wave functions to be square integrable (and thereby bring them into the generalized Hilbert space). As we will show later, sometimes the resonance wave functions become square integrable when applying a shifted version of the conventional scaling transformation so that i
r → r0 + 共r − r0兲e ,
autoionization phenomena. First, the one-electron molecular potentials are not dilation analytic potentials and in order to overcome this problem one should use complex scaled basis functions rather than a complex scaled Hamiltonian (see Refs. [10,1]). By scaling the molecular basis functions by a complex factor one shifts the center of the electronic Gaussian basis functions from its natural center associated with the position of the positively charged nuclei of the molecule in a given geometry. The price one pays for this kind of shifting and, in particular, in the calculations of the inner shell molecular orbitals, is in a slow convergence of the numerical calculations. The second disadvantage in using CS and SCS is in the fact that one cannot use the available commercial quantum chemistry packages for molecular electronic structure calculations without a major revision of the codes. In an attempt to avoid the need to scale the entire potential it was proposed [7] to add to the Hamiltonian a FR-CAP which has been rigorously derived by applying the smooth exterior scaling transformation. Under the specific conditions discussed below, by applying the FR-CAP one can avoid the need to scale the entire potential. Let us first consider the one-particle three-dimensional (3D) Hamiltonian with a spherical potential when the total angular momentum is equal to zero (i.e., can be reduced into a 1D problem), 2 2 ˆ 共兲 = − ប + V共兲. H 2M 2
By applying the SES transformation
共2兲
where 0 艋 r 艋 ⬁. We denote this kind of transformation as the shifted complex scaling (SCS) transformation (see Fig. 1). It is important to notice that by using the two methods mentioned above (CS and SCS) the entire potential is complex scaled and not only at the asymptote of the potential [assuming that V共r0兲 ⫽ 0] as in the SES method that will be described later, where is a function of x when 共x兲 → const as x → ⬁. The requirement of using complex scaled potentials as imposed by CS and SCS introduces a serious numerical difficulty in the study of many electron molecular
共3兲
= F共r兲 → rei
as
r → ⬁,
共4兲
one obtains the transformed Hamiltonian ˆ + Vˆ ˆ =H H f FR-CAP,
共5兲
ˆ is the original untransformed Hamiltonian and the where H FR-CAP is defined as 共f兲 VˆFR-CAP = ⌬V f 共r兲 + VˆCAP ,
where ⌬V f 共r兲 is defined as
012709-2
共6兲
PHYSICAL REVIEW A 70, 012709 (2004)
RESONANCES FOR COULOMBIC POTENTIALS BY…
⌬V f 共r兲 = V„F共r兲… − V共r兲
共7兲
共f兲 and VˆCAP is given by
2 共f兲 VˆCAP = V0共r兲 + V1共r兲 + V2共r兲 2 r r and V0共r兲 =
冉 冊
2 f ប2 −3 2 f 5ប2 −4 f 共r兲 2 − f 共r兲 4M r 8M r2
2
共9兲
,
ប2 −3 f f 共r兲 , M r
共10兲
ប2 关1 − f −2共r兲兴, 2M
共11兲
V1共r兲 =
V2共r兲 =
共8兲
FIG. 2. SFR-CAP using r⬘ = 15, r0 = 13. Note that the FR-CAP is simply the SFR-CAP, where r⬘ = 0.
where f=
F共r兲. r
共12兲
One can consider any contour F共r兲 (provided [F共r兲 → r exp共i兲 as r → ⬁] to be placed in Eq. (6) and calculate the resonances by solving the transformed time-independent Schrödinger equation. From now on we refer here to the FR-CAP [7] as the CAP constructed with the following contour:
冋
冉
1 cosh关共r − r0兲兴 F共r兲 = r + 共e − 1兲 r + ln 2 cosh关共r + r0兲兴 i
冊册
. 共13兲
A simple variation of this contour provides a shifted FR-CAP (SFR-CAP) which will be shown later to be the preferable one in our numerical calculations of resonances. The SFRCAP is defined by using (see Fig. 2)
冋
F共r兲 = r + 共ei − 1兲 r − r⬘ +
冉
1 cosh关共r − r⬘ − r0兲兴 ln 2 cosh关共r − r⬘ + r0兲兴
冊册
⌬V f 共r兲 (i.e., 1 / 关兩f共rជi兲 − f共rជ j兲兩兴 − 1 / 兩rជi − rជ j兩) unless one assumes that the electronic repulsion vanishes when the molecule is ionized.
III. TEST-CASE STUDY OF A 1D ANALYTICAL SOLUBLE PROBLEM
To the best of our knowledge the only known onedimensional LRP potential for which the resonances have been analytically calculated is V共r兲 = 共1 / r − ␥ / r2兲. The expression for the resonances derived by Doolen [19] is given by
Eres =
2m2 ប2
. ⫻
共14兲 One can see that the SFR-CAP is reduced to the FR-CAP when r⬘ = 0. Let us return to the discussion of Eq. (6) while defining as a perturbation the term ⌬V f 共r兲. If ⌬V f vanishes there is no need to consider the complex transformation of the potential 共f兲 . This condition is easy to satisfy when and VˆFR-CAP = VˆCAP short-range potentials (SRP’s) are considered. For piecewise potentials (SRP’s by definition) VSRP共r兲 = 0 when 兩r兩 ⬎ r0 for given values of r0 and therefore ⌬V f = 0. For more general SRP’s, ⌬V is not strictly equal to zero in the neighborhood of r = r0 and therefore we denote the case where ⌬V f = 0 as “FRCAP” (i.e., approximation to the FR-CAP). However, when dealing with LRP the situation isn’t so simple and now we must consider the case where ⌬V f = ⑀ where due to numerical difficulties (e.g., due to the desire to use as small a number of grid points or basis functions as possible) ⑀ cannot be taken as small as one wishes. It should be stressed that when solving a many body problem, it is hard to take into consideration the electronic repulsion in
8m␥/ប2 − 4n2 − 4n − 2 − i共4n + 2兲共8m␥/ប2 − 1兲1/2 . 共4n2 + 4n + 8m␥/ប2兲2 共15兲
Even though Doolen claims that CS is used in the solution of the resonance energies in Ref. [19], it is not the CS of the Hamiltonian that makes the resonance eigenfunctions square integrable but rather the use of the hyper confluent geometric function with a contour which has not been uniquely defined. As discussed by Moiseyev and Hirschfelder the only condition that this contour should satisfy is that its asymptote is rotated into the complex coordinate plane by the angle which has a sufficiently large value (as defined by Balslev and Combes). We will show later that resonances cannot be calculated in this case if CS (i.e., a rotated coordinate) is chosen as a contour. As has been shown before [1], there is a connection between the resonances of the Doolen potential and the bound states of the potential ˜V共r兲 = −V共r兲 that was given by Landau and Lifshitz [20]. This inverted potential supports bound states but no resonances. The eigenvalues for the bound states are given by
012709-3
PHYSICAL REVIEW A 70, 012709 (2004)
KLAIMAN, GILARY, AND MOISEYEV
TABLE I. Exact values (in a.u.) for the positions and widths of the two narrowest Doolen potential resonances as obtained from Eq. (15). Resonance position Eexact
Resonance width ⌫exact
First resonance 0.246093750000000 8.804240366863004E − 2 Second resonance 0.204152249134948 0.233967356461964
FIG. 3. Eigenstates of the Doolen potential V共r兲 = 4共1 / r − 4 / r2兲 as a function of the box size. The number of basis functions is kept constant, N = 400. In the inset the ground state is portrayed at a constant box size of 50 a.u. as a function of the number of basis functions N.
Ebound = −
2m2 关2n + 1 + 8 m␥/ប2兴1/2 . ប2
共16兲
As was shown in Ref. [1] by carrying out the transformation → − one can get the resonances for the Doolen potential from the expression derived by Landau and Lifshitz for the bound state eigenvalues of the inverted potential (and vice versa). Therefore, there is a one-to-one mapping between the bound states and resonance states of the two problems and we can get the energies of the resonances without specifically stating the contour which has been used to impose square integrable boundary conditions on the resonances solutions. The Doolen potential has two numerically “challenging” regions. It has a long-range Coulombic tail as r → ⬁ but also a strong singularity of −1 / r2 at the origin. The −1 / r2 potential term absorbed the scattered particles (see Ref. [9]) and therefore the spectrum of such a potential is unbounded from below as can be seen in Fig. 3. Since the variational energy is bounded from below by the real energy, as the box size increases the energy of the bound states increases due to the fact that the number of basis functions is kept constant. As the box size increases, more basis functions are required to describe this state accurately. The −1 / r2 potential term is such that the ground state that is located in the potential well has an energy of −⬁ (see inset in Fig. 3). In fact as the number of basis functions increases, more bound states appear which have the same tendency to −⬁. We chose for our calculations the parameters ␥ = 4 and = 4. From the analytically derived expression for the resonance positions and widths, given in Eq. (15), we calculated the values for the first two narrowest resonances (see Table I). Attempting to numerically calculate the resonance positions and widths by applying the CS method and use of the complex-variational principle [21] completely failed. This is not a technical problem that can be overcome by increasing
the numerical effort. This is due to the singularity of the Doolen potential at the origin as discussed above. Because of this singularity the particles escape from the top of the potential barrier also to the origin r = 0 where they are absorbed and not only to r → ⬁. Therefore we should avoid the use of a contour that does not scale the potential at r = 0 and does not impose absorbing boundary condition at the origin. Therefore the CS method is not applicable in this special case. However, using the SCS transformations, while taking r0 to be the distance of the top of the potential barrier from the origin (i.e., r0 = 2␥ / ␦), we have calculated the resonance positions and widths. The results obtained by the SCS method are presented in Fig. 4. Calculating the resonance wave functions (see, for example, Fig. 5) we observed that the resonances of the Doolen potential are not associated with the tunneling through the potential barrier as one might expect but are localized inside or over the potential barrier rather than inside a potential well, similarly to the behavior of the resonances of the Eckart potential barrier [1]. The location of the Doolen resonance wave functions at and over the potential barrier can be explained by recognizing the oneto-one mapping of the Doolen resonances with the bound state functions of the Landau-Lifshitz potential. As mentioned before by varying only one parameter in the potential the bound states (which are exponentially localized inside and over the potential well) in the Landau-Lifshitz potential are turned into resonances which are exponentially localized at the top of the Doolen’s potential barrier (i.e., the inverted potential well in Landau-Lifshitz problem). Since the reso-
FIG. 4. trajectory of the first resonance position and width relative error ⌬E = E / Eexact − 1 where Eexact is given in Table I. ⌬⌫ = ⌫ / ⌫exact − 1, where ⌫exact is given in Table I.
012709-4
PHYSICAL REVIEW A 70, 012709 (2004)
RESONANCES FOR COULOMBIC POTENTIALS BY…
FIG. 5. First (dot dashed) and second (dashed) resonance wave functions shown on a base line at their respective real part of the energy.
FIG. 7. First resonance wave function as r0 changes. Showing the exponential divergent nature of the scaled broad resonance solution.
nances are localized at the top of the Doolen potential barrier there is a decay of temporarily trapped particles into two opposite directions. The particles can escape in both directions r → ⬁ and r → 0, and the resonance wave functions exponentially diverge at these boundaries. This is the reason that in order to make the corresponding resonance wave functions square integrable functions and bring them back into the generalized Hilbert space, one has to scale the potential by a complex factor (and thereby impose absorbing boundary conditions) both as r → ⬁ and as r → 0. Applying the FR-CAP we encounter the same difficulty as in the CS because we do not scale the origin as r → 0. In order to scale the potential at both boundaries we used the SFR-CAP instead. We calculated the resonance energies and widths as a function of r0 while keeping the value of r⬘ at r0 − 2. As demonstrated in Fig. 6, the cusp (stationary solution in the complex variational space) has been obtained at
r0 = 2.6. This sensitivity to the value of r0 can be reduced by increasing the number of the basis functions or grid points. The resonance wave functions, as was discussed before, diverge exponentially as r → ⬁ and by rotating the asymptote into the complex plane we make them square integrable. The SES transformation only scales the Hamiltonian from a certain r0 and therefore the wave function decays only beyond the point where the contour exits the real axis and goes into the complex plane, i.e., for r 艌 r0. This is the explanation to the behavior of the resonances wave functions which are presented in Fig. 7. Our calculations also show very clearly that in order to get accurate values for the resonance positions and widths for the Doolen potential using a small number of grid points or basis functions (i.e., in the numerical calculations r is varied from 0 to L where 1 / L cannot be neglected and taken as equal to zero) one cannot avoid the scaling of the Coulomb potential and the ⌬V f 共r兲 term cannot be neglected. However, it is clear that when the Coulomb potential is scaled accurate results have been obtained. IV. THE BOUNDED PROBLEM WITH THE INVERTED DOOLEN POTENTIAL
The Landau and Lifshitz [20] potential mentioned above V共r兲 = 共␥ / r2 − 1 / r兲 has only bound states which can be calculated analytically, see Eq. (16). Using the same parameters as for the resonance problem, ␥ = 4 and = 4 (see Table II for the energies), we attempted to better understand the role of ⌬V f in LRP. TABLE II. Exact energies (in a.u.) for the ground and first excited bound states of the Landau and Lifshitz potential as obtained from the analytical expression given in Eq. (16). FIG. 6. Convergence of the SFR-CAP calculations by changing r0. The results presented are for the first resonance position and width relative error ⌬E = E / Eexact − 1, where Eexact is given in Table I. ⌬⌫ = ⌫ / ⌫exact − 1, ⌫exact is given in Table I.
Bound state energy First bound state Second bound state
012709-5
−0.209539778548435 −0.155229110778659
PHYSICAL REVIEW A 70, 012709 (2004)
KLAIMAN, GILARY, AND MOISEYEV
FIG. 8. r0 trajectory of the first bound state. FR-CAP refers to ⌬V f ⫽ 0 and “FR-CAP” refers to ⌬V f = 0. ⌬E = E / Eexact − 1 where Eexact is given in Table II. ⌫ refers to ⌫exact given in Table II.
Two possibilities were considered in our calculations. 共1兲 ⌬V f ⫽ 0 throughout the entire length of the box (i.e., 0 艋 r 艋 L) used in our numerical calculations. In our calculation L = 100. 共2兲 ⌬V f = 0 throughout the entire length of the box. The results as plotted in Fig. 8 clearly show the significance of the value of r0. One can choose ⌬V f = 0 only if one takes r0 to be far enough from the origin such that the variation of the potential (scaled or unscaled) no longer affects the temporarily trapped particle, i.e., where the resonance wave function almost vanishes. Exploring the effect of the numerical effort on the results as a function of r0 (see Fig. 9) one can see that for the case of the FR-CAP, any r0 can be used as long as the basis set is sufficiently large.
FIG. 9. r0 trajectory of the first bound state using the FR-CAP and increasing the number of basis functions used. ⌬E = E / Eexact − 1 where Eexact is given in Table II. ⌫ refers to ⌫exact given in Table II.
Our calculations of the resonances embedded in the Doolen potential and of the bound states in the Landau-
Lifshitz potential (both are analytically soluble eigenvalue problems) clearly show that the FR-CAP can be used for calculating the poles of the scattering matrix without the need to modify the LRP (e.g., scale the potential by a complex factor). This can be done if the FR-CAP is introduced in a region where the resonance or bound wave functions get exponentially small values. This result is important as it suggests that the FR-CAP can be improved when used in longrange calculations (i.e., we can introduce the cap closer to the interaction region of our system) by defining ⌬V f 共r兲 differently. Instead of the entire approximation noted above as “FR-CAP” we can define a better one by taking ⌬V f 共r兲 as 1 / F共r兲 − 1 / r. This implies that the potential reduces approximately to a Coulombic interaction in the region where the scaling is introduced. It is expected that this proposed CAP will be found to be useful in the calculations of molecular resonances of neutral molecules.
[1] N. Moiseyev, Phys. Rep. 302, 211 (1998). [2] E. Balslev and J. M. Combes, Commun. Math. Phys. 22, 280 (1971). [3] B. Simon, Commun. Math. Phys. 27, 1 (1992); Ann. Math. 97, 247 (1973). [4] G. Jolicard and E. J. Austin, Chem. Phys. 103, 295 (1986). [5] U. V. Riss and H.-D. Meyer, J. Phys. B 26, 4503 (1993). [6] U. V. Riss and H.-D. Meyer, J. Phys. B 31, 2279 (1998). [7] N. Moiseyev, J. Phys. B 31, 1431 (1998). [8] R. Zavin, I. Vorobeichik, and N. Moiseyev, Chem. Phys. Lett. 288, 413 (1998). [9] J. R. Taylor, Scattering Theory (Wiley, New York, 1972), Chap. 2. [10] C. T. Corcoran and N. Moiseyev, Phys. Rev. A 20, 814 (1979). [11] T. N. Rescigno, A. E. Orel, and C. W. McCurdy, J. Chem. Phys. 73, 6347 (1980).
[12] T. Sommerfeld, U. V. Riss, H.-D. Meyer, and L. S. Cederbaum, Phys. Rev. Lett. 79, 1237 (1997). [13] T. Sommerfeld, F. Tarantelli, H.-D. Meyer, and L. S. Cederbaum, J. Chem. Phys. 112, 6635 (2000). [14] T. Sommerfeld, J. Phys. Chem. A 104, 8806 (2000). [15] I. B. Müller, R. Santra, and L. S. Cederbaum, Int. J. Quantum Chem. 94, 75 (2003). [16] U. Peskin and N. Moiseyev, J. Chem. Phys. 97, 6443 (1992). [17] T. N. Rescigno, M. Baertschy, D. Byrum, and C. W. McCurdy, Phys. Rev. A 55, 4253 (1997). [18] U. Peskin and N. Moiseyev, J. Chem. Phys. 99, 4590 (1993). [19] G. D. Doolen, Int. J. Quantum Chem. 14, 523 (1978). [20] L. D. Landau and E. M. Lifshitz, Quantum Mechanics (Pergamon, Oxford, 1965), problem 3 in Sec. 36. [21] N. Moiseyev, P. R. Certain, and F. Weinhold, Int. J. Quantum Chem. 14, 727 (1978).
V. CONCLUDING REMARKS
012709-6
