Nuclear Physics A463 (1987) 139c - 144c North-ltolland, Amsterdam
139c
THE LANCZOS ALGORITHMAND SPLINE EXPANSIONS G.L. PAYNE Department of Physics and Astronomy, The U n i v e r s i t y of lowa, lowa C i t y , lowa 52242, U.S.A.*
1. INTRODUCTION There e x i s t many techniques to solve the Schrodinger equation for three nucleons with r e a l i s t i c two- and three-body i n t e r a c t i o n s .
The solution of the
Faddeev equations in configuration space has the advantage that the equations have a simple form, which can be programmed in a straightforward manner. difficulty
The
with t h i s method is that the standard techniques for solving the
p a r t i a l d i f f e r e n t i a l equations can lead to numerical c a l c u l a t i o n s which are p r o h i b i t i v e even on modern supercomputers.
In t h i s paper we review two tech-
niques which have been used to reduce the problem to one which can be solved e a s i l y on a CRAY or s i m i l a r computer.
T r a d i t i o n a l l y , Faddeev c a l c u l a t i o n s ,
both in configuration space and momentum space, are categorized by numbers of channels.
These channels specify the angular momentum quantum numbers of the
i n t e r a c t i n g pair of nucleons, and the corresponding quantum numbers of the remaining spectator nucleon with respect to the i n t e r a c t i n g pair.
To obtain a
converged solution to the bound-state problem, one must include enough channels to accurately represent the three Faddeev amplitudes. t i o n s I in configuration space used 5 or 18 channels.
The i n i t i a l
calcula-
Recently, i t has been
shown2 that i t is necessary to use 34 channels to obtain a converged r e s u l t for the case where there is a three-body force as well as the usual two-body force. 2. CONFIGURATIONSPACE FADDEEV EQUATIONS The Faddeev equations in c o n f i g u r a t i o n space were f i r s t f o r the three-body scattering problem.
derived by Noyes3
These equations were f i r s t
used for the
bound-state problem by the Grenoble group 4 who used f i n i t e difference techniques to solve the equations.
To derive the equations one introduces the
Jacobi variables for three equal mass p a r t i c l e s with position vectors ~ i , ~ j , and ~k : ~
;
Yi = ~ ( r j + ~k) - ~i
,
*This work was supported, in part, by the U.S. Department of Energy. 0 3 7 5 - 9 4 7 4 / 8 7 / $ 0 3 . 5 0 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)
(1)
G.L. Payne / The Lanczos algorithm and spline expansions
140c
where i , j ,
k imply c y c l i c permutation.
Using these v a r i a b l e s the Schrodinger
wave f u n c t i o n is w r i t t e n as the sum of the three Faddeev amplitudes: +
= (1 + P- + P+)@(~z,Yl) =
e(~1,Yl) + ~(xZ,Y2) + ~(X3 ,Y3)
= ~i
+ ~b2 + ~3
where P- and P+ are the c y c l i c permutation operators.
(2)
,
For a system with a
Hamiltonian of the form: H = T + V(~I) + V(~2) + V(;~3) = T + VI + V2 + V3
,
(3)
where T is the k i n e t i c energy operator and V(~i) is the two-body p o t e n t i a l acting between p a r t i c l e s j and k, the Schrodinger equation can be replaced by the three coupled Faddeev equations (T - E)~ i + V ( ~ i ) ( ~ i + ~j + ~k) = 0
(4)
For the case of i d e n t i c a l p a r t i c l e s a l l three Faddeev amplitudes have the same f u n c t i o n a l form, and i t
is only necessary to solve one of the equations.
A f t e r expanding the Faddeev amplitude in a complete set of channel f u n c t i o n s , one obtains a set of coupled e l l i p t i c
partial differential
equations. The
standard numerical methods f o r solving these equations f o r the eigenvalue, E, and the e i g e n f u n c t i o n lead to p r o h i b i t i v e matrix equations to solve.
The size
of the numerical c a l c u l a t i o n can be considerably reduced by using a scheme s i m i l a r to one used in momentum space c a l c u l a t i o n s . 5
This method consists of
r e w r i t i n g the Faddeev equation in the form: [E - (T + V1)]~l : XVl(@2 + @3)
(5)
Now one assumes a value f o r E and solves t h i s equation f o r the eigenvalue x. The value of E is varied u n t i l
the eigenvalue ~ is u n i t y .
Using the Lanczos
method described below, one must i n v e r t a matrix corresponding to the l e f t - h a n d side of equation (5).
Since the operator (T + V - E) couples at most two chan-
nels ( f o r the case with a tensor f o r c e ) , the r e s u l t i n g matrix equation never i n v o l v e s more than two channels at a time.
Therefore, the computer time
required to solve the equations varies l i n e a r l y with the number of channels, r a t h e r than as the cube. 3. SPLINE EXPANSION To solve equation (5), we f i r s t
expand the Faddeev amplitude in a complete
set of channel functions
,Yl) = Z
*~(xI'Yl) xlYl
(6)
G.L. Payne / The Lanczos algorithm and spline expansions
141 c
where ~ l a b e l s the various angular momentum and s p i n - i s o s p i n f u n c t i o n s used to d e s c r i b e the three-body system.
A l s o , we have introduced the reduced channel
f u n c t i o n s @ a ( x l , y l ) which have the boundary c o n d i t i o n s t h a t the f u n c t i o n s are zero when xz or Yl is zero.
M u l t i p l y i n g t h i s equation by x l Y l and t a k i n g the
i n n e r product w i t h <~I gives a set of coupled e l l i p t i c
differential
equations.
In a d d i t i o n , we i n t r o d u c e the h y p e r s p h e r i c a l v a r i a b l e s x l = p cos o
,
This change of v a r i a b l e s s i m p l i f i e s
Yl = ~
P sin e
(7)
the p r o j e c t i o n i n t e g r a l s on the r i g h t - h a n d
side of the e q u a t i o n , as they depend only on the v a r i a b l e e and not p.
Thus,
we o b t a i n the set of coupled equations
(E - T)@a(p,O) - Z V B(p COS 8)@8(p,8) O+ V¥(p COS O) f _ Kyg(e,o')@6(p,O)de' By o
,
= ~ Z
(8)
where V~B(p cos e) is the projection of the two-body potential onto the channel states, and the kernel Ky~(e,o'), which depends upon the angular momentum and isospin variables, can be evaluated by standard techniques. 6 To solve this equation we expand the channel functions using a basis set of bicubic Hermite s p l i n e s on a r e c t a n g u l a r g r i d in the p-O coordinates. We choose an expansion o f the form: M
@a(P,O) = [ Z
m e-KP Z amn Sm(P)Sn(O)] p~/Z
N
m=l n=l
,
(g)
where ~ K 2 / M = -E, and the Sm and s n are the cubic Hermite s p l i n e s . t o have a smoother f u n c t i o n to be f i t
out the asymptotic behavior of the bound-state f u n c t i o n . t i o n s f o r the d i f f e r e n t i a l
The boundary c o n d i -
equation are i n c o r p o r a t e d i n t o the basis f u n c t i o n s .
The Hermite s p l i n e s are l o c a l f u n c t i o n s w i t h a continuous f i r s t These f u n c t i o n s are defined by d i v i d i n g the region to be f i t intervals,
In o r d e r
by the s p l i n e expansion, we have f a c t o r e d
derivative.
i n t o a number of
where the ends o f these i n t e r v a l s are c a l l e d the b r e a k p o i n t s .
At
each b r e a k p o i n t t h e r e are two s p l i n e s which are nonzero in the a d j a c e n t i n t e r vals.
In Figure I we i l l u s t r a t e
the s t r u c t u r e of these f u n c t i o n s f o r the case
w i t h f i v e i n t e r v a l s and the boundary c o n d i t i o n t h a t the f u n c t i o n be zero at each end of the r e g i o n .
Since the s p l i n e s are l o c a l f u n c t i o n s , t h e r e are at
most f o u r nonzero s p l i n e s on each i n t e r v a l
f o r the cubic s p l i n e s .
We use the
method of orthogonal c o l l o c a t i o n 7 t o determine the expansion c o e f f i c i e n t s
asmn.
G.L Payne / The Lanezos algorithm and spline expansions
142c
The collocation procedure requires that one determine the values of as mn for which the d i f f e r e n t i a l equation
c- c.B6 584
HERMITESPLINES l
-
-
is s a t i s f i e d at M d i s t i n c t values of Pi and N d i s t i n c t values of ej. I f
-
one chooses these points to be the ×
5
two point Gauss quadrature points on
6 7
-
each i n t e r v a l , the procedure is known as orthogonal collocation. The col-
-
8
location points for the example shown
9
in Figure i are marked by x's along
IO x
x
O
I
5
x
x~
x
x
IO
; x
15
x
I x
20
25
x
the x-axis. The orthogonal collocation tnethods leads to the set of simultaneous equations for the expansion coefficients a(% mn
FIGURE 1 Hermite spline functions for f i v e integral s. M
N
~ (~ +1) m
L (L +1)
m
m=IZ n=l ~ {Sm(Pi)Sn(°~)J + P-~I [¼_ c ° s ~ e J
(%
- ~
m
]
Sm(Pi)Sn(ej)
i (% + --2 Pi Sm(Pi)Sn(O~) J - 2~ s'm(Pi)Sn(Sj)}amn M
N
Z
Z
M
N
Z
Z
Nc
m=l n=l ~=I
:x
[v 6(p i cos ej)Sm(Pi)Sm(ej)]a~n 8+
Nc
m=l n=l ~=1
{~¥ v y(p i cos e j ) [ f 0 1 Ky6(ej,e')Sn(O')de']Sm(Pi)}aBmn "(i0)
In these equations ~m is the orbital angular momentum of the interacting pair 2-3, L(% is the orbital angular momentum of the spectator p a r t i c l e i , and vm6 is (M/~Y~)Vm6. This equation can be written symbolically as
M N M c Aijm,mn6 a~ ~ ~c Bijm,mn~ a6 mn Z N Z mn = ~ ~ m=l n=l 6=i m=l n=l B=I
(11)
which can be written as the generalized matrix eigenvalue problem Aa = ~Ba
,
(12)
where the eigenvector is the array of coefficients for the spline expansion. For a r e a l i s t i c calculation the size of this matrix can be quite large. A typical 34-channel (Nc = 34) calculation may require 30 splines for the p
G.L. Payne / The Lanczos algorithm and spline expansions variable (M = 30) and 30 splines for the 9 variable (N = 30).
143c
Thus, the total
number of expansion coefficients would be 34 x 30 x 30 or 30,600 unknowns. Consequently, this matrix is too large to be solved by straightforward techniques, and we use a variation of the Lanczos algorithm to find the solution. 4. LANCZOS METHOD Since we wish t o f i n d the e i g e n v a l u e s c l o s e s t t o u n i t y , e q u a t i o n (12) in the form A- i
B a = (I/~)a
we f i r s t
rewrite
or
H a = A a To s o l v e t h i s
(13)
e q u a t i o n , we use a m o d i f i e d Lanczos a l g o r i t h m t o g e n e r a t e a small
basis set which can be s o l v e d by standard t e c h n i q u e s .
Since the m a t r i c e s A and
B are not symmetric, we have t o generate a b i o r t h o g o n a l basis s e t . done using the method o f Saad. 8
Starting
w i t h an i n i t i a l
This i s
vector ai,
e r a t e the basis set a i as w e l l as the b i o r t h o g o n a l basis set bi f o r . . o ~.
These basis v e c t o r s are generated by f i r s t
we geni = 1,2,
using the r e c u r s i o n
relations N
a i + l = Hai - ~ i a i
(14)
- ~iai-1
and bi+l
= HTbi - ~ i b i
where a i = b~Hai , Bi = b T a i / Y i , basis v e c t o r s .
and ¥i =
" Yibi-1
'
(15)
I T iL I/2,
t o o b t a i n the unnormalized
The n o r m a l i z e d v e c t o r s are d e f i n e d as a i = a i / ~ i
bi = ~ i / B i . Thus, one can generate u v e c t o r s ai and t h e i r
and
b i o r t h o g o n a l v e c t o r s b i such
t h a t bT aj = a i j . Using t h i s
small basis s e t , we use t h e expansion V
a = Z ~iai
(16)
i=1
Substituting this expansion into equation (13) and taking the inner product with bj leads to the matrix equation V
v L
~ibTHai = Z Hji~i = A~j
i =1
The e i g e n v a l u e s and e i g e n v e c t o r s o f t h i s t o the l a r g e m a t r i x e q u a t i o n .
(17)
i=1 m a t r i x e q u a t i o n are a p p r o x i m a t i o n s
In p r a c t i c e one i n c r e a s e s u u n t i l
A i s l e s s than the d e s i r e d accuracy o f s i x s i g n i f i c a n t
figures.
the change in
G.L. Payne / The Lanezos algorithm and spline expansions
144c
The most d i f f i c u l t
numerical procedure of this method is solving the recur-
sion r e l a t i o n s in equations (14) and (15). sion of the matrix A.
This requires the numerical inver-
However, the choice of the equations to be solved leads
to a reduced matrix whose blocks never have more than two channels, and each of these blocks can be inverted separately.
F i n a l l y , the inversion of these
blocks is considerably s i m p l i f i e d because of the local nature of the splines. Since there are at most four nonzero splines on each i n t e r v a l , the r es u lt in g matrix is a sparse banded matrix.
By the proper ordering of the equations and
the unknowns, one obtains a block t r i d i a g o n a l matrix whose blocks are banded matrices.
This matrix can be inverted by using an algorithm which works on a
few blocks at a time.
Thus, the amount of computer memory needed does not
become excessive, and the inversion of A can be performed in an e f f i c i e n t manner. 6. CONCLUSIONS The use of spline expansions provides an e f f i c i e n t method for obtaining accurate solutions of the configuration-space Faddeev equations. i s t i c c a l c u l a t i o n one must use a modified form of these equations.
For a r e a l The Lanczos
algorithm f o r nonsymmetric matrics can be used to find the solution f o r the bound-state problem. REFERENCES i) G.L. Payne, J.L. F r i a r , B.F. Gibson, I.R. Afnan, Phys. Rev. C 22 (1980) 823. 2) C.R. Chen, G°L. Payne, JoL. F r i a r , B.F. Gibson, Phys. Rev. C 33 (1986) 1740. 3) H.P. Noyes, in Three Body Problem in Nuclear and P a r t i c l e Physics, ads. J.S.C. McKee and P.M. Rolph (North-Holland, Amsterdam, 1970), p. 2. 4) S.P. Merkuriev, C. Gignoux, and A. Laverne, Ann. Phys. (NY) 99 (1976) 30. 5) Ch. Hajduk and P.U. Sauer, Nucl. Phys. A369 (1981) 321. 6) C.R. Chen, Effects of three-nucleon forces on the trinucleon system, Ph.D. d i s s e r t a t i o n , University of lowa, 1985. 7) P.M. Prenter, Splines and V a r i a t i o n a l Methods (Wiley, New York, 1975). 8) Y. Saad, SIAM J. Num. Anal. 19 (1982) 485.