3 Wave Functions for Scattering and Bound States

In order to describe details of the nuclear transitions realistically, it is necessary to specify in sufficient detail the initial and final states of the nuclei involved. To start with, the excitation energies, spins and parities of all the states in each mass partition need to be specified, along with the nuclear masses, charges and relative Q-values of the partitions. Each pair projectile and target excited states is then a distinct channel with its own scattering wave function and boundary conditions. The initial projectile and target states will select one such channel as the `incoming channel', with its boundary conditions specifying an incoming plane wave. All channels (including the incoming channel) will have outgoing spherical waves. Particular attention must be given to the consistent placement of i^L factors in these definitions.

The individual nuclear states are then specified in sufficient detail for the particular reaction mechanisms involved. It is not necessary to specify the full quantum mechanical states of all the nucleons in the nucleus, but rather, only the states of those changed in the reactions being considered. In particular, one and two-nucleon wave functions will have to be described, if those nucleons are to be transferred to other nuclei. If a nuclear state consists of a particle of spin s bound outside a nucleus with possible core states $\phi_I$, then the bound state radial wave functions $u_{\ell s j I} (r)$ will have to be found by solving a coupled-channels set of equations for negative energy eigen-solutions. If the particle is not bound, in the other hand, then its continuum range of energies must be discretised into a finite collection of `bin' states which can be scaled to unit normalisation. If the nuclear state consists of two particles of intrinsic spins s₁ and s₂ outside a core, then it is usually specified by a shell-model or by a Sturmian-basis calculation in terms of the independent coordinates ${\bf r}_1$ and ${\bf r}_2$. To calculate transfer rates, however, the two-particle wave functions need to be given in terms of the collective coordinates (usually ${\bf r}=\frac{1}{2}( {\bf r}_1 + {\bf r}_2 )$ and ${{\mbox{\boldmath$\rho$}}} = {\bf r}_1 - {\bf r}_2$). In order to use the states in a reaction calculation, therefore, equations are given for the transformation from the independent coordinates.

When we have calculated the scattering wave functions, or at least their asymptotic parts in terms of their S-matrix elements, we can find the cross sections for each outgoing pair of projectile and target states in each partition. Furthermore, if the initial projectile has non-zero spin J_p, then the effect on these cross sections of polarisation of the projectile is specified by the tensor analysing powers T_kq (for $1 \leq k \leq 2J_p$ and $0 \leq q \leq k$). Integrated cross sections and fusion polarisations can also be found using the S-matrix elements.

3.1 Total wave function

In each partition $\kappa$ of the system into a projectile of mass $A_{\kappa p}$ and a target of mass $A_{\kappa t}$, the coupling order is

$${\bf L}+ {\bf J}_p = {\bf J};~~~ {\bf J}+ {\bf J}_t = {\bf J}_T ,$$ (24)

which may be defined by writing

$$\Psi_{\kappa J_T}^{M_T} = \left \vert ~(L~J_p )J, J_t ; J_T \right\rangle$$ (25)

where J_p = projectile spin, J_t = target spin, L = orbital partial wave, and J_T = total system angular momentum.

The set $\{\kappa pt,~(L~J_p )J, J_t ; J_T \}$ will be abbreviated by the single variable $\alpha$. Thus, in each partition,

$$\Psi_{\kappa J_T}^{M_T} \left ( {\bf R}_\kappa , \xi_p , \xi_t \right )= \sum_{\begin{array}{cc}L J_p J J_t\\ M \mu_p M_J \mu_t\end{array}} \phi_{J_p} ( \xi_p )~\phi_{J_t} ( \xi_t )~ i^L Y_L^M ({\hat {\bf R}_\kappa} ) {1 \over {R_\kappa}} ~ f^{~\kappa J_T}_{(LJ_p )J,J_t} (R_\kappa)$$

$$\left\langle LM J_p \mu_p \vert JM_J \right\rangle \left\langle JM_J J_t \mu_t \vert J_T M_T \right\rangle$$ (26)

here ${\bf R}_\kappa$ = radial coordinate from the target to the projectile in partition $\kappa$, $\xi_p$ = internal coordinates of projectile, $\xi_t$ = internal coordinates of target, and

$$f^{~\kappa J_T}_{(LJ_p )J,J_t} (R) \equiv f _\alpha (R)$$ (27)

are the radial wave functions. The i^L factors arise from the spherical Bessel expansion of the incoming plane wave. Some formalisms include extra powers of i in the equation (26), in order to make the coupling interactions $V_{\alpha :\alpha'}$ real. Inelastic coupling interactions can be made real (for integer-valued spins, at least), by including a factor i^{J_p + J_t} in the definition of $\Psi$, and transfer couplings can be made real by including a factor $i^ \ell$ for orbital angular momentum $\ell$ of the bound particle state in this partition. In a general purpose code [34], however, there may be clashes between these different conventions. The ground state of ⁷Li, for example, would have a factor of i3/2 on the rotational model convention, but a factor of i¹ if the state were regarded as a $\ell = 1$ bound state of a $\alpha$ core and a triton cluster. It seems simplest, therefore, to omit these addition phase factors completely. The coupling interactions can very often be made real, nevertheless, if the i L factors are included explicitly in the CRC equations, as in the next section.

The wave function $\Psi$ could also have been defined using the `channel spin' representation (as in [43]) $\Psi = \vert L, (J_p J_t)S; J_T\rangle$, which is symmetric upon projectile $\rightarrow$ target interchange except for a phase factor (-1)^{S - J_p - J_t}. This would simplify the subsequent description of the coupling elements in section 4, as the formulae for projectile mechanisms and target mechanisms would differ only by this phase factor. However, the channel spin representation has the disadvantage that the projectile spin-orbit force is not diagonal in this basis. This would not matter if coupled-channels solutions were always sought, but one of the advantages of sometimes solving the CRC equations iteratively is that the DWBA solutions of first and second order (etc.) may be obtained. In order for the partially-iterated CRC solutions to reproduce the results of DWBA codes, it is necessary to treat spin-orbit forces without approximation, and since spin-orbit forces almost always are those of the projectile, the asymmetric representation of channel (24) is advisable.

3.1.0.1 Identical Nuclei

If one partition ($\kappa$ say) is identical to another $\kappa$ except that the projectile and target nuclei are exchanged, then the total wave function should be formed from $(1 + \pi P_{\kappa\kappa'})$ times the above expression, where $\pi =\pm1$ is the intrinsic parity of the two nuclei under exchange. A simple method of dealing with this exchange is to first form the wave function of equation (26), and then operate with $(1 + \pi P_{\kappa\kappa'})$ on both the wave functions and the S-matrix elements, before cross sections are calculated. This is equivalent to the replacement

$${f _\alpha (R) \leftarrow f _\alpha (R) + c_{\alpha ,\alpha'} f_{\alpha'} (R) }$$ (28)

where

$$c_{\alpha ,\alpha'} = \pi (-1)^L \delta_{L,L'} (-1)^{J_T + L - J - J'} \hat{J} \hat{J'} W(J_p L' J_T J_t ; J J')$$ (29)

with $\alpha' = \vert (L' J_t)J' , J_p ; J_T\rangle$.

3.2 Coupled equations

The CRC equations are in many cases of the form

$$\left [ E_{\kappa pt } - T_{\kappa L } {(R_\kappa)}- U_\kappa {(R_\kappa)}\right ] f_\alpha {(R_\kappa)} \sum_{\alpha' , \Gamma> 0 } { i ^ {L' - L} ~ V^\Gamma_{\alpha :\alpha' } {(R_{\kappa'})}f_{\alpha'} {(R_{\kappa'})}}$$ =

$$\sum_{\alpha' ,\kappa' \neq\kappa } i ^ {L' - L} ~ {\int_ 0 ^ {R_... ...ha' }( {(R_\kappa)}, {R_{\kappa'}} ) f_{\alpha'} (R_{\kappa'}) d R_{\kappa'} }}$$ + (30)

where the kinetic energy term is

$$T_{\kappa L} (R) = - {\hbar^2 \over {2 \mu_\kappa}} ~ \left ( {d^2 \over dR^2} - {{L(L+1)} \over R^2} \right ) ,$$ (31)

$U_\kappa (R_\kappa)$ is the diagonal optical potential with nuclear and Coulomb components, and R_m is a radius limit larger than the ranges of $U_\kappa (R_\kappa)$ and of the coupling terms. The equations (30) are in their most common form: they become more complicated when non-orthogonalities are included by the method of section 2.3. The $V^\Gamma_{\alpha :\alpha' } {(R_{\kappa'})}$ are the local coupling interactions of multipolarity $\Gamma$, and the $V_{\alpha :\alpha' }( R_\kappa , {R_{\kappa'}} )$ are the non-local couplings between mass partitions that arise from particle transfers.

For incoming channel $\alpha_0$, the $f_\alpha {(R_\kappa)}$ satisfy the boundary conditions

$$f_\alpha {(R_{\kappa'})}= _ {R_\kappa > R_m } {i \over 2} \left [... ..._{\alpha_0\alpha } H^{(+)}_{L {\eta _\alpha}} ( K_\alpha {(R_\kappa)}) \right ]$$ (32)

where $H^{(-)}_{L \eta}$ and $H^{(+)}_{L \eta}$ are the Coulomb functions [44] with incoming and outgoing boundary conditions respectively. The asymptotic kinetic energies are

$$E_{\kappa pt} = E + Q_\kappa - \epsilon_p - \epsilon_t$$ (33)

for excited state energies $\epsilon_p , \epsilon_t$ and Q-value $Q_\kappa$ in partition $\kappa$, and

$$K_\alpha = \sqrt {\left [ {2 \mu_\kappa} \over {\hbar^2} E_{\kappa pt} \right ] }$$ (34)

where $\mu_\kappa = A_{\kappa p} A_{\kappa t} / (A_{\kappa p} + A_{\kappa t})$ is the reduced mass in the channel with partition $\kappa$, and

$$\eta _\alpha = {{2 \mu_\kappa} \over \hbar^2} ~ { Z_{\kappa p} Z_{\kappa t} e^2 \over 2 K _\alpha }$$ (35)

is the Sommerfeld parameter for the Coulomb wave functions.

3.3 Single-nucleon states

If $\phi_{JM} ( \xi )$ is a core+particle bound state, then for coupling order $\left \vert ( \ell s)j,I;~JM \right\rangle ,$ the wave function is

$$\phi_{JM} ( \xi_c , {\bf r}) = \sum_{\ell jI} A_{\ell sj}^{jIJ} ~ \left [\phi_I (\xi_c) \varphi_{\ell sj} ({\bf r}) \right ] _{JM}$$ (36)

$$\sum_{ {\ell jI} ,{m \mu m_s m_\ell}} A_{\ell sj}^{jIJ} \left\lan... ...rangle Y_\ell^{m_\ell} ( \hat {\bf r})\phi_s^{m_s} \frac{1}{r} u_{\ell sjI} (r)$$ =

where $\xi_c$ = core internal coordinates, $\phi_{I mu} (\xi_c )$ = core internal state, $\phi_s^{m_s}$ = particle internal spin state, $u_{\ell s j I} (r)$ = particle core radial wave function, and $A_{\ell sj}^{jIJ}$ is the coefficient of fractional parentage.

3.3.1 Bound States

If the single-particle is bound at negative energy E around the core, then its wave function may be found as the eigen-solution of a given potential:

$$[ T_\ell (r) + V (r) + \epsilon_I - E ] u_{\ell sjI} (r) + \sum_{... ... ~\Gamma>0 } V^\Gamma_{\ell sjI : \ell' s j' I' } (r) u_{\ell' s j' I'} (r) = 0$$ (37)

with boundary conditions $u_{\ell sjI} (0) = 0$ and, as $r \geq R_m$, of $u_{\ell sjI} (r) \propto W_{\ell} (k_I r)$ where $W_\ell (\rho)$ is the Whittaker function and $k_I^2 = { 2 \mu \vert E - \epsilon_I \vert / \hbar^2 }$ is the asymptotic wave number.

If the core cannot be excited, then these coupled equations reduce to one uncoupled equation, but solving this equation can still be regarded as a special case of the coupled bound state problem. Eigen-solutions can be found by solving either for the bound state energy E, or by varying the depth of the binding potential. In general, however, we can choose to vary any multipole of any part of the binding potentials (except the Coulomb component), so one method of solving the full coupled bound-state problem will be outlined below.

To define the phase ($\pm 1$) of the overall wave function, some convention has to be adopted. One component (say that around a core I=0 state) can be set to either positive towards the origin ( $r\rightarrow 0$), or positive towards large distances ( $r\rightarrow \infty$). The former choice is made in the FRESCO code, following the Mayer-Jensen phase convention, which is also used for harmonic oscillator wave functions in many structure calculations.

3.3.2 Solution of the Coupled-Channels Eigenvalue Problem

When, for example, the problem is to find the bound state of a particle in a deformed potential, then several channels with different angular momenta will be coupled together. There are various techniques for calculating the wave functions of these bound states: for a review see ref. [10]. The Sturmian expansion method [45] can be used, or the coupled equations can be solved iteratively. The Sturmian method has the advantage that all solutions in the deformed potential are found, where sometimes the iterative method has difficulty in converging to a particular solution if there are other permitted solutions near in energy. The iterative method has the advantage that the radial wave functions (once found) are subject only to the discretisation error for the Schrödinger's equation, and are not dependent on the (time-consuming) diagonalisation of large matrices, often of the order of 100 or more. As they satisfy the correct boundary conditions independently of the size of a basis-state set, the radial wave functions of the iterative method therefore more accurately reflect the details of the coupling potentials and of the core excitation energies. As nuclear reactions are often confined to the surface region, it is important to satisfy the exterior boundary conditions as accurately as possible.

A method for solving the uncoupled eigenstate problem has to be included in a reaction code in any case, and since it can be generalised as described in this section to solving the coupled problem, it seems a worthwhile facility to have available. Bound states from a previous Sturmian solution can still be included as explicit linear combinations of the single particle (uncoupled) basis states used in the Sturmian expansion.

The general problem of finding eigen-solutions of a set M coupled-channels equations can be represented as the problem of finding $\lambda$ such that the equations

$$\left [ {d^2 \over dr^2} - {\ell_i (\ell_i + 1) \over r^2} \right... ...\sum _ {j=1} ^M \left [ U_{ij} (r) + \lambda V_{ij} (r) \right ] \psi_j (r) = 0$$ (38)

with boundary conditions

$$\psi_i (R) a_i W_{\ell_i \eta_i} (k_i R)$$ = (39)

$$\psi_i (R + \delta R) a_i W_{\ell_i \eta_i} (k_i (R + \delta R))$$ = (40)

$$\psi_i (0)$$ = 0 (41)

(with $k_i^2 \equiv\kappa^2_i + \theta \lambda$ and $\eta_i \equiv {nu_i} / (2 k_i)$) for given partial waves $\ell_i ,$ fixed potentials U_ij (r), variable potentials V_ij (r), matching radius R, and Coulomb proportionality constants $\nu_i$. The energy constants $\kappa^2_i$ are the asymptotic components of the diagonal U_ii (r), and $\theta$ is the asymptotic component of the diagonal V_ii (r) (assumed all equal).

The solution begins by constructing the trial integration functions for a trial value of $\lambda$, on either side of an intermediate matching point $r = \rho$:

$f_{i;j}^{\rm in} (r)$ by integrating r from h to $\rho$,
starting with $f_{i;j}^{\rm in} (h) = \delta_{i,j} ~ h^{\ell_i + 1} / (2 \ell_i + 1)!!$, and

$f_{i;j}^{\rm out} (r)$ by integrating r from R in to $\rho$,
starting with $f_{i;j}^{\rm in} (R) = \delta_{i,j} ~ W_{\ell_i \eta_i} (k_i (R + \delta R)).$

The intermediate point $r = \rho$ should be chosen where the wave functions are oscillatory, to avoid having to integrate outwards in the classically forbidden region.

The solution is therefore

$$\psi_i (r) = \left \{ \begin{array}{ll} \sum _ j b_j f_{i;j}^{\rm... ... _ j c_j f_{i;j}^{\rm out} (r) \mbox{ for } r \geq \rho , \end{array} \right. ~$$ (42)

and the matching conditions are the equality of the two expressions and their derivatives at $r = \rho$. The normalisation is still arbitrary, so we may fix c₁ = 1. In general the equations (38) have no solution as $\lambda$ is not exactly an eigenvalue. The method therefore uses the discrepancy in the matching conditions to estimate how $\lambda$ should be changed to $\lambda+\delta\lambda$ to reduce that discrepancy, and iterates this process to reduce $\delta\lambda$.

Thus at each iteration we first solve as simultaneous equations the 2M-1 of the matching conditions

$$\sum _ j b_j f_{i;j}^{\rm in} (\rho) \sum _ j c_j f_{i;j}^{\rm out} (\rho) \mbox{ for all } i$$ = (43)

$$\sum _ j b_j f_{i;j}^{\rm in} (\rho)' \sum _ j c_j f_{i;j}^{\rm out} (\rho) \mbox{ for all } i \neq 1$$ = (44)

along with c₁ =1 for the 2M unknowns b_i , c_i. If the function $\psi_i (r)$ is then constructed using equation (42), there will be a discrepancy as

$$\psi'_ {\rm in} \equiv\psi_1 (r) \vert _ {r< \rho} ~~~ \neq ~~~ \psi'_ {\rm out} \equiv\psi_1 (r) \vert _ {r> \rho} ,$$ (45)

and this difference will generate $\delta\lambda$ via

$$\delta \lambda \sum _ {ij} ^ M \int_ 0 ^ R \psi_i (r) V_{ij} (r)\psi_j (r) dr =\psi_1 (\rho) [ \psi'_ {\rm out} -\psi'_ {\rm in} ].$$ (46)

It is necessary while iterating in this manner to monitor the number of nodes in one or more selected components of the wave function, as in general a given potential will have different eigensolutions with different numbers of radial nodes. When the iterations have converged to some accuracy criterion on the size of $\delta\lambda$, the set of wave functions can be normalised in the usual manner:

$$\sum _ i ^ M \int _ 0 ^ \infty \vert \psi_i (r)\vert ^2 dr = 1$$ (47)

and perhaps some of the components i omitted if their contribution to this norm is below some preset threshold.

3.3.3 Continuum States

If the initial and/or final single-particle states of a transfer step are unbound $E - \epsilon> 0$, the use of single energy eigenstates $\phi_k (r)$ will result in calculations of the transfer form factors which will not converge, as the continuum wave functions do not decay to zero as $r\rightarrow \infty$ sufficiently fast as to have square norms. One way [39], [40] of dealing with this divergence is to take continuum states not at a single energy, but averaged over a range of energies. These `bin' states that result are square integrable, and if defined as

$$\Phi (r) \sqrt {2 \over \pi N} ~~ \int_ {k_1} ^ {k_2} w(k)\phi_k (r) dk$$ = (48)

$$\mbox{ with } N \int_ {k_1} ^ {k_2} \vert w(k)\vert ^ 2 dk$$ = (49)

for some weight function w(k), then they are normalised $\langle\Phi \vert\Phi\rangle = 1$ provided a sufficiently large maximum radius for r is taken, and that the $\phi_k$ are eigensolutions of a potential which is energy-independent. They are orthogonal to any bound states, and are orthogonal to other bin states if their energy ranges do not overlap. The construction can be easily generalised to give coupled-channels bin wave functions.

The weight function w(k) is best chosen ([40] p. 148) to include some of the effects known to be caused by the variation of $\phi_k (r)$ within the bin range $k_1 \leq k \leq k_2$. If $w(k) = \exp ( - i \delta_k )$, where delta_k is the scattering phase shift for $\phi_k (r)$, then it includes the effects of the overall phase variations of $\phi_k$, at least in the DWBA limit. If, however, $w(k) = \exp ( - i \delta_k ) \sin {\delta_k} \equiv T_k^*$, where T_k is the T-matrix element for $\phi_k (r)$, then it includes in addition a scale factor which is useful if the | T_k | varies significantly, as it does, for example, over resonances. Both choices result in a real-valued wave function $\Phi (r)$ (for real potentials), which is computationally advantageous.

If the maximum radius (R_m say) is not sufficiently large, then the wave functions $\Phi$ will not be normalised to unity by the factors given in equation (48). The rms radius of a bin wave function increases as the bin width k₂ - k₁ decreases, approximately as 1/(k₂ - k₁). These bin constructions can be used to describe the narrow resonant wave functions of say the 3⁺ state in ⁶Li, or the 7/2^- state in ⁷Li, but these states will require a large limiting radius R_m unless the w(k) = T_k^* weighting factor is used to emphasise the increase in the interior wave function over the resonance. The 3⁺ state in ⁶Li at 0.71 MeV, for example, for which the resonance width is approximately 40 keV, yields the normalisations shown in 3.3.3. It can be seen that without a scale factor which emphasises the resonance peak, very large radii R_m will be needed to obtain unit normalisation.

Figure 1: Normalisations of a continuum bin state. For this 3⁺ state in ⁶Li at 0.71 MeV, Saxon-Woods potentials were used with V = 77.05 MeV, R = R_c = 1.2 * 4^1/3 fm., and a = 0.65 fm.

3.4 Two-particle bound states

3.4.1 Centre-of-mass coordinates

If $\phi_{JT} (\xi_c , {\bf r}, {{\mbox{\boldmath$\rho$}}} )$ is a two-particle bound state with total spin J and isospin T, then for coupling order $\left \vert \{~L , (\ell , (s_1 s_2) S ) j \} J_{12} , I ;~J \right\rangle$ we have

$$\phi_{J M} = \sum_{\begin{array}{ccc} L \ell S\\ j J_{12} I\end{array}} A_{L j J_{12}}^{J_{12} I J} ~ \phi_{I \mu_I} (\xi_c) . \phi_{s_1}... ...ell^\mu ( \hat {{\mbox{\boldmath$\rho$}}} ) ~ {1 \over r \rho} u_{12} (r, \rho)$$

$$\langle J_{12} M_{12} I \mu_I \vert J M\rangle \langle L \Lambda ... ...a \vert j m_{12}\rangle \langle s_1 \sigma_1 s_2 \sigma_2 \vert S \Sigma\rangle$$ (50)

where A_{L j J12}^{J₁₂ I J} is the coefficients of fractional parentage, and $\phi_{s_1}^{\sigma_1} ~\phi_{s_2}^{\sigma_2}$ are the intrinsic spins of the two particles.

Note that two neutron transfer can be viewed as the transfer of a `structured particle' $( \ell , (s_1 s_2) S ) j$, and then becomes similar to single-particle transfers of above.

The radial wave function $u_{12} (r, \rho)$ can be given either as a cluster product of single-particle wave functions $u_{12} (r, \rho) =\Phi_L (r)\phi_\ell (\rho),$ or input directly as a two-dimensional distribution e.g. from a Faddeev bound-sate calculation, or calculated from the correlated sum of products of single-particle states, as in the next section.

3.4.2 Independent Coordinates

Two-particle states from shell-model calculations or from Sturmian-basis calculations [11], and are then usually described by means of the $\vert {\bf r}_1 , {\bf r}_2\rangle$ coordinates, and then transformed internally into the centre-of-mass coordinates $\vert {\bf r}, {{\mbox{\boldmath$\rho$}}} \rangle$ of equation (50) using ${\bf r}_i = x_i {\bf r}+ y_i {{\mbox{\boldmath$\rho$}}}$. For equal mass particles, x₁ = x₂ = 1, and $y_1 = - y_2 =\frac{1}{2}$. The second description is as

$$\varphi_{12} ( {\bf r}_1 , {\bf r}_2 ) \sum _ i c_i ~ \left \vert ( \ell_1 (i), s_1 )j_1 (i), ( \ell_2 (i), s_2 )j_2 (i);~ J_{12} T \right\rangle$$ = (51)

$$\rightarrow %%\vert _{\mbox{ coordinate transformation}} \sum _... ...s_1 s_2 )S )j ; J_{12} T \right\rangle \phi^{J_{12} T,i}_{L(\ell S)j} (r, \rho)$$ (52)

The transformation of the i'th component in the cluster basis is

$$\phi^{J_{12} T,i}_{L(\ell S)j} (r, \rho) \left\langle L , ( \ell , (s_1 s_2 )S )j ;J_{12} T \right \vert \... ...( \ell_1 (i), s_1 )j_1 (i), ( \ell_2 (i), s_2 )j_2 (i);~ J_{12} T \right\rangle$$ = (53)

$$\times\langle \left [ Y_L (\hat {\bf r}) Y_\ell (\hat \rho) \righ... ...} ({\bf r}_1)\varphi_{\ell_2 s_2 j_2} ({\bf r}_2) \right ] _ {J_{12} T} \rangle$$ (54)

where (suppressing the i indices for clarity)

$$\langle L , ( \ell , (s_1 s_2 )S )j ;J_{12} T \vert ( \ell_1 , s_1 )j_1 , ( \ell_2 , s_2 )j_2 ;~ J_{12} T\rangle =$$

$$\sum _ \lambda \hat \lambda \hat S \hat {j_1} \hat {j_2} \left ( \begin{array}{c... ... }} ~ \hat j \hat \lambda W(L \ell J_{12} S; \lambda j) (-1)^{\ell+L-\lambda} .$$ (55)

The radial overlap integral can be derived by means of harmonic-oscillator expansions [12], with the Bayman-Kallio expansion [13] or using the Moshinsky solid-harmonic expansion[71]. This last method gives

$$K^\lambda_{\ell L: \ell_1 \ell_2} (r, \rho) \langle \left [ Y_L (\hat{\bf r}) Y_\ell (\hat{{\mbox{\boldmath$\... ...phi_{\ell_1} ({\bf r}_1)\varphi_{\ell_2} ({\bf r}_2) \right ] ^ \lambda \rangle$$ =

$$\sum_{n_1 n_2} ~ \left ( \begin{array}{c} 2 \ell_1 +1\\ 2n_1\end... ...(x_1 r)^{\ell_1 - n_1} (y_1 \rho)^{n_1} (x_2 r)^{n_2} (y_2 \rho)^{\ell_2 - n_2}$$ =

$$\times \sum _ Q {\bf q}_{\ell_1 \ell_2}^Q (r, \rho) ~ (2Q+1)~ \ha... ...ell_1} \hat {\ell_2} \hat {\ell_1 - n_1} \hat {\ell_2 - n_2} ~ \hat L \hat \ell$$

$$\times \sum_{\Lambda_1 \Lambda_2} \left ( \begin{array}{ccc}\ell_... ...right ) \left ( \begin{array}{ccc}\Lambda_2&\ell&Q\\ 0&0&0\end{array} \right )$$

$$\times (-1)^{\ell_1 + \ell_2 +L+ \Lambda_2} (2 \Lambda_1 + 1) (2 ... ...n_1 & \ell_2-n_2& \Lambda_2\\ \ell_1 & \ell_2 & \lambda \end{array} \right ) .$$ (56)

where $\left ( \begin{array}{c}a \\ b \end{array}\right )$ is the binomial coefficient (see Appendix A).

The kernel function ${\bf q}_{\ell_1 \ell_2}Q (r, \rho)$ which appears in this expression is the Legendre expansion of the product of the two radial wavefunctions in terms of u, the cosine of the angle between ${\bf r}$ and ${{\mbox{\boldmath$\rho$}}}$:

$${\bf q}^Q_ {\ell_1 , \ell_2} (r, \rho) =\frac{1}{2} \int_ {-1}^{+... ..._1) \over r_1}^{\ell_1 +1} {u_{\ell_2} (r_2)\over r_2}^{\ell_2 +1} ~ P_Q (u) du$$ (57)

3.5 Scattering Amplitudes

The Rutherford amplitude for pure Coulomb scattering (with no $e^{2i \sigma_0}$ factor) is

$$F_c ( \theta ) = - {\eta \over 2k} ~~ { \exp (-2 i \eta \ln(\sin \theta /2)) \over \sin^2 \theta /2}$$ (58)

The Legendre coefficients for the scattering to the projectile state J'p and target state J't from initial projectile state J_p and target state J_t are given by

$$\sum_{L,J,J' ,J_T } \langle L0 J_p m \vert Jm\rangle\langle Jm J_t M \vert J_T M_T\rangle$$ A^L'_{m' M' ;mM} =

$$\langle L' M_{L'} J'_ p m' \vert J' M_{L'} + m'\rangle \langle J' M_{L'} +m' J'_ t M' \vert J_T M_T\rangle$$

$${4 \pi \over k} \sqrt {\frac{k'}{\mu'}\frac{\mu}{k}} e^{i( \sigma_L - \sigma_0 )} e^{i( \sigma'_ {L'} - \sigma'_ 0 )}$$

$$\left ( {i \over 2} \right ) \left [ \delta_{\alpha ,\alpha'} - S^{J_T}_{\alpha ,\alpha'} \right ] \sqrt {{2L+1 \over {4 \pi}}}~ Y_c (L' ,M_{L'})$$ (59)

where Y_c (L,M) is the coefficient of $P_L^{\vert M\vert} (\cos \theta ) e^{iM\phi}$ in $Y_L^M (\theta ,\phi )$, $\sigma_L = \arg \Gamma (1 + L + i \eta )$ is the Coulomb phase shift, $\alpha'$ refers to the primed values $L' J'_ p J'_ t k' \mu'$ etc., and $\alpha$ refers to the unprimed values $L J_p J_t k \mu$.

For each outgoing channel J'p , J't, we may then calculate the angular-dependent scattering amplitudes

$$f_{m' M' : mM} (\theta) = \delta_{J_p , J'_ p} \delta_{J_t , J'_ t} F_c (\theta) + \sum_{L'} A_{m' M' : mM}^{L'} P_{L'}^{m' +M' -m-M} (\cos \theta)$$ (60)

in terms of which the differential cross section is

$${d \sigma(\theta) \over d \Omega} = {1 \over (2J_p + 1)(2J_t + 1) } \sum_{m' M' m M} \left \vert f_{m' M' : mM} (\theta) \right \vert ^2 .$$ (61)

The near-side and far-side decompositions [14] of this cross section are defined by the same process, with P^M_L (u) replaced by $\frac{1}{2}[P^M_L (u)\pm 2i / \pi Q^M_L (u)]$ respectively. The Coulomb scattering of equation (58) is included in the near-side component [15].

The spherical tensor analysing powers T_kq describe how the outgoing cross section depends on the incoming polarisation state of the projectile. If the spherical tensor $\tau_{kq}$ is an operator with matrix elements

$$(\tau_{kq})_{m m''} = \sqrt{2k+1}\langle J_p m k q \vert J_p m''\rangle,$$

we have

$$T_{kq} (\theta) {Tr ({\bf f} \tau_{kq} {\bf f}^+)} \over {Tr ({\bf f} {\bf f}^+)}$$ = (62)

$$\hat{k} { \sum_{m' M' m M} f_{m' M' : m M} (\theta)^* \langle J_p... ...m'' M} (\theta) \over \sum_{m' M' m M} \vert f_{m' M' : mM} (\theta) \vert ^2 }$$ = (63)

The polarisations in the `transversity frame' [16] are then

$$\sqrt 2 i T_{11}$$ ^T T₁₀ = (64)

$$-\frac{1}{2}(T_{20} + \sqrt 6 T_{22})$$ ^T T₂₀ = (65)

$$-\frac{1}{2}( \sqrt 3 i T_{31} + \sqrt 5 i T_{33} ) .$$ ^T T₃₀ = (66)

The S-matrix elements can also be used to directly calculate the integrated cross sections

$$\sigma = \int_ {4 \pi} {~ {d \sigma(\theta) \over d \Omega} d \Omega }$$ (67)

to give

$$\sigma = {1 \over k^2} {k' \over \mu'} {\mu \over k} ~ {4 \pi \ov... ...ha ~\alpha'} (2J_T + 1) \left \vert S^{J_T}_{\alpha ,\alpha'} \right \vert ^2 .$$ (68)

The fusion cross section is defined as that amount of flux which leaves the coupled-channels set because of the imaginary parts of the optical potentials. If the incoming projectile is not spherical, then the fusion rate will depend on its orientation, and hence on the magnetic substate quantum number m. One can therefore define the fusion polarisation as the distribution $\sigma^{\rm fus}_m$, which can be calculated from the S-matrix elements as

$$\sigma^{\rm fus}_m {\pi \over k^2} \sum_{J_T \geq m} (2J_T + 1)$$ =

$$\times { \left [ 1 - {1 \over 2J_t + 1} \sum_{J'_ p L' M \omega} ... ...gle e^{i \sigma_L} S^{J_T \omega}_{\alpha ,\alpha'} \right \vert ^ 2 \right ] }$$ (69)

where $\omega$ is the parity ($\pm 1$) of the coupled-channels set for each total angular momentum J_T

Partial Wave Interpolation: Heavy ion reactions typically involve a range of partial waves L up to several hundred or more, especially when Coulomb excitations dominate the highest partial waves. In such cases it is often advantageous to solve the coupled channels sets (30) for, say, every n'th value of J_T, and interpolate the intermediate values. Different values of n can be used in different reaction regions: n can be small (1 or 2) for the grazing partial waves, and up to 5 or 10 for the Coulomb-dominated peripheral processes, and can be adjusted for the required balance between speed and accuracy.

This interpolation may be performed on the S-matrix elements themselves, or on the Legendre amplitudes of equation (59) In this second method (that used in ref. [34]), cubic spline interpolations are used. The main factor limiting the accuracy of this process is that the rate of change with J_T of the Coulomb phase shifts $\exp i(\sigma_L + \sigma'_ {L'})$ does not diminish as J_T increases. For that reason, it is advisable to interpolate not on the A^L' of equation (59), but on a $\tilde A^{L'}$ defined with a revised phase shift factor $\exp i(\sigma_L - \sigma'_ {L'}) .$ Since L and L' both tend to be near J_T, it is only the difference the phase shifts which limits the accuracy of the interpolation. It will therefore be more accurate for smaller projectile and target spins, and incoming and outgoing channels with similar Sommerfeld parameter $\eta$ (equation 35).