- Email: [email protected]
The application of exterior complex scaling in calculations on resonances in nuclear motion in molecular systems
Chemical ?hysics 71(1982) 127-133 Nortix-HollandPublishingCompany
THE APPLICATION OF EXTERIOR CMPLEX SCAL%G IN CALCULATIONS ON RESONANCES IN NUCLEA’? MOTION IN MOLECULAR SYSTEMS Julia TURNER and C. William McCURDY * DepuRmenr of Che?r’Gssny,Ohio Srore University. Columbus, Ohio 43210. USA Received 18 January 1962
The application of the method of exterior coniplex scaling to find complex resonance energies by direct numerical integration is presented. Details ofthe numerical procedures required are &cussed and results are presented for rotational predissociation resonances in Hz. The advantage of the exterior-scaling approach is that it requires analytic continuation of the potential to complex values of the coordinates only in the asymptotic or near-asymptotic region. Therefore any repiesentation of the potential is adequate, including piecewise-contiuuous representations such as cubic splices. in contrast, application of the usual complex-scaling (rotated coordites) procedure generally requires that the potential be represented by a single analytic function.
1. Introduction The method of complex coordinates has found most of its applications in finite-basis-set calculations of resonance parameters for electronic resonances. (See for some early work on mathematical aspects ref. [I ] and for some applications in atomic resonance calculations ref. [2] _) The usual complexcoordinate transformation is the one appropriate for dilatation analyiic hamiltoinians, e.g. atomic systems, and consists of scaling the radial coordinates r_ of each electron by a complex scaling parameter, rj + ei@$.. For a range of the angle Q in the scaling parameter e@ this transformation converts the wavefunction whose
ener,gy eigenvalue is the complex resonance ener,y, and which increases exponentially at large real values Of 3, to a square integrable Function. Thus the complex-coordinate transformation allows the resonance wavefunction of the transformed problem to be expanded in a basis of square integrable functions. Consequently the complex resonance energy corresponding to a pole of the S-matrix can be computed (approximately) as the eigenvalue of a finite matrix. It
would seem that basis-set expansion is the ideal way to take advantage of the properties of the complex* Alfred P. Sloan Research Fellow. O30~-olo4/82/OOOO-OOao/S 02.75 0 1982 North-Holland
coordinate transformation [3]. Recently however, Atabek et al. [4] have demonstrated that the complex-coordinate transformation can also be applied in numerical integration of the Schrijdinger equation to fmd resonance solutions. For resonances in nuclear motion, for example molecular predissociation resonances, it is likely that direct numerical integration will be more efficient than basis-set methods for resonance calculations because of the short wavelengths involved. However, there is a considerable drawback in the use of complex coordinates for nuclear motion because the potential V(r ei@) is required at complex values of the nuclear coordinates. Unlike the case of electronic resonances in atoms where the potentials are all Coulomb potentials and hence .known analytically for all values of the electronic coordinates, the potential for nuclear motion is generally known oniy from electronicstructure calculations at discrete values of the nuclear coordinates. The usual procedure is to fit these discrete values with some functional form, but in most cases the functional form is only piecewise continuous. For example, for a diatomic molecule cubic splines or power series might be used for small values of the internuclear separation and reciprocal powers at large values. Although such a functional form might be necessary to represent the potential accurately for
J. Turner. C IV. McCwdy
i28
/ ExtenM complex scaiing in molecularsystems
real values of the nuclear coordinates, it does not defme 2 useful andytic continuation to complex values, 2nd such an andytic continuation is exactly wP;at is required in the complexcoordinate approach. In this paper we propose a solution to this difficulty by demonstrating that Simon’s [S] “exteriorscaling” theory can be applied to numerical calcul2tions on resonances in nuckar motion. Exterior sc2ling requires the use of complex coordinates only in ihe asymptotic or near-asymptotic region, 2nd thus does not require azdyticity of the potential except in that region. The piecewise-analytic functions necessary for accurate representations cf the potential are therefore easily accommodnted in this approach. The outline of this paper is 2s follows. Section 2 briefly describes Simon’s exterior scaling theory and section 3 presents tie details of its application in numerical integration of the Schriidinger equation. In section 4 we present the results of the application of this approach to fiiding rotational predissociation resonances in H7 and in section 5 we comment on possible extensions of the method to other types of problems.
7,. Exterior
scaling
Simon proposed his theory of exterior scaling to treat the problem of ?he definition of molecuIar electronic resonance energics in the Born-Oppenheimer approximation [S] _In that problem the potential which electrons experience due to the fLved nuclei is not dilatation 21Aytic aboct any center, 2nd therefore ordinary complex scaling fails to provide 2 method for finding resonance states. The basic idea of exterior scaling is that, if the non-analytic region of the potential can be enclosed in a sphere of radius R,, one c2n make coordinates complex only outside of that sphere and still obtain the Same transformation of the spectrum of the hamiltoni2n which ordinary complex scaling achieves in the dilatation analytic case. For simplicity let us specialize our discussion to the case of 2 single particle. If r denotes the coordinates of the particle, Sinon chooses the contourR(r), parameterized by the magnitude T of r
R(r) = i- ,
O,cra?,
=RofeiO(r--Ro),
~
R,,
(1)
shown in fig. 1, along which to continue analytic2Ily the hamiltonian. He also introduces the operator uR,,@) which, when it acts on a function of r, carries out the exterior-scaling transformation according to uRo(d)*(r)
=J(f@@(rjf)
,
(2)
where J(r) is the square root of the jacobian
_
J(r) = [dR(r)[dr]
‘12R(r)/r
The hamiltonian
is thus transformed
= uR,(@H@i;(@)
&,f@‘)
so that the hamiltonian fi = (-@/2J7z)Vz becomes, H&(Q)
(3) to (4)
> of our one-particle
system
+ V(r)
(9
under the exterior-scaling = (-fi2/k)V’
according
tranfiformation
+ v(r),
= (-fiz/2m)e-ziQVa
+ V(R(r)?),
OGrGRg, R, < r. (6)
Under *Jle assumptions that the non-analytic behavior in the potential occurs for r < R, Simon has shown that the spectrum of the tmnsformed hamilton& fiR,,(Q) Consists Of discrete bound st2tes 2t red energies, and continuum branch cuts (one in our one-particle case) rotated into the complex plane by an 2ngIe equal to --2Re(Q). In other words, the spectrum of the electronic hamiltonian of z molecule in the Born-Oppenheimer approximation is transformed by the exterior-scr&ng transformation in precisely the same w2y 2s the spectrum of an atomic h2miltorkn is transformed ucder the complex-scaling transformation.
r plane
Fig. 1. The comples r plane showbg
the radial coordinate under the complex-scali% transformation (dotted) and under the exterior-scaling transformation (solid).
J. l’hmer. C. IV. McCimiy /Exterior
Briefly stated, the reason the exterior-scaling transformation works in this way is that it is the same as ordinary complex scaling asymptotically, and it is the asymptotic behavior of the wavefunction which is ivportant in *&e determination of the spectrum. Of course, any hamiltonian which is dilatation analytic is also analytic under the exterior-scaling transformation. Interaction potentials for nuclear motion generally fall off at larger as fast or faster than r-1 and have non-analytic behavior at smaller values of the coordinates which are induced only by the functional fomrs with which we choose to represent them. In the following section we show how to use the transformed hamiltonian of eq. (6) in a numerical-integration appreach to finding its resonance eigenvalues and wavefunctions.
3. Numerical method
and a computational
example
Again we restrict the discussion to a one-particle system; extension to the case of several channels will be discussed briefly in section 5. If we consider a case for which the potential is local and central (these are not essential restrictions) we obtain, on passing to partial waves, an equation for r times the radial function, denoted by P(r), d%‘(r)/dr*
= [U(r) - e]P(r)
3
(7)
where U(r) = (2m/@) [ V(r) + j(j i- 1)/2vzr2] ,
j is the angular-momentum
W
quantum number, m is the mass and E = 2mE/fi2. To fmd bound states of eq. (7) by numerical inte,oration it is often convenient to employ Cooley’s algori’&m [6], an iterative procedure which involves outward integration to some point and inward integration to the same point. Cooley’s contribution was to fmd an efficient algorithm for computing the step D(E) to take in the search for the bound-state energy at which the tivo integrations produce the same logarithmic derivative at the matching point. The entire procedure is based on the Nurnerov [7] integrator. If Pi denotes the value ofP(r) at the ith integration step, the Numerov integration is carried out with the threepoint recursion formula
complex
scaling in moleadarsystems
Yi,i + Yj_i - 2Yj
124
=h’(Uj - E)Pj
I
(9)
where 5
= [1-
(@112)(~j
- f)]Pi )
(10)
h is the step size and Vi is V(r) at the ith step. The error is of order h6 in each step. The results of the inward and outward integrations for Y are scaled to make their values equal at the matching point. Cooley’s expression for D(E) which results,
01) where Y, is the value of Y at the matching point, leads to second-order convergence in AE, the distance from the bound-state energy. Note that the error in the numerator of eq. (1 i) is also of order h6. Atabek et al. [4] demonstrated that eqs. (9)-(11) can be used to find resonance states in a complexscaling calculation. If ordinary complex scaling, r -+reib, is employed, Cooley’s algorithm is essentially unchanged except for char@ng the step h to heio. The values ofr at each step are then just multiples of hei*. The error at each step of the integration is still of order h6. Unfortunately, the Cooley procedure apparently cannor be similarly generalized to the exterior-scaling contour without changing the order in h of the error accumulated at a critical point in the integration. The three-term recursion formula in eq. (9) has an error of order h6 only if the integration
steps are equal, but they must be equal in directio~z as well as in magnitude if the cancellations which cause terms of order 123through h5 to vanish are to occur. Therefore where the exterior-scaling contour R(r) changes direction at R. the recursion formula in eq. (9) accumulates an error of order h3_ However, this point is a crucial one in the integration
because generally it will
be where exponential
decay of P(r) begins. One of the tests of the numerical calculation of the resonance ener,v is the extent to which the result is independent of the scaling angle 9. The leading error in the integration step withRu at its center is (@/3!)eiO(ezie
- l)(d3P/dr3)R,
,
and can be responsible for sensitivity of the resonance energy to 4 when all other components of the calcula-
J. ZU?ner, C IV. MXkrdy
130
/Exterior
complex scalitig in mdecular
tion show negligible dependence on & This sensitivity to Q is zpparent if we choose R0 as the matching point in Coolmy’s algorithm, in which case eq. (11) becomes’
’
D(E) = [(eio + I)/21 {(urn - ~)p, + [-SOY,
gxtems
b-
\p(Ro)
_1 + (ei@ + i)Y,, - Yrr+l ]/S@eib(&+l)] _
+ 1 @(R(r))[dR(r)ldr]
dr])-l
.
04)
Ro
The numerator of this expression is just the generalization of eq. (9) (divided by h2) whenRo is the mth point. The error in D(E) is thus of order h3 and if the criterion fo: convergence to the resonance energy is that Ll(E) vanish, the error in the computed resonance ener,v is of order /23 even though the error in every other component of the calculation is of order h6. However, usingRO as the matching point suggests a way of avoiding this problem altoge’rher. Instead of Cooley’s algorithm for D(E) we can use the corrector due to Douglas and Ridley [S]
me derivatives in the numerator are evaluated at R. and P(r) and P(R(r)) are obtained by Numerov integration in the ouhvard and inward directions respectively and scaled to match at R, with the (arbitrary) value P(R,-J. In the inward integration the step h in eqs. (9) and (10) becomes he@‘_ To computeB{E) in eq. (14) so that it will vanish up to order h6 we merely need to compute the derivatives in the numerator to that order with six-point differentiation formulas using the results of the six integration steps adjacent toRo_ For example, for the outward derivative we use [9] [rlPold%,t
=(1/120h)(-24Pm_S
- 400 P,I_3
’
[pj&[P2(d d,!-’ ,
+ SOP,_4
+ 6C0 P,B_2 - 600 P,,_1 f 274 Pm) , 05)
03)
where ‘;n is the matching point, and Ph,(,,) and P;,(rrn) are the derivatives ofP(r) at the matching point from the outward and inward integrations respectively. This form forD(E) does not lead to quadratic convergence in the energ iteration, but we can cause it to vanish up to order h6 at convergence. The Douglas-Ridley corrector in eq. (13) can be generaiized for the‘ exterior sc&ng contour R(r), and ifR, is chosen as the matching point we obtain
* Derivation of eq. (12) requires that Pi and Yi be redefined to contain jacobian factors from eq. (3). Pi becomes e*/zpi for i > m md [@ t 1)/211” 1Dn2 at the matching point. With these changes and some simple manipulations, the kinetic-rnergy contribution to Cooley’s matrix F.4in his equsrion (3.3) czm be made symmetric.
with F, denoting P(R& and a similar formula obtains for the other derivative. We have found L3at it is adequate to use extended Simpson’s rule [9] to do the two integrals in the denominator of eq. (14). Then the numerator vanishes with error sixth order in 12at convergence while the denominator contains errors of order h5. We have found that eq. (12) and sirrar modikation of Cooley’s algorithm are inferior to our modification of the Douglas-Ridley co.qector and can lead to instability with respect to the scaling angle 9. Our more stable alternative, the numerical exterior-scaling algorithm, can be summarized as follows: (1) Choose a first g!JessE to the resonance ener,qy. (2) Integrate outward from P = 0 to r =R, with Numerov integration [eqs. (9) and (lo)]. (3) Starting with a value ofR(r) in the asymptotic region and the initial condition P(R(r)) = t$~&R(r)), where i,k is the Ricatti Hankel function (see for example ref. [lo]), integrate inward using eqs.(9) and (10) with h changed to heiS. The scaling angle Q must be large enough to expose the resonance
by rotating the
Jl Turner, C PJ.McCurdy /Exterior complex scaring in molecule? systems continuum
and the potential is extended dependence for r > 9 .S a0 ,
branch cut pest it;The choice ofRO is not critical but the least numerical sensitivity to its value occurs withR, chosen larger than the location of any classical turning points in r. (4) Scale the results of the integrations to match with the value P(Ra) at R, . Compute the numerical derivatives in eq. (14) with six-point formulas [eq. (15)] ami compute the integrals with extended Simpson’s rule. (5) Change the value of E to E +D(E) [eq. (14)] and return to step (2). In the following section we will show how this method can be used with an accurate representation of a potential which does not provide a useful analytic continuation to complex7 for small lrl. Note that no complex arithmetic isnecessary to compute U(r) for the outward integration.
V(r) = --c,/fi
4.A numericalexample Waech and Bernstein [ 111 have given a representation of the highly accurate KoTos and Wolniewicz potential [12] for Hz. For values ofr between 0.4 and 9.5 a0 they use a 33-term power series of the form (atomic units)
06)
- c,/r* - c&*
131 with reciprocal
.
power
(17)
The coefficients an, C', ,C8 and Cl0 are given in reference 1121 and references therein. There is a small discontinuity at r = 9.5 of 0.21 cm-l. Attempts to use ordinary complex scaling for this problem run into grave trouble immediately, because, although the interior and exterior parts of the potential match well at real I, for complex r there is no simple boundary in the complex plane along which (16) and (17) give the same value. This probfem is more acute at large values of Q because the contour rei@ can probe the edges of the region of convergence of eq. (IS) where oscillations occur. Exterior scaling is ideally suited to this problem because the matching point R, can be chosen larger than 9.5 a0 and then for a wide range of Q (e.g. 9 < 90’) there is nc difficulty because eq. (17) gives an analytic function in the exterior region. In table 1 we compare the results of the exterior-scaling procedure for the resonance energies of a few rotational predissociation resonances with the results obtained from analysis of phase shifts for this potential by Waech end Bernstein [ 1l]. The values of the parameters used in the calculations in table 1 were: R. = 15 ao, h = 1 0e3a0 and Q = 37”. The results are insensitive to the choices of Q from 12’ to 87’, and changes in R, to essentially any value larger than 15 a0 do not
Table 1 Comparison of results of exterior-scaling calculations with the phase-shift aaaiysis of Waech and Bernstein [ 121. u and j denote vibrational and rotational quantum numbers and allenergies are given in atomicunitstimeslo4 u
i
0 1 2 4 6 8 9 10 10 11 12 12 14
37 36 33 29 25 21 19
17 16 14 12 11 5
ReiE) 296.90662 293.75604
213.67048 148.31509 99.11053 64.10670 51.096GY
40.44963 36.70338 21.87402 17.54129 9.82197 2.08740
-0.13745 -3.78848 -0.47782 -0.57767 -0.61843 -0.92279 -1.36003 -2.19593 -0.06672 -0.40946 -1.63406 -0.05982 -0.40269
Re(E)
ImU
296.8 293.5 f 0.4 213.6 148.4 99.0 63.6 f 0.4 51.1 40.4 26.7 21.8 17.6 9.80 2.08
-0.14 f 0.01 -3.42 i 0.2 -0.48 r 0.05 -0.43 f.0.2 -0.66 c 0.05 -1.05 f 0.05 -1.3 to.02 -8.56 2 0.05 -0.07 i 0.05 -0.4310.02 -1.73 * 0.2 -0.07 2 O.Oi -0.45 IO.01
3....
.
.
’
I
I
1
.’
I
0
IO
20 r (00)
i’
I
I
I
30
.
I
I!
40
Fig. 2. Plot of RejP(R(r))] for the u = 10, j = 17 resonancr in II2 on the Waech-Bernstein potential. The scaling angle 0 is 37” and Ro is 22.75 an.
It is necessary to abandon Cooley’s algorithm which provides a conver?ient approach in the usual complexscaling case, but we have found th2t this is not a se& ous drawback. Beginning fmm a reasonable guess it generally requires less than ten iterations to fmd the resonance energies to the accuracy reported in table 1. A great advantage of this approach is that it requires analytic continuation of the potential to complex values of the coordinates only in the asymptotic or near-asymptotic region. This feature of the method allows it to be applied to essentially any form for th,e potential V(r). In fact, if the potential is of finite range, or can be approximated as such, no analytic continuation of V(r) is necessary at all. We have found that greatly increased stability with respect to the scaling angle Q results in exterior-scaling calculations for this reason. Atabek and Lefebvre [ 131 have shown how to extend the numerical complex-scaling approach to multichannel systems. This generalization can also be made for the exterior-scaling approach along similar lines. The advantages of tie exterior-scaling approach should obtain in that case as well.
Acknowledgement change the results. The maximum value of jR(r)l for the outer contour for all but the lowest-energy resonance in table 1 can be taken at any valce 20 a0 and must be 335 Q,, for the last entry in table 1. A plot of the real part ofP(r) for one of the Hz resonances is given in fig. 2 which reveals precisely how the exterior-scaling algorithm works. In this case the value of& was chosen to be 22.75 no_ The rapid oscillation at small r occurs in the effective potential well of V(r) and is similar to the expected behavior for a bound vibrational state. At larger values ofr the exponentially increasing, oscillatory behavior is essentially that of the outgoing Hankel function I;i’(,&) for complex k. At r =RO we “turn ‘Lhecorner” of the R(r) contour and exponentially decreasing oscillatory Sehavior follows.
5. Discussion We have shown how exterior scaling can be used in a numerica! calculation of complex resonance energies.
This work was supported by U.S. National Science Foundation Grant No. CHE-7907787.
References
[l] J. AguiIar and J. Combes, Common. Math. Phys. 22 (i971) 269; E. Balslev and J. Combes, Commun hfath. Phys. 22 (19%) 280. B. Simon, dommun. hfath Phys 27 (1972) 1. [2] G. Doolen, I. Nuttall and R. Stag& Phys. Rev. A10 (1974) 1612; G.D. Doolen, J. Phys. B8 (1975) 525; A.P. Hickman, A.D. Isascson and W.H. hlilier, Chem. Phys. Letters 37 (1976) 63: TX. RescIgno. C.W. McCurdy and A.E. Orel, Phys. Rev. -417 (1978) 1931; B.R. Junker and C.L. Huang, Phys. Rev. A18 (i978) 313; Y.K. Ho, Phys. Rev. A23 (1981) 2137; C.W. McCurdy, J.G. Lauderdale 2nd R.C. Mowrey, I. Chem. Phys 75 (1981) 1835. [3] Intern. J. Quan:cm Chem. 14 (1978) No. 4. [4] 0. Atabek, R. Lefebvre and A. Requena, hid. Phys. 40 (1980) 1107.
J. ?km@, C W McCurdy /Exterior complex scalingin moiecukvsystems [51 B. Shnon, Phys Letters 7IA f1979) 111. Iti] J-W. Cocfey, Math. Cornpitt. 15 (1961) 363. I71 3. Numerov, Publ. Observatoire Central Astrophys. Russ. 2 (1933) 188_ 181 I). Harttee, CaIcuIation of atomic sttucture (H&y, xeii York, 19573 p. as. [9! Af. Abmmowitz and I. Stegun, Haadbook of matbematicef functions (US Govt. Printing Offtce, Washington, 1964) cit. 25.
133
[iOl J.R. TayIor, Scattering theory CWiley,New York, 1972) p. 184 Iill T-G. Wiechand R.B. Bernstein, J. Chem. Phys. 46 (1967) $905. IlIZ] W. Koios and L. Wolniewicz, J. Chem. Phys. 43 (1965) 2429.
[i3] Cl. Atabek and R. Lefe@te, Chem. Phys. 56 (1981) 195.
THE APPLICATION OF EXTERIOR CMPLEX SCAL%G IN CALCULATIONS ON RESONANCES IN NUCLEA’? MOTION IN MOLECULAR SYSTEMS Julia TURNER and C. William McCURDY * DepuRmenr of Che?r’Gssny,Ohio Srore University. Columbus, Ohio 43210. USA Received 18 January 1962
The application of the method of exterior coniplex scaling to find complex resonance energies by direct numerical integration is presented. Details ofthe numerical procedures required are &cussed and results are presented for rotational predissociation resonances in Hz. The advantage of the exterior-scaling approach is that it requires analytic continuation of the potential to complex values of the coordinates only in the asymptotic or near-asymptotic region. Therefore any repiesentation of the potential is adequate, including piecewise-contiuuous representations such as cubic splices. in contrast, application of the usual complex-scaling (rotated coordites) procedure generally requires that the potential be represented by a single analytic function.
1. Introduction The method of complex coordinates has found most of its applications in finite-basis-set calculations of resonance parameters for electronic resonances. (See for some early work on mathematical aspects ref. [I ] and for some applications in atomic resonance calculations ref. [2] _) The usual complexcoordinate transformation is the one appropriate for dilatation analyiic hamiltoinians, e.g. atomic systems, and consists of scaling the radial coordinates r_ of each electron by a complex scaling parameter, rj + ei@$.. For a range of the angle Q in the scaling parameter e@ this transformation converts the wavefunction whose
ener,gy eigenvalue is the complex resonance ener,y, and which increases exponentially at large real values Of 3, to a square integrable Function. Thus the complex-coordinate transformation allows the resonance wavefunction of the transformed problem to be expanded in a basis of square integrable functions. Consequently the complex resonance energy corresponding to a pole of the S-matrix can be computed (approximately) as the eigenvalue of a finite matrix. It
would seem that basis-set expansion is the ideal way to take advantage of the properties of the complex* Alfred P. Sloan Research Fellow. O30~-olo4/82/OOOO-OOao/S 02.75 0 1982 North-Holland
coordinate transformation [3]. Recently however, Atabek et al. [4] have demonstrated that the complex-coordinate transformation can also be applied in numerical integration of the Schrijdinger equation to fmd resonance solutions. For resonances in nuclear motion, for example molecular predissociation resonances, it is likely that direct numerical integration will be more efficient than basis-set methods for resonance calculations because of the short wavelengths involved. However, there is a considerable drawback in the use of complex coordinates for nuclear motion because the potential V(r ei@) is required at complex values of the nuclear coordinates. Unlike the case of electronic resonances in atoms where the potentials are all Coulomb potentials and hence .known analytically for all values of the electronic coordinates, the potential for nuclear motion is generally known oniy from electronicstructure calculations at discrete values of the nuclear coordinates. The usual procedure is to fit these discrete values with some functional form, but in most cases the functional form is only piecewise continuous. For example, for a diatomic molecule cubic splines or power series might be used for small values of the internuclear separation and reciprocal powers at large values. Although such a functional form might be necessary to represent the potential accurately for
J. Turner. C IV. McCwdy
i28
/ ExtenM complex scaiing in molecularsystems
real values of the nuclear coordinates, it does not defme 2 useful andytic continuation to complex values, 2nd such an andytic continuation is exactly wP;at is required in the complexcoordinate approach. In this paper we propose a solution to this difficulty by demonstrating that Simon’s [S] “exteriorscaling” theory can be applied to numerical calcul2tions on resonances in nuckar motion. Exterior sc2ling requires the use of complex coordinates only in ihe asymptotic or near-asymptotic region, 2nd thus does not require azdyticity of the potential except in that region. The piecewise-analytic functions necessary for accurate representations cf the potential are therefore easily accommodnted in this approach. The outline of this paper is 2s follows. Section 2 briefly describes Simon’s exterior scaling theory and section 3 presents tie details of its application in numerical integration of the Schriidinger equation. In section 4 we present the results of the application of this approach to fiiding rotational predissociation resonances in H7 and in section 5 we comment on possible extensions of the method to other types of problems.
7,. Exterior
scaling
Simon proposed his theory of exterior scaling to treat the problem of ?he definition of molecuIar electronic resonance energics in the Born-Oppenheimer approximation [S] _In that problem the potential which electrons experience due to the fLved nuclei is not dilatation 21Aytic aboct any center, 2nd therefore ordinary complex scaling fails to provide 2 method for finding resonance states. The basic idea of exterior scaling is that, if the non-analytic region of the potential can be enclosed in a sphere of radius R,, one c2n make coordinates complex only outside of that sphere and still obtain the Same transformation of the spectrum of the hamiltoni2n which ordinary complex scaling achieves in the dilatation analytic case. For simplicity let us specialize our discussion to the case of 2 single particle. If r denotes the coordinates of the particle, Sinon chooses the contourR(r), parameterized by the magnitude T of r
R(r) = i- ,
O,cra?,
=RofeiO(r--Ro),
~
R,,
(1)
shown in fig. 1, along which to continue analytic2Ily the hamiltonian. He also introduces the operator uR,,@) which, when it acts on a function of r, carries out the exterior-scaling transformation according to uRo(d)*(r)
=J(f@@(rjf)
,
(2)
where J(r) is the square root of the jacobian
_
J(r) = [dR(r)[dr]
‘12R(r)/r
The hamiltonian
is thus transformed
= uR,(@H@i;(@)
&,f@‘)
so that the hamiltonian fi = (-@/2J7z)Vz becomes, H&(Q)
(3) to (4)
> of our one-particle
system
+ V(r)
(9
under the exterior-scaling = (-fi2/k)V’
according
tranfiformation
+ v(r),
= (-fiz/2m)e-ziQVa
+ V(R(r)?),
OGrGRg, R, < r. (6)
Under *Jle assumptions that the non-analytic behavior in the potential occurs for r < R, Simon has shown that the spectrum of the tmnsformed hamilton& fiR,,(Q) Consists Of discrete bound st2tes 2t red energies, and continuum branch cuts (one in our one-particle case) rotated into the complex plane by an 2ngIe equal to --2Re(Q). In other words, the spectrum of the electronic hamiltonian of z molecule in the Born-Oppenheimer approximation is transformed by the exterior-scr&ng transformation in precisely the same w2y 2s the spectrum of an atomic h2miltorkn is transformed ucder the complex-scaling transformation.
r plane
Fig. 1. The comples r plane showbg
the radial coordinate under the complex-scali% transformation (dotted) and under the exterior-scaling transformation (solid).
J. l’hmer. C. IV. McCimiy /Exterior
Briefly stated, the reason the exterior-scaling transformation works in this way is that it is the same as ordinary complex scaling asymptotically, and it is the asymptotic behavior of the wavefunction which is ivportant in *&e determination of the spectrum. Of course, any hamiltonian which is dilatation analytic is also analytic under the exterior-scaling transformation. Interaction potentials for nuclear motion generally fall off at larger as fast or faster than r-1 and have non-analytic behavior at smaller values of the coordinates which are induced only by the functional fomrs with which we choose to represent them. In the following section we show how to use the transformed hamiltonian of eq. (6) in a numerical-integration appreach to finding its resonance eigenvalues and wavefunctions.
3. Numerical method
and a computational
example
Again we restrict the discussion to a one-particle system; extension to the case of several channels will be discussed briefly in section 5. If we consider a case for which the potential is local and central (these are not essential restrictions) we obtain, on passing to partial waves, an equation for r times the radial function, denoted by P(r), d%‘(r)/dr*
= [U(r) - e]P(r)
3
(7)
where U(r) = (2m/@) [ V(r) + j(j i- 1)/2vzr2] ,
j is the angular-momentum
W
quantum number, m is the mass and E = 2mE/fi2. To fmd bound states of eq. (7) by numerical inte,oration it is often convenient to employ Cooley’s algori’&m [6], an iterative procedure which involves outward integration to some point and inward integration to the same point. Cooley’s contribution was to fmd an efficient algorithm for computing the step D(E) to take in the search for the bound-state energy at which the tivo integrations produce the same logarithmic derivative at the matching point. The entire procedure is based on the Nurnerov [7] integrator. If Pi denotes the value ofP(r) at the ith integration step, the Numerov integration is carried out with the threepoint recursion formula
complex
scaling in moleadarsystems
Yi,i + Yj_i - 2Yj
124
=h’(Uj - E)Pj
I
(9)
where 5
= [1-
(@112)(~j
- f)]Pi )
(10)
h is the step size and Vi is V(r) at the ith step. The error is of order h6 in each step. The results of the inward and outward integrations for Y are scaled to make their values equal at the matching point. Cooley’s expression for D(E) which results,
01) where Y, is the value of Y at the matching point, leads to second-order convergence in AE, the distance from the bound-state energy. Note that the error in the numerator of eq. (1 i) is also of order h6. Atabek et al. [4] demonstrated that eqs. (9)-(11) can be used to find resonance states in a complexscaling calculation. If ordinary complex scaling, r -+reib, is employed, Cooley’s algorithm is essentially unchanged except for char@ng the step h to heio. The values ofr at each step are then just multiples of hei*. The error at each step of the integration is still of order h6. Unfortunately, the Cooley procedure apparently cannor be similarly generalized to the exterior-scaling contour without changing the order in h of the error accumulated at a critical point in the integration. The three-term recursion formula in eq. (9) has an error of order h6 only if the integration
steps are equal, but they must be equal in directio~z as well as in magnitude if the cancellations which cause terms of order 123through h5 to vanish are to occur. Therefore where the exterior-scaling contour R(r) changes direction at R. the recursion formula in eq. (9) accumulates an error of order h3_ However, this point is a crucial one in the integration
because generally it will
be where exponential
decay of P(r) begins. One of the tests of the numerical calculation of the resonance ener,v is the extent to which the result is independent of the scaling angle 9. The leading error in the integration step withRu at its center is (@/3!)eiO(ezie
- l)(d3P/dr3)R,
,
and can be responsible for sensitivity of the resonance energy to 4 when all other components of the calcula-
J. ZU?ner, C IV. MXkrdy
130
/Exterior
complex scalitig in mdecular
tion show negligible dependence on & This sensitivity to Q is zpparent if we choose R0 as the matching point in Coolmy’s algorithm, in which case eq. (11) becomes’
’
D(E) = [(eio + I)/21 {(urn - ~)p, + [-SOY,
gxtems
b-
\p(Ro)
_1 + (ei@ + i)Y,, - Yrr+l ]/S@eib(&+l)] _
+ 1 @(R(r))[dR(r)ldr]
dr])-l
.
04)
Ro
The numerator of this expression is just the generalization of eq. (9) (divided by h2) whenRo is the mth point. The error in D(E) is thus of order h3 and if the criterion fo: convergence to the resonance energy is that Ll(E) vanish, the error in the computed resonance ener,v is of order /23 even though the error in every other component of the calculation is of order h6. However, usingRO as the matching point suggests a way of avoiding this problem altoge’rher. Instead of Cooley’s algorithm for D(E) we can use the corrector due to Douglas and Ridley [S]
me derivatives in the numerator are evaluated at R. and P(r) and P(R(r)) are obtained by Numerov integration in the ouhvard and inward directions respectively and scaled to match at R, with the (arbitrary) value P(R,-J. In the inward integration the step h in eqs. (9) and (10) becomes he@‘_ To computeB{E) in eq. (14) so that it will vanish up to order h6 we merely need to compute the derivatives in the numerator to that order with six-point differentiation formulas using the results of the six integration steps adjacent toRo_ For example, for the outward derivative we use [9] [rlPold%,t
=(1/120h)(-24Pm_S
- 400 P,I_3
’
[pj&[P2(d d,!-’ ,
+ SOP,_4
+ 6C0 P,B_2 - 600 P,,_1 f 274 Pm) , 05)
03)
where ‘;n is the matching point, and Ph,(,,) and P;,(rrn) are the derivatives ofP(r) at the matching point from the outward and inward integrations respectively. This form forD(E) does not lead to quadratic convergence in the energ iteration, but we can cause it to vanish up to order h6 at convergence. The Douglas-Ridley corrector in eq. (13) can be generaiized for the‘ exterior sc&ng contour R(r), and ifR, is chosen as the matching point we obtain
* Derivation of eq. (12) requires that Pi and Yi be redefined to contain jacobian factors from eq. (3). Pi becomes e*/zpi for i > m md [@ t 1)/211” 1Dn2 at the matching point. With these changes and some simple manipulations, the kinetic-rnergy contribution to Cooley’s matrix F.4in his equsrion (3.3) czm be made symmetric.
with F, denoting P(R& and a similar formula obtains for the other derivative. We have found L3at it is adequate to use extended Simpson’s rule [9] to do the two integrals in the denominator of eq. (14). Then the numerator vanishes with error sixth order in 12at convergence while the denominator contains errors of order h5. We have found that eq. (12) and sirrar modikation of Cooley’s algorithm are inferior to our modification of the Douglas-Ridley co.qector and can lead to instability with respect to the scaling angle 9. Our more stable alternative, the numerical exterior-scaling algorithm, can be summarized as follows: (1) Choose a first g!JessE to the resonance ener,qy. (2) Integrate outward from P = 0 to r =R, with Numerov integration [eqs. (9) and (lo)]. (3) Starting with a value ofR(r) in the asymptotic region and the initial condition P(R(r)) = t$~&R(r)), where i,k is the Ricatti Hankel function (see for example ref. [lo]), integrate inward using eqs.(9) and (10) with h changed to heiS. The scaling angle Q must be large enough to expose the resonance
by rotating the
Jl Turner, C PJ.McCurdy /Exterior complex scaring in molecule? systems continuum
and the potential is extended dependence for r > 9 .S a0 ,
branch cut pest it;The choice ofRO is not critical but the least numerical sensitivity to its value occurs withR, chosen larger than the location of any classical turning points in r. (4) Scale the results of the integrations to match with the value P(Ra) at R, . Compute the numerical derivatives in eq. (14) with six-point formulas [eq. (15)] ami compute the integrals with extended Simpson’s rule. (5) Change the value of E to E +D(E) [eq. (14)] and return to step (2). In the following section we will show how this method can be used with an accurate representation of a potential which does not provide a useful analytic continuation to complex7 for small lrl. Note that no complex arithmetic isnecessary to compute U(r) for the outward integration.
V(r) = --c,/fi
4.A numericalexample Waech and Bernstein [ 111 have given a representation of the highly accurate KoTos and Wolniewicz potential [12] for Hz. For values ofr between 0.4 and 9.5 a0 they use a 33-term power series of the form (atomic units)
06)
- c,/r* - c&*
131 with reciprocal
.
power
(17)
The coefficients an, C', ,C8 and Cl0 are given in reference 1121 and references therein. There is a small discontinuity at r = 9.5 of 0.21 cm-l. Attempts to use ordinary complex scaling for this problem run into grave trouble immediately, because, although the interior and exterior parts of the potential match well at real I, for complex r there is no simple boundary in the complex plane along which (16) and (17) give the same value. This probfem is more acute at large values of Q because the contour rei@ can probe the edges of the region of convergence of eq. (IS) where oscillations occur. Exterior scaling is ideally suited to this problem because the matching point R, can be chosen larger than 9.5 a0 and then for a wide range of Q (e.g. 9 < 90’) there is nc difficulty because eq. (17) gives an analytic function in the exterior region. In table 1 we compare the results of the exterior-scaling procedure for the resonance energies of a few rotational predissociation resonances with the results obtained from analysis of phase shifts for this potential by Waech end Bernstein [ 1l]. The values of the parameters used in the calculations in table 1 were: R. = 15 ao, h = 1 0e3a0 and Q = 37”. The results are insensitive to the choices of Q from 12’ to 87’, and changes in R, to essentially any value larger than 15 a0 do not
Table 1 Comparison of results of exterior-scaling calculations with the phase-shift aaaiysis of Waech and Bernstein [ 121. u and j denote vibrational and rotational quantum numbers and allenergies are given in atomicunitstimeslo4 u
i
0 1 2 4 6 8 9 10 10 11 12 12 14
37 36 33 29 25 21 19
17 16 14 12 11 5
ReiE) 296.90662 293.75604
213.67048 148.31509 99.11053 64.10670 51.096GY
40.44963 36.70338 21.87402 17.54129 9.82197 2.08740
-0.13745 -3.78848 -0.47782 -0.57767 -0.61843 -0.92279 -1.36003 -2.19593 -0.06672 -0.40946 -1.63406 -0.05982 -0.40269
Re(E)
ImU
296.8 293.5 f 0.4 213.6 148.4 99.0 63.6 f 0.4 51.1 40.4 26.7 21.8 17.6 9.80 2.08
-0.14 f 0.01 -3.42 i 0.2 -0.48 r 0.05 -0.43 f.0.2 -0.66 c 0.05 -1.05 f 0.05 -1.3 to.02 -8.56 2 0.05 -0.07 i 0.05 -0.4310.02 -1.73 * 0.2 -0.07 2 O.Oi -0.45 IO.01
3....
.
.
’
I
I
1
.’
I
0
IO
20 r (00)
i’
I
I
I
30
.
I
I!
40
Fig. 2. Plot of RejP(R(r))] for the u = 10, j = 17 resonancr in II2 on the Waech-Bernstein potential. The scaling angle 0 is 37” and Ro is 22.75 an.
It is necessary to abandon Cooley’s algorithm which provides a conver?ient approach in the usual complexscaling case, but we have found th2t this is not a se& ous drawback. Beginning fmm a reasonable guess it generally requires less than ten iterations to fmd the resonance energies to the accuracy reported in table 1. A great advantage of this approach is that it requires analytic continuation of the potential to complex values of the coordinates only in the asymptotic or near-asymptotic region. This feature of the method allows it to be applied to essentially any form for th,e potential V(r). In fact, if the potential is of finite range, or can be approximated as such, no analytic continuation of V(r) is necessary at all. We have found that greatly increased stability with respect to the scaling angle Q results in exterior-scaling calculations for this reason. Atabek and Lefebvre [ 131 have shown how to extend the numerical complex-scaling approach to multichannel systems. This generalization can also be made for the exterior-scaling approach along similar lines. The advantages of tie exterior-scaling approach should obtain in that case as well.
Acknowledgement change the results. The maximum value of jR(r)l for the outer contour for all but the lowest-energy resonance in table 1 can be taken at any valce 20 a0 and must be 335 Q,, for the last entry in table 1. A plot of the real part ofP(r) for one of the Hz resonances is given in fig. 2 which reveals precisely how the exterior-scaling algorithm works. In this case the value of& was chosen to be 22.75 no_ The rapid oscillation at small r occurs in the effective potential well of V(r) and is similar to the expected behavior for a bound vibrational state. At larger values ofr the exponentially increasing, oscillatory behavior is essentially that of the outgoing Hankel function I;i’(,&) for complex k. At r =RO we “turn ‘Lhecorner” of the R(r) contour and exponentially decreasing oscillatory Sehavior follows.
5. Discussion We have shown how exterior scaling can be used in a numerica! calculation of complex resonance energies.
This work was supported by U.S. National Science Foundation Grant No. CHE-7907787.
References
[l] J. AguiIar and J. Combes, Common. Math. Phys. 22 (i971) 269; E. Balslev and J. Combes, Commun hfath. Phys. 22 (19%) 280. B. Simon, dommun. hfath Phys 27 (1972) 1. [2] G. Doolen, I. Nuttall and R. Stag& Phys. Rev. A10 (1974) 1612; G.D. Doolen, J. Phys. B8 (1975) 525; A.P. Hickman, A.D. Isascson and W.H. hlilier, Chem. Phys. Letters 37 (1976) 63: TX. RescIgno. C.W. McCurdy and A.E. Orel, Phys. Rev. -417 (1978) 1931; B.R. Junker and C.L. Huang, Phys. Rev. A18 (i978) 313; Y.K. Ho, Phys. Rev. A23 (1981) 2137; C.W. McCurdy, J.G. Lauderdale 2nd R.C. Mowrey, I. Chem. Phys 75 (1981) 1835. [3] Intern. J. Quan:cm Chem. 14 (1978) No. 4. [4] 0. Atabek, R. Lefebvre and A. Requena, hid. Phys. 40 (1980) 1107.
J. ?km@, C W McCurdy /Exterior complex scalingin moiecukvsystems [51 B. Shnon, Phys Letters 7IA f1979) 111. Iti] J-W. Cocfey, Math. Cornpitt. 15 (1961) 363. I71 3. Numerov, Publ. Observatoire Central Astrophys. Russ. 2 (1933) 188_ 181 I). Harttee, CaIcuIation of atomic sttucture (H&y, xeii York, 19573 p. as. [9! Af. Abmmowitz and I. Stegun, Haadbook of matbematicef functions (US Govt. Printing Offtce, Washington, 1964) cit. 25.
133
[iOl J.R. TayIor, Scattering theory CWiley,New York, 1972) p. 184 Iill T-G. Wiechand R.B. Bernstein, J. Chem. Phys. 46 (1967) $905. IlIZ] W. Koios and L. Wolniewicz, J. Chem. Phys. 43 (1965) 2429.
[i3] Cl. Atabek and R. Lefe@te, Chem. Phys. 56 (1981) 195.