2.2 A program for calculating cross sections 19
calculate differential cross sections, we need Legendre polynomials too. In
Appendix A2, iterative methods for evaluating special functions are discussed.
- Finally, we complete the program with a routine for calculating the cross
sections from the phase shifts.
2.2.1 Numerov’s algorithm for the radial Schrödinger equation
The radial Schrödinger equation is given inEq. (2.3). We define
F(l,r,E)=V(r)+
^2 l(l+ 1 )
2 mr^2
−E (2.10)
so that the radial Schrödinger equation now reads:
^2
2 m
d^2
dr^2
u(r)=F(l,r,E)u(r). (2.11)
Units are chosen in which^2 /( 2 m)assumes a reasonable value, that is, not
extremely large and not extremely small (see below). You can choose a lib-
rary routine for integrating this equation but if you prefer to write one yourself,
Numerov’s method is a good choice because it combines the simplicity of a regular
mesh with good efficiency. The Runge–Kutta method can be used if you want to
have the freedom of varying the integration step when the potential changes rapidly
(see Problem 2.1).
Numerov’s algorithm is described inAppendix A7.1. It makes use of the special
structure of this equation to solve it with an error of orderh^6 (his the discretisation
interval) using only a three-point method. For^2 / 2 m≡1 it reads:
w(r+h)= 2 w(r)−w(r−h)+h^2 F(l,r,E)u(r) (2.12)
and
u(r)=
[
1 −
h^2
12
F(l,r,E)
]− 1
w(r). (2.13)
It is useful to keep several things in mind when coding this algorithm.
- The functionF(l,r,E), consisting of the energy, potential and centrifugal
barrier, given in Eq. (2.10), is coded into a functionF(L,R,E), withLan
integer andRandEbeing real variables.
- As you can see fromEq. (2.9a), the value of the wave function is needed for two
values of the radial coordinater, both beyondrmax. We can taker 1 equal to the
first integration point beyondrmax(if the grid constanthfor the integration fits
an integer number of times intormax, it is natural to taker 1 =rmax). The value
ofr 2 is larger thanr 1 and it is advisable to take it roughly half a wavelength
beyond the latter. The wavelength is given byλ= 2 π/k= 2 π/
√
2 mE.As