aa r X i v : . [ nu c l - t h ] A p r Adiabatic hyperspherical approach to large-scale nuclear dynamics
Y. Suzuki
1, 2 Department of Physics, Niigata University, Niigata 950-2181, Japan RIKEN Nishina Center, Wako 351-0198, Japan
We formulate a fully microscopic approach to large-scale nuclear dynamics using a hyperradiusas a collective coordinate. An adiabatic potential is defined by taking account of all possible con-figurations at a fixed hyperradius, and its hyperradius dependence plays a key role in governing theglobal nuclear motion. In order to go to larger systems beyond few-body systems, we suggest basisfunctions of a microscopic multicluster model, propose a method for calculating matrix elements ofan adiabatic Hamiltonian with use of Fourier transforms, and test its effectiveness.
I. INTRODUCTION
Atomic nuclei present a unique example of self-bound, finite quantum many-body systems. They not only exhibit avariety of excitation modes but also decay or fission into two or a few fragments. Exploring the excitation mechanismbased on single-particle, collective and clustering degrees of freedom is an interesting subject. Intrinsically differentshapes such as prolate-oblate may coexist or mix at close energies, leading to the so-called large-amplitude collectivemotion [1]. Spontaneous fission and sub-barrier fusion are also typical examples of the collective motion that involvesa large-scale change of the nuclear size [2, 3]. A fully microscopic description of their dynamics is still a long-standingchallenging problem.All of the above phenomena should in principle be described starting from a Hamiltonian of the system. What isoften performed is, however, to solve an equation of motion with some constraints [4] or to calculate energy surfacesassuming different shapes in order to look for a path along which the collective motion proceeds. In the case of adeep sub-barrier fusion an initial fragment decomposition is maintained for the whole fusion process and the relevantfusion potential is calculated as a function of the relative distance of the fragments. It is hard for that approach totake into account couplings with configurations corresponding to a different mass distribution of the fragments. Sincethe phenomena are very complicated, those approaches sound reasonable. However, neither the geometrical shapenor the relative distance between the fragments is a nuclear collective coordinate in a strict sense. A question thusarises of whether or not we can describe the large-scale dynamics by employing a true collective coordinate.The purpose of this paper is to make use of a hyperradius as a collective coordinate, and to step forward for aconsistent formulation of the large-scale dynamics together with the underlying collective potential. Most of theneeded ingredients are available in the literature. The hyperradius is a global coordinate that measures matter size,and it is widely used in three-body problems [5–8]. After pioneering work with the hyperspherical approach [9], itsextension to N -body systems has been proposed for solving various problems [10–14]. A common foundation of allthe hyperspherical approaches is that a total wave function of the system is expanded in terms of a product of thehyperradial and hyperangular functions. There are two types of realization for describing the hyperangular functions.One is to use hyperspherical harmonics [10, 11], and the other called an adiabatic hyperspherical approach is toemploy channel wave functions that are defined by diagonalizing the hyperangular part of the Hamiltonian [13, 14].The former has the advantage that the hyperspherical harmonics are well-known eigenfunctions of the angular partof the multi-dimensional Laplacian, but its use in a real problem is fairly complicated and so far limited to few-particle systems. Moreover, a convergence with that expansion is rather slow. See the first paper in Ref. [11] fordetails on the development and difficulty. The latter is widely used in atomic and molecular physics. Since anadiabatic hyperspherical potential defined there reflects the large-scale change of the system, we adopt the adiabatichyperspherical approach in what follows.The equation of motion in the hyperspherical method is the same independent of the number of particles in thesystem, which is an appealing feature of the hyperspherical approach. In spite of the various efforts, only smallsystems have so far been investigated mainly because calculating the matrix element of an adiabatic Hamiltonian isstill not trivial and solving its eigenvalue problem is hard for general N -body systems. Correlated Gauss functions(CG) are employed for studying cold atom physics and electron-positron systems in Refs. [13, 14], but their applicationis limited to few-body systems with the total orbital angular momentum L = 0 and 1. The use of harmonic-oscillatorshell-model wave functions is discussed in Ref. [11] for calculating the matrix element needed in the hypersphericalapproach. The oscillator basis is convenient for representing such one-centered configurations that are not highlyexcited from the ground state, but it is not flexible enough to cope with a description of large-scale change suchas, for example, the clustering and fragmentation. We instead need basis functions of cluster type to describe suchconfigurations. We attempt here an extension of microscopic multicluster wave functions [15, 16] used to describethe structure of light nuclei. In the multicluster model the intrinsic fragment wave functions are described withshell-model type configurations, while the relative motion among the fragments is described with the CG [17, 18]. Weapply a Fourier integral for evaluating the matrix element as it is applicable for any type of many-body basis function.The structure of the present paper is as follows. We define in Sect. II the hyperspherical coordinates and separatethe kinetic energy of the system into hyperradial and hyperangular parts. In Sect. III we define the eigenvalue problemof the adiabatic Hamiltonian and present the equation of motion for hyperradial functions. In Sect. IV we discussqualitative features of the adiabatic potential together with a separation of active and inactive degrees of freedom.In Sect. V we define basis functions of the multicluster model and give a method for calculating the matrix elementintegrated over the hyperangles together with examples for the overlap and kinetic energy. In Sect. VI we show howto extract the evolution of intrinsic shapes of the system as a function of the hyperradius. In Sect. VII we touch onan eigenvalue problem of the full Hamiltonian with a constraint of the mean-square matter radius in comparison withthe present approach. Conclusions are drawn in Sect. VIII. II. HYPERSPHERICAL COORDINATES
We start from defining the hyperradius for a general case consisting of K particles. The mass of the i th particle is A i in units of a suitable mass m . By denoting its position coordinate by R i , we define a set of Jacobi coordinates by X i = √ µ i (cid:16) R i +1 − A ...i i X j =1 A j R j (cid:17) , (1)where A ...i = P ij =1 A j and µ i is the reduced mass factor, µ i = A ...i A i +1 /A ...i +1 . The square of the hyperradius ρ is defined by ρ = K − X i =1 X i , (2)which is also rewritten in several ways as ρ = K X i =1 A i ( R i − R cm ) = K X i =1 A i R i − A ...K R = 1 A ...K K X j>i =1 A i A j ( R i − R j ) , (3)where R cm is the center-of-mass (cm) coordinate of the system, R cm = P Ki =1 A i R i /A ...K . Note that mρ is equalto the trace of the moment of inertia tensor of the system.It is straightforward to extend the above definition to an N -nucleon system. Protons and neutrons are assumed tohave an equal mass, the nucleon mass, which is taken as m . By denoting the nucleon’s position coordinate by r i , wedefine Jacobi coordinates as x i = √ µ i (cid:16) r i +1 − i i X j =1 r j (cid:17) (4)with µ i = i/ ( i + 1). Then ρ reads ρ = N − X i =1 x i = N X i =1 ( r i − R cm ) = N X i =1 r i − N R = 1 N N X j>i =1 ( r i − r j ) . (5)We often use a matrix notation. For example, x = ( x i ) stands for an N − N − × x stands for its row vector. The ρ is simply written as a scalar product, ρ = ˜ xx . It isclear that ρ is equally defined by any coordinates that are related to x by an orthogonal transformation. In fact ρ is independent of any choice of such coordinates.A measure of the nuclear size, ρ /N is an operator for the mean-square matter radius. Symmetric with respect tothe nucleons’ coordinates, ρ is a collective coordinate that has a unit of length. The other 3 N − d x = d x d x . . . d x N − = ρ d − dρd Ω , (6)where d is the dimension of the spatial coordinates excluding the cm coordinate d = 3( N − . (7)It is well known that the volume V d of a d -dimensional hypersphere with radius ρ = ρ is given by V d ≡ R ˜ xx ≤ ρ d x =( ρ √ π ) d / Γ( d/ V d is equal to R ρ ρ d − dρ R d Ω, the surface area of thehypersphere is Z d Ω = 2 √ π d Γ( d/ . (8)The volume element in the single-particle coordinates reads d r d r . . . d r N = N / d x d R cm .Let us introduce dimensionless coordinates ξ i by x i = ρ ξ i . They are subject to the constraint P N − i =1 ξ i = ˜ ξξ = 1.An explicit form of Ω may be constructed from the N − ξ i , but it is not needed in what follows. It shouldbe noted, however, that a variety of configurations or shapes of the nucleus correspond to different functions of Ω.The total kinetic energy T of the N -nucleon system, with its cm kinetic energy T cm being subtracted, is separatedinto hyperradial ( T ρ ) and hyperangular ( T Ω ) parts: T = − ~ m N X i =1 ∂ ∂ r i − T cm = − ~ m N − X i =1 ∂ ∂ x i = T ρ + T Ω (9)with T ρ = − ~ m (cid:18) ∂ ∂ρ + d − ρ ∂∂ρ (cid:19) = − ~ m ρ − ( d − ∂∂ρ ρ d − ∂∂ρ . (10)The hyperangular kinetic energy T Ω may be expressed as T Ω = ~ K (Ω) mρ , (11)where K (Ω) is the square of the grand angular momentum. An explicit form of K (Ω) is available in a recursive waytogether with the definition of Ω [6, 13].Suppose that the N -nucleon system develops into K fragments or clusters each of which has N i nucleons ( P Ki =1 N i = N ). It is convenient to divide ρ of Eq. (5) into two groups: ρ = ρ + ρ = K X i =1 ρ i + K X i =1 N i ( R i − R cm ) , (12)where R i is the cm coordinate of the i th fragment and ρ i is its squared hyperradius, ρ i = N i X j =1 ( r N ...i − + j − R i ) , ( N = 0) . (13)The first term ρ of Eq. (12) gives a measure of the sum of the squared matter radii of the fragments (each ρ i /N i is themean-square matter radius of the i th fragment), while the second term ρ is exactly the same as that of Eq. (3) with A i = N i , giving a measure of the spatial extension of the relative motion of the fragments. It is natural to arrange thecoordinates into cluster-internal and cluster-relative to describe the motion of the K fragments. The cluster-internalcoordinates, denoted by ( x , x , . . . , x N − K ), consist of a collection of Jacobi coordinates of each fragment, and thecluster-relative coordinates denoted by ( x N − K +1 , x N − K +2 , . . . , x N − ) are Jacobi coordinates as defined by Eq. (1).Clearly ρ is independent of the number of fragments into which the N -nucleon system develops. III. EQUATION OF MOTION IN ADIABATIC HYPERSPHERICAL EXPANSION
To solve a Schr¨odinger equation for the system in the hyperspherical method, a total wave function Ψ isusually expanded in terms of a complete set of the hyperspherical harmonics or K -harmonics Y λ (Ω): Ψ = ρ − ( d − / P λ χ λ ( ρ ) Y λ (Ω) [10, 11, 19, 20]. Here Y λ (Ω) is an eigenfunction of T Ω labeled by λ . The hyperradialfunctions χ λ ( ρ ) are determined from a set of coupled-channels equations. This method is successfully used in nuclearfew-body systems [5]. However, the number of hyperspherical harmonics needed to reach a converged solution be-comes very large at large ρ values [21]. Moreover, the coupling matrix elements between different Y λ (Ω) are of thesame orders of magnitude as the diagonal matrix elements especially for the Coulomb interaction, which also makesthe convergence slow.In low-energy phenomena, the hyperradial motion is expected to be slow compared to the hyperangular motion.Thus an adiabatic potential that takes account of all possible hyperangular motion at a fixed hyperradius gives insightinto the dynamics of the system’s evolution [22]. We adopt the adiabatic hyperspherical expansion method [6, 8, 12–14]used extensively in atomic and molecular physics.We define an adiabatic Hamiltonian H ad by H ad = T Ω + V + ~ ( d − d − mρ = H − T ρ + ~ ( d − d − mρ . (14)Here V is the total potential energy and H = T + V is the total Hamiltonian of the system. The nucleon-nucleoninteraction of V is assumed to be an effective interaction that contains no strong short-ranged repulsion. Assumingthat V contains no derivative of ρ , we solve an eigenvalue problem of the Hermitian operator H ad , H ad Φ ν ( ρ, Ω) = U ν ( ρ )Φ ν ( ρ, Ω) , (15)to obtain channel wave functions Φ ν ( ρ, Ω) and real adiabatic potentials U ν ( ρ ) that are labeled by ν . The quantumnumbers of H such as spin-parity J π are preserved as those of H ad , and the antisymmetry requirement on Ψ applieson Φ ν ( ρ, Ω) as well. Note that ρ appears parametrically in Eq. (15). At fixed ρ , all possible couplings amongvarious hyperangular configurations are taken into account to obtain U ν ( ρ ) and Φ ν ( ρ, Ω). The Φ ν ( ρ, Ω) form a set oforthonormal functions at each ρ , h Φ ν ′ ( ρ, Ω) | Φ ν ( ρ, Ω) i Ω = δ ν,ν ′ , (16)where h . . . i Ω indicates that the integration is carried out over Ω with ρ being fixed. Actually, the Φ ν ( ρ, Ω) also containspin and isospin coordinates that have to be integrated, but they are omitted for the sake of simplicity. Apparently U ν ( ρ ) contain the minimum ‘centrifugal potential’, ~ ( d − d − / mρ , for N ≥ K (Ω) vanishes.The Schr¨odinger equation, H Ψ = E Ψ, is solved by expanding Ψ in terms of Φ ν ( ρ, Ω):Ψ = ρ − ( d − / X ν f ν ( ρ )Φ ν ( ρ, Ω) . (17)The normalization of Ψ is P ν R ∞ | f ν ( ρ ) | dρ = 1 for a bound state. The hyperradial functions f ν ( ρ ) are determinedby solving a set of coupled-channels equations, (cid:20) − ~ m d dρ + U ν ( ρ ) − E (cid:21) f ν ( ρ ) − ~ m X ν ′ (cid:20) P νν ′ ( ρ ) ddρ + Q νν ′ ( ρ ) (cid:21) f ν ′ ( ρ ) = 0 , (18)with non-adiabatic coupling terms P νν ′ ( ρ ) = h Φ ν ( ρ, Ω) | ∂∂ρ Φ ν ′ ( ρ, Ω) i Ω , Q νν ′ ( ρ ) = h Φ ν ( ρ, Ω) | ∂ ∂ρ Φ ν ′ ( ρ, Ω) i Ω . (19)Equations (15), (18), and (19) give a microscopic description of the large-scale dynamics.A unique advantage of the adiabatic hyperspherical approach is that both lower and upper bounds to the exactlowest energy of H are readily obtained [23, 24]. As shown in Appendix A, we have P νν ( ρ ) = 0 . (20)Its differentiation with respect to ρ leads to ¯ Q νν ( ρ ) + Q νν ( ρ ) = 0 , (21)where ¯ Q νν ( ρ ) = h ∂∂ρ Φ ν ( ρ, Ω) | ∂∂ρ Φ ν ( ρ, Ω) i Ω (22)is non-negative, and consequently Q νν ( ρ ) ≤
0. The potential defined by W ν ( ρ ) = U ν ( ρ ) − ~ m Q νν ( ρ ) (23)always satisfies W ν ( ρ ) ≥ U ν ( ρ ). The lowest eigenvalue E obtained by truncating Eq. (18) to a single-channel equationwith the lowest adiabatic potential U ( ρ ) or W ( ρ ) gives a lower or upper bound to the exact lowest energy of H . SeeAppendix A for details. Convergence of the solution of Eq. (18) is checked by increasing the number ν of channels.A time-dependent Schr¨odinger equation is convenient for studying how final configurations in e.g. few-body decayand sub-barrier fusion evolve from their initial states. The wave function at time t is assumed asΨ( t ) = ρ − ( d − / X ν f ν ( ρ, t )Φ ν ( ρ, Ω) . (24)Once f ν ( ρ,
0) are given, f ν ( ρ, t ) for t > i ~ ∂∂t f ν ( ρ, t ) = (cid:20) − ~ m ∂ ∂ρ + U ν ( ρ ) (cid:21) f ν ( ρ, t ) − ~ m X ν ′ (cid:20) P νν ′ ( ρ ) ∂∂ρ + Q νν ′ ( ρ ) (cid:21) f ν ′ ( ρ, t ) . (25) IV. HYPERRADIUS DEPENDENCE OF ADIABATIC POTENTIAL
The ρ -dependence of U ν ( ρ ) or W ν ( ρ ) governs how the nucleus responds to its change of size. The kinetic energy andthe Coulomb potential respectively give 1 /ρ and 1 /ρ contributions to U ν ( ρ ) at large ρ values. Short-range pairwisenuclear interactions give a ρ − n ( n ≥
3) contribution [25]. Let us focus on the lowest adiabatic potential with the samespin-parity J π as that of the ground state. U ( ρ ) has a minimum at ρ ≈ ρ min corresponding to the matter size ofthe ground state. As ρ decreases from ρ min , U ( ρ ) rises because of a loss of nuclear potential energy as well as anincrease in the kinetic energy. As ρ increases from ρ min , various configurations contribute to determining Φ ( ρ, Ω).Here, deformations, shell effects, couplings with different modes and so on participate in determining the adiabaticpotential. U ( ρ ) reaches a peak at some ρ value or may even have a couple of local peaks at different ρ values. As ρ increases further, U ( ρ ) approaches the lowest decay threshold of the nucleus.The above global feature of the adiabatic potential well corresponds to the decomposition of ρ in conformity witha formation of fragments or clusters. As shown in Eq. (12), the different decomposition of the fragments can betreated on an equal footing in the hyperspherical approach, which makes it possible to assess what configurationsplay an important role in determining the adiabatic potential. If one instead calculates a sort of adiabatic potentialor potential energy surface as a function of the relative distance between two fragments, there is no way to comparesuch potentials for different fragment decompositions because their relative distances have a different meaning.What fragment decompositions or configurations contribute to the lowest adiabatic potential clearly depends on ρ .The expectation value of H is a major contribution to the adiabatic potential (see Eq. (14)). We rewrite H accordingto the fragment decomposition: H = K X i =1 ( T i + V i ) + T rel + V rel , (26)where T i + V i is the intrinsic Hamiltonian of the i th fragment, T rel the kinetic energy of the relative motion among thefragments, and V rel denotes the potential energies acting between the nucleons belonging to the different fragments. V rel depends on both cluster-internal and cluster-relative coordinates, thus causing a coupling of the relative motionamong the fragments with their intrinsic motion. When ρ rel is so large compared to ρ in that the nucleon-nucleoninteractions of V rel can be neglected and only the leading term of the Coulomb potentials of V rel is retained, V rel reduces to V rel → V C rel = K X j>i =1 Z i Z j e | R i − R j | = e ρ rel C (Ω rel ) , (27)where Z i e is the charge of the i th fragment and ρ rel and Ω rel are the hyperradius and hyperangles constructed from thecluster-relative coordinates ( x N − K +1 , x N − K +2 , . . . , x N − ). With increasing ρ the intrinsic motion of each fragmentis stabilized toward its own ground state, while the configurations responsible for the relative motion are decoupledfrom the intrinsic motion. Both the coupling and decoupling of various degrees of freedom are naturally taken intoaccount in the hyperspherical approach. U ν ( R ) ( M e V )
100 200 300R (fm)00.10.20.3 U ν ( R ) ( M e V ) (a) (b) FIG. 1: The 10 lowest adiabatic potential curves of the three- α system with J π = 0 + that are taken from Ref. [26]. Thehyperradius denoted R here is defined by R = √ P i =1 ( R i − R cm ) , and the energy is measured from the three- α threshold.The solid line denotes the adiabatic potential dominated by the Be+ α channel, while the dashed lines denote the potentialsdominated by the three- α continuum channel. Panel (b) is an enlarged view of the potentials at large R where a number ofsharp avoided crossings successively appear. When there are several thresholds corresponding to different fragment decompositions, avoided crossings of theadiabatic potential energy curves may occur. As an example, we show the case of C that is described with a clustermodel of three α -particles [26]. The eigenvalue problem (15) for H ad is solved accurately, and an analysis of theadiabatic potentials clarifies how the contributions of the hyperangular kinetic energy, the nuclear potential and theCoulomb potential change as a function of ρ . Figure 1, taken from Fig. 2 of Ref. [26], displays the 10 lowest adiabaticpotential curves for J π = 0 + . The lowest potential U ( R ) has a minimum at R ≈ . C. Furthermore, the lowest potential reaches a broad peak around12 fm, corresponding to the second 0 + state of C, the Hoyle resonance state. The adiabatic potential indicated bythe solid line is dominated by the two-body Be+ α state and approaches the Be+ α threshold at large R , while theother potentials indicated by the dashed lines are all dominated by the three- α continuum states. As seen in Fig. 1(b), an avoided crossing begins to occur at R ≈
140 fm, which is because the three- α continuum state comes downclosely to the two-body Be+ α state. Since the avoided crossing actually occurs within a small range of R , it maybe hard to see it in the figure. Refer to Fig. 3 of Ref. [26] to confirm the crossing clearly. Since the Be+ α thresholdis higher than the three- α threshold, a number of avoided crossings successively appear below the adiabatic potentialindicated by the solid line. As is well known, the non-adiabatic coupling terms (19) may be singular especially whenthe avoided crossing is sharp, namely it occurs within a small range of ρ . In that case, a diabatic procedure is proposedfor accurately solving Eq. (18) [6, 27, 28]. The slow variable discretization method combined with a complex absorbingpotential makes it possible to solve Eq. (18) and to reproduce the energy and width of the Hoyle resonance in goodagreement with experiment [26].Let us speculate concerning the adiabatic potential curves of Cf that are crucially important for determining itsdecay mode. The ground state of
Cf decays mostly by an α -particle emission. The rest is a spontaneous fission(SF), emitting 3.7 neutrons on average. To make things simple, we approximate the SF as occurring through a singlechannel of Xe +
Ru + 4 n . The two decay modes contain different numbers of fragments, two in α + Cmand six in the SF, but the hyperspherical approach can treat both in a unified way. The threshold of α + Cm is6.2 MeV below the ground state of
Cf, whereas that of the SF is 200.4 MeV lower than the ground state. See theschematic diagram of Fig. 2. The lowest adiabatic potential U ( ρ ) approaches the SF threshold at large ρ . Abovethat threshold many U ν ( ρ ) curves, not drawn in Fig. 2, show up corresponding to the continuum states of the SFmode. A unique U ν ( ρ ) with the two-body α + Cm character appears high above the SF threshold. When movinginward from this asymptotic region, the Coulomb potential (27) produces a distinct difference between the two decaymodes. The charge factor Z Z of the SF mode is more than ten times larger than that of the α channel. Thus those U ν ( ρ ) curves that are dominantly contributed by the SF configurations rise up rapidly, while the U ν ( ρ ) curve of the α channel increases much more slowly. At the avoided crossing point ρ ac , the lowest curve U ( ρ ) comes very close tothat of the α curve, and for ρ < ρ ac the α channel makes a dominant contribution to U ( ρ ). With further decrease of ρ many different configurations begin to mix due to an increasing role of the nuclear interaction V . The U ( ρ ) reachesa barrier top around some point and reaches its minimum at ρ min corresponding to the matter radius of the groundstate of Cf. Though much more complicated than the C case, the gross feature of the adiabatic potential curvesof
Cf should have some similarity to those of C, and the decay branch of
Cf will be determined by solvingEq. (18). ρ min ρ ac A d i aba t i c po t en t i a l ρ FIG. 2: A schematic diagram of the adiabatic potential curves of
Cf as a function of hyperradius ρ . Energy is measuredfrom the ground state of Cf. The potential dominated by the α + Cm channel (solid line) goes to − . ρ values, while that by the Xe +
Ru + 4 n channel (dashed line) to − . ρ ≈ ρ ac . Thelowest adiabatic potential U ( ρ ) changes its dominant character from α + Cm to the
Xe +
Ru + 4 n channel around ρ ac . V. MULTICLUSTER APPROXIMATION AND INTEGRATION OVER HYPERANGLES
Solving Eq. (15) is of vital importance in the adiabatic hyperspherical approach. Its accurate solution is obviouslyvery hard except for few-body system. The difficulty is enhanced by the fact that the matrix element has to becalculated by integrating over Ω only. Some efforts have been made for extending to larger systems [11, 13, 14]. Wetake up this problem assuming the use of many-body wave functions that contain all the coordinates.Before discussing the eigenvalue problem (15), we note that a usual approach defines an adiabatic potential barrieror energy surface at a given ‘collective’ coordinate by searching for a minimum of V for various parameters thatcharacterize the nuclear density or shape [29]. This makes sense in that V is a major part of H ad , and because,since V is a function of ρ and Ω, its minimum gives information on the most important Ω values contributing to thelowest adiabatic potential. As mentioned before, the adiabatic hyperspherical approach can go beyond that by takingaccount of various couplings with different degrees of freedom.Let us assume that the channel wave function Φ ν ( ρ, Ω) at a given ρ is expanded in terms of suitable basis functions φ i ( x ): Φ ν ( ρ, Ω) = X i C νi ( ρ ) φ i ( x ) . (28)Equation (15) is then reduced to the following generalized eigenvalue equation for determining the coefficients C νi ( ρ )and the adiabatic potential U ν ( ρ ): X j [ H ij ( ρ ) − U ν ( ρ ) B ij ( ρ )] C νj ( ρ ) = 0 , (29)where H ij ( ρ ) and B ij ( ρ ) are adiabatic Hamiltonian and overlap matrices defined by H ij ( ρ ) = h φ i ( x ) | H ad | φ j ( x ) i Ω , B ij ( ρ ) = h φ i ( x ) | φ j ( x ) i Ω . (30)We include only those basis functions that give a c -number ρ for the expectation value of the squared hyperradiusoperator ˜ xx : h φ i ( x ) | ˜ xx | φ i ( x ) ih φ i ( x ) | φ i ( x ) i ≈ ρ . (31)We face two problems. One is what basis functions we use for φ i ( x ). The other is how to calculate the matrixelement in Eq. (30). The first one is crucially important for assessing the quality of Φ ν ( ρ, Ω) and U ν ( ρ ). Though it isdifficult to give a general answer, our ansatz is to employ a microscopic multicluster approximation [15, 16]. This isbecause, as mentioned in Sects. I and IV, the structure change we are interested in includes a variety of configurationsranging from one-centered shell-model wave functions to those with a few fragments or subsystems. A general formof the multicluster wave function containing K fragments reads φ ( K ) ( x ) = A{ Ψ ( z )Ψ ( z ) · · · Ψ K ( z K ) χ ( x N − K +1 , . . . , x N − ) } , (32)where A is an antisymmetrizer, Ψ i an antisymmetrized intrinsic state of the i th fragment containing N i nucleonsand χ is the relative motion function for the fragments. The cluster-internal coordinates ( x , x , . . . , x N − K ) areabbreviated as ( z , z , . . . , z K ), where e.g. z stands for the first N − x , x , . . . , x N − ). Thespin-isospin coordinates are again suppressed. In general Ψ i may represent not only the ground state of the fragmentbut also its excited state. The quantum numbers for characterizing Ψ i are omitted. The coupling of the angularmomenta of the Ψ i s and χ to a total angular momentum JM is implicitly understood in Eq. (32). We presume φ i ( x )to belong to the space spanned by { φ (1) ( x ) } + { φ (2) ( x ) } + . . . (33)Note that any states in { φ ( K ) ( x ) } are in general nonorthogonal to each other even when they belong to different K subspaces. The questions of what intrinsic states of the fragments are important and what K subspaces haveto be included depend on a given system and energy range of interest. To proceed further, we assume that Ψ i isapproximated by harmonic-oscillator shell-model wave functions, while χ is described well with a superposition ofGauss functions [17, 18, 30] as developed in few-body problems.We have to calculate a matrix element for some operator O ( x ), O ( ρ ) = h φ i ( x ) | O ( x ) | φ j ( x ) i Ω , ρ = ρ , (34)by integrating over Ω at fixed ρ , say ρ . The calculation of the matrix element of T ρ in H ad can be aided with use ofthe identity ∂∂ρ φ j ( x ) = 1 ρ (cid:16) N − X i =1 x i · ∂∂ x i (cid:17) φ j ( x ) . (35)See Appendix B for an example. In calculating the matrix element of H , the cluster-intrinsic term P Ki =1 ( T i + V i ) (seeEq. (26)) may be replaced by Hφ ( K ) ( x ) → (cid:16) K X i =1 E i (cid:17) φ ( K ) ( x )+ A{ ( T rel + V rel )Ψ ( z )Ψ ( z ) · · · Ψ K ( z K ) χ ( x N − K +1 , . . . , x N − ) } , (36)using the observed energy E i of Ψ i . This approximation looks reasonable and practically useful because any nuclearinteraction can not satisfactorily reproduce the saturation property of nuclear binding energies despite the fact thatreproducing the threshold energy for the fragment decomposition is important in the present approach.The second problem has so far been examined using integral transform techniques [14, 31]. We use a δ functiontechnique as in Ref. [14]. Using the expression for Dirac δ function δ ( ρ − ρ ) = 1 πρ Z ∞−∞ e iω (1 − ρ /ρ ) dω, (37)we can express O ( ρ ) as a Fourier transform of F ρ ( ω ): O ( ρ ) = 1 π Z ∞−∞ e iω F ρ ( ω ) dω, (38) F ρ ( ω ) = Z e − iω ˜ ξξ (cid:0) φ i ( ρ ξ ) (cid:1) ∗ O ( ρ ξ ) φ j ( ρ ξ ) d ξ . (39)Note that x i is changed to ρ ξ i with a dimensionless variable ξ i . In Eq. (39) d ξ stands for d ξ d ξ . . . d ξ N − , wherethe integration range of each ξ i covers the whole three-dimensional space. Since e − iω ˜ ξξ = Q N − k =1 e − iω ξ k results in asimple modification of the basis function, F ρ ( ω ) can be calculated with a technique developed in microscopic clustermodels [32, 33].In some cases the Fourier integral (38) can easily be obtained by Cauchy’s integral formula that reduces to aresidue calculation. Whether or not we have a practical means for evaluating Eq. (34) for a general case depends onhow fast and accurately the Fourier integral is computed. For this aim we test the Whittaker cardinal series or theWhittaker-Shannon interpolation formula [34]: F ρ ( ω ) = ∞ X n = −∞ F ρ ( ω n )sinc πh ( ω − ω n ) , (40)where sinc x is the sinc function, sin x/x , and ω n = nh ( n = 0 , ± , ± , . . . ) is the grid of the sampling points. Theseries (40) is known to converge if F ρ ( ω ) is a band-limited function. Because sinc nπ = δ n, , the series is exact atall the sampling points. It is in fact an expansion in terms of orthogonal functions { sinc πh ( ω − ω n ) } that have theproperties: Z ∞−∞ sinc πh ( ω − ω n ) dω = h, Z ∞−∞ sinc πh ( ω − ω m ) sinc πh ( ω − ω n ) dω = hδ m,n , Z ∞−∞ e iω sinc πh ( ω − ω n ) dω = he iω n ( hπ < . (41)The third equation called the Dirichlet integral leads to an approximation for O ( ρ ): O ( ρ ) ≈ hπ M X n = − M F ρ ( ω n ) e iω n , (42)which is nothing but a trapezoidal rule for the integration. This result is due to the fact that the Fourier transformof the sinc function is the rectangular function and vice versa. To determine M , we need to know how fast F ρ ( ω )decreases as a function of ω . The mesh size h ( h < π ) is determined by examining how accurate the expansion is at,e.g. ω = ( n + ) h , the midpoint of ω n and ω n +1 .Other interpolations, e.g. a cubic spline interpolation may also be worthwhile testing because it leads to a simpleexpression for Eq. (38) and in addition the mesh size can be taken as piecewise variable. Once dF ρ ( ω ) /dω valuesat both boundaries of the interpolation are calculated, we can completely fix the interpolating function of the cubicspline.Since ω -dependence of F ρ ( ω ) is of practical importance, we examine it for the diagonal matrix elements ( φ i ( x ) = φ j ( x )) of O ( x ) = 1 and T Ω in a very schematic model. As the model, we employ CG ignoring the antisymmetryrequirement of the wave function and focus only on its spatial part. See Appendix B for some basic matrix elementswith the CG. For a spherical CG, exp( − ˜ x A x ), the positive-definite symmetric matrix A is set to Tr A − = ρ because of Eqs. (31) and (B2). We may choose A to be diagonal, A = ( a i δ ij ), as far as the diagonal matrix elementof O ( ρ ) is concerned.Our first choice for A is a uniform nuclear expansion, a i = a , leading to a hyperscalar Gaussian, exp( − ˜ x A x ) =exp( − aρ ). This function is totally symmetric and Ω-independent. By taking a as ( N − /a = ρ , the overlapmatrix element is (see Eq. (B3)) F ρ ( ω ) = (cid:18) (2 π ) N − (2 aρ + 2 iω ) N − (cid:19) / = (cid:18) πd/ iω (cid:19) d/ . (43)Clearly | F ρ ( ω ) | becomes very small if ω is significantly larger than d/
2. The Fourier transform (38) can be rigorouslycomputed in this case. If d/ F ρ ( ω ) has a pole of order d/ ω = id/
2, so that the integral is reducedto a residue calculation, yielding O ( ρ ) = 2 √ π d Γ( d/ e − d/ . (44)Even when d/ d/ iω = − t , O ( ρ ) is reduced to O ( ρ ) = i √ π d π e − d/ Z − d/ − i ∞− d/ i ∞ e − t ( − t ) − d/ dt. (45)By changing the integration path to the Hankel contour and using Hankel’s integral representation and Euler’sreflection formula for the gamma function, we find the above integral to be 2 π/i Γ( d/ ρ = ρ is exp( − aρ ) = e − d/ and hence O ( ρ )must be e − d/ R d Ω. We note that O ( ρ ) for O ( x ) = T Ω vanishes because the hyperscalar Gaussian is Ω-independentand, when acted on by T Ω , vanishes. This is also confirmed by using Eq. (B6).0 -14-12-10-8-6-4-2 0 0 400 800 1200 1600 2000 l og r( ω ) N = 4 η = 135 0 400 800 1200 1600 2000 ω N = 40
0 400 800 1200 1600 2000
N = 240
FIG. 3: The ratio r ( ω ) = | F ρ ( ω ) /F ρ (0) | for the overlap corresponding to N nucleons’ symmetric fission as a function of ω . ρ is set to η √ NR rms ( N ): Solid, dashed, and dotted lines correspond to η = 1 ,
3, and 5, respectively. See text for details. -14-12-10-8-6-4-2 0 0 400 800 1200 1600 2000 l og r( ω ) N = 4 η = 135 0 400 800 1200 1600 2000 ω N = 40
0 400 800 1200 1600 2000
N = 240
FIG. 4: The same as Fig. 3 but for the hyperangular kinetic energy.
The next example is a ‘symmetric fission’, that is, the nucleus fissions into two identical fragments with massnumber N/ ρ . Let R rms ( N ) denote the root-mean-square radius of a nucleus with mass number N , and set it equal to p / r N / ( r = 1 . A for the symmetric fission is chosen as a = a = . . . = a N − = a , and a and a N − are determined by the condition N − a + 1 a N − = 23 ρ , N/ − a = 23 ρ f , (46)where ρ f is fixed to p N/ R rms ( N/ N is changed to 4, 40, and 240, and for each N ρ istaken as ρ = η √ NR rms ( N ) ( η = 1 , , r ( ω ) for the overlap, where r ( ω ) = | F ρ ( ω ) /F ρ (0) | .In the case of N = 4, the fall-off of r ( ω ) is slow with increasing ω and η . For example, log r ( ω ) at ω = 5000 is − . , − . , − . η = 1 , ,
5, respectively. For N = 40, r ( ω ) rapidly drops to 10 − as a function of ω , butits decrease becomes slower for η = 5. This behavior is also valid for N = 240, and the decrease in ω becomes evenslower with increasing η . As shown in Fig. 4, the ratio log r ( ω ) for T Ω is very similar to that of the overlap. VI. EVOLUTION OF INTRINSIC SHAPES
It is interesting to know how an intrinsic shape of the nucleus changes with increasing ρ . When a decay or an SFis considered as a tunneling through a barrier, the shape will give insight into where the fragments are formed andhow they evolve during the passing through the barrier. The barrier is conventionally calculated by assuming somedensity distribution constrained with shape or deformation parameters such as quadrupole and octupole [2, 29]. Suchdeformations are not observable, however. Our view is to reverse this approach. Since the nucleus should in principlepreserve the total angular momentum, it is not trivial to imagine the intrinsic shape in the space-fixed frame. Forexample, any state with L = 0 is spherical in that frame, but it can happen that such state is intrinsically deformedand rotates. As shown in Ref. [35], the intrinsic two- α structure of the rotational state of Be emerges from the wavefunction obtained by a quantum Monte Carlo calculation.Following the procedure of Ref. [35], we can get the intrinsic density or deformation indicated by e.g. the lowestchannel wave function Φ ( ρ, Ω). Since Φ ( ρ, Ω) is normalized as in Eq. (16), its square, P ρ (Ω) = | Φ ( ρ, Ω) | , gives the1probability density as a function of Ω at a given ρ . First, we generate many sampling points Ω , Ω , . . . , Ω M accordingto the distribution of P ρ (Ω) using the Metropolis-Hastings algorithm. Secondly, we define a body-fixed intrinsicframe for each Ω j = ( ξ j , ξ j , . . . , ξ jN − ) as follows. By using Eq. (4) together with R cm = 0, Jacobi coordinates x ji = ρ ξ ji ( i = 1 , , . . . , N −
1) specify the positions of N nucleons ( r j , r j , . . . , r jN ) in the space-fixed frame. Fromthese position vectors, we calculate the moment of inertia tensor I jαβ = N X i =1 r ji α r ji β , (47)where r ji α ( α = x, y, z ) is the Cartesian component of r ji . Diagonalizing the 3 × I j determinesthe principal moments of inertia, which define the axes of the intrinsic frame. For example, the axis is called x ′ , y ′ , z ′ in increasing order of the principal moment of inertia. The direction of the axis also has to be chosen consistently.By reading ( r j , r j , . . . , r jN ) as ( r ′ j , r ′ j , . . . , r ′ jN ) in reference to the intrinsic frame, we obtain the desired positioncoordinates of N nucleons in the intrinsic frame. Finally, accumulating these position coordinates over j = 1 , , . . . , M leads to the intrinsic single-particle density at ρ . Once the intrinsic density is obtained, it is easy to extract multipoledeformations. VII. EIGENVALUE PROBLEM OF HAMILTONIAN WITH RADIUS CONSTRAINT
It looks as though the adiabatic hyperspherical approach has some relationship to an eigenvalue problem of theHamiltonian with a constraint [4]. Let us attempt to find a solution of the Schr¨odinger equation by first constrainingthe expectation value of the squared hyperradius to a fixed ρ value. Suppose that the solution is expanded in termsof some basis functions: Ψ ρκ ( x ) = X i G κi ( ρ ) φ i ( x ) . (48)The constraint (31) is not necessarily imposed on φ i ( x ) itself, but we demand the solution we are looking for to satisfythe condition h Ψ ρκ ( x ) | ˜ xx | Ψ ρκ ( x ) ih Ψ ρκ ( x ) | Ψ ρκ ( x ) i = ρ . (49)Here κ is a label related to the Lagrange multiplier. The unknown coefficients G κi ( ρ ) and the energy eigenvalue E κ ( ρ )are determined from the following equation X j (cid:2) H ij − κ ( Q ij − ρ B ij ) − E κ ( ρ ) B ij (cid:3) G κj ( ρ ) = 0 , (50)where H , B , and Q are matrices defined by H ij = h φ i ( x ) | H | φ j ( x ) i , B ij = h φ i ( x ) | φ j ( x ) i , Q ij = h φ i ( x ) | ˜ xx | φ j ( x ) i . (51)Unlike Eq. (30), the above matrices are obtained by integrating over the whole coordinates. To determine thecoefficients G κj ( ρ ) from Eq. (50), the value of κ has to be given. Actually κ should be such that both Eqs. (49) and(50) are simultaneously met. Apparently Ψ ρ ′ κ ′ ( x ) and Ψ ρκ ( x ) are not orthogonal to each other even for ρ ′ = ρ .The next step is to use the generator coordinate method in which a solution Ψ for the Schr¨odinger equation isassumed as Ψ = X κ Z C κ ( ρ )Ψ ρκ ( x ) dρ. (52)The coefficients C κ ( ρ ) are determined from the Hill-Wheeler equation X κ Z h Ψ ρ ′ κ ′ ( x ) | H − E | Ψ ρκ ( x ) i C κ ( ρ ) dρ = 0 , (53)which should be satisfied for any ρ ′ and κ ′ values. An approximate solution to the Hill-Wheeler equation gives anupper bound to the ground-state energy. Note that the adiabatic hyperspherical approach gives both lower and upperbounds as discussed in Sect. III.2We refer to two interesting calculations with a constraint in comparison to the adiabatic hyperspherical approach.One is a Hartree-Fock-Bogoliubov calculation performed by constraining the mean-square radius, P Ni =1 r i /N , tostudy how self-conjugate nuclei fragment into α clusters [36]. As Eq. (5) indicates, this constraint is equivalent tothat of ρ provided the contribution of R to the squared radius remains a constant. The treatment of the cmmotion in Ref. [36] does not satisfy this condition as usual in a mean-filed model. It would be a challenge for themean-field approximation to cope with such diverse structure at large distances that is composed of different numbersof fragments. What should be further pursued at this moment is to establish the essential relationship between theadiabatic hyperspherical approach and ‘beyond mean-field’ calculations or configuration interaction calculations thatconstrain the mean-square matter radius.Another is a simultaneous study of both α + He reactions and the structure change of Be in a microscopic α + α + n + n model [37], in which a distance parameter between the two α -clusters is constrained. Since the motion of thetwo neutrons is restricted to either molecular or atomic orbits around the α -clusters, the main configurations includedare α + He and He+ He two-body types. The adiabatic energy surfaces are calculated within that approximation.An avoided crossing is treated by the generator coordinate method. As noted in Sect. IV, the relative distance of thefragments is not a collective coordinate. If one constrains ρ as the generator coordinate, it would be possible in thesame four-body model to take account of possible couplings with the Be+ n channel that is the lowest threshold of Be as well as the three- and four-body channels, Be+ n + n and α + α + n + n , that are open in the energy regiontreated in Ref. [37]. VIII. CONCLUSION
Stressing that the hyperradius is a collective coordinate, we have formulated a fully microscopic adiabatic hyper-spherical approach to large-scale nuclear dynamics. The equation of motion for hyperradial functions is universal,independent of the number of nucleons, and enables one to consistently treat the dynamics from confined nuclearmotion to relative motion among fragments in their asymptotic region. It is possible to describe in a unified waycases where the nucleus fragments into several channels. No spurious center-of-mass motion appears and couplingswith different degrees of freedom can naturally be taken into account. These properties are due to the fact thatboth the squared hyperradius and the kinetic energy are flexibly decomposed into cluster-internal and cluster-relativequantities responding to the fragment formation.The adiabatic potential as a function of the hyperradius plays a key role in the present approach. It is unambigu-ously defined solely by the Hamiltonian of the system, and there is no need to assume specific geometrical shapesor deformations to compute it. Conversely the shape or intrinsic density, if necessary, comes out after the adiabaticpotential is obtained or the equation of motion for the hyperradial functions is solved. The calculation of the adia-batic potential involves the integration over all the coordinates but the hyperradius. Expecting that a microscopicmulticluster model is a promising candidate for applying the present approach to larger systems, we have discussedthe use of Fourier transforms for evaluating the matrix elements needed to obtain the adiabatic potential. The ma-trix elements can be obtained in exactly the same way as the usual matrix elements needed in nuclear many-bodycalculations. A merit of the Fourier transform technique is its simplicity, and test calculations indicate that accurateevaluations of the matrix elements are feasible.Although the calculation of the adiabatic potential still requires much computer time for large systems, a realchallenge is whether we can provide large enough basis functions to cover important configurations for fixed ρ . Furtherdevelopments are certainly indispensable for a microscopic, realistic description of large-scale nuclear dynamics. Acknowledgment
The author is greatly indebted to H. Suno for many instructive discussions. He also thanks W. Horiuchi, K. M.Daily, and C. H. Greene for valuable communications. This work is supported in part by JSPS KAKENHI Grant No.24540261.
Appendix A: Lower and upper bounds
In this appendix we rigorously prove Eq. (20) and show that both lower and upper bounds to the ground-stateenergy are respectively obtained by solving single-channel equations.3The ground-state wave function may be expressed in the hyperspherical coordinates asΨ = ρ − ( d − / f ( ρ )Φ( ρ, Ω) (A1)with the normalization condition Z ∞ | f ( ρ ) | dρ = 1 , h Φ( ρ, Ω) | Φ( ρ, Ω) i Ω = 1 . (A2)The hyperradial function f ( ρ ) has to vanish at ρ = 0. The ground-state energy reads E exact = − ~ m Z ∞ f ∗ ( ρ ) h d f ( ρ ) dρ + 2 df ( ρ ) dρ P ( ρ ) + f ( ρ ) Q ( ρ ) i dρ + Z ∞ | f ( ρ ) | U ( ρ ) dρ, (A3)where P ( ρ ) = h Φ( ρ, Ω) | ∂∂ρ Φ( ρ, ω ) i Ω , Q ( ρ ) = h Φ( ρ, Ω) | ∂ ∂ρ Φ( ρ, Ω) i Ω ,U ( ρ ) = h Φ( ρ, Ω) | H ad | Φ( ρ, Ω) i Ω . (A4)From the normalization condition of Φ( ρ, Ω), we obtain ddρ h Φ( ρ, Ω) | Φ( ρ, Ω) i Ω = P ( ρ ) ∗ + P ( ρ ) = 0 . (A5)Thus P ( ρ ) must be pure imaginary or zero. If P ( ρ ) is not zero but pure imaginary, f ∗ ( ρ ) df ( ρ ) /dρ in Eq. (A3) mustalso be pure imaginary because E exact is real. With f ( ρ ) = g ( ρ ) + ih ( ρ ), where g ( ρ ) and h ( ρ ) are real functions, f ∗ ( ρ ) df ( ρ ) /dρ reads f ∗ ( ρ ) df ( ρ ) dρ = 12 ddρ { g ( ρ ) + h ( ρ ) } + i { g ( ρ ) h ′ ( ρ ) − h ( ρ ) g ′ ( ρ ) } , (A6)which leads to d { g ( ρ ) + h ( ρ ) } /dρ = 0. Thus g ( ρ ) + h ( ρ ) is a constant, and it must be zero because of f (0) = 0.Namely, f ( ρ ) vanishes identically, which can not be accepted. Using P ( ρ ) = 0 in Eq. (A3) leads to E exact = Z ∞ f ∗ ( ρ ) (cid:18) − ~ m d dρ + W ( ρ ) (cid:19) f ( ρ ) dρ (A7)with W ( ρ ) = U ( ρ ) − ~ m Q ( ρ ) . (A8)Suppose that for Φ( ρ, Ω) we take the Φ ( ρ, Ω) that gives the lowest adiabatic potential. The corresponding quantities Q ( ρ ) and U ( ρ ) in Eq. (A4) are denoted by Q ( ρ ) and U ( ρ ), respectively. It follows from the Ritz variational principlethat E exact ≤ Z ∞ f ∗ ( ρ ) (cid:18) − ~ m d dρ + W ( ρ ) (cid:19) f ( ρ ) dρ. (A9)If f ( ρ ) is chosen to be the solution of the equation (the adiabatic approximation), (cid:18) − ~ m d dρ + W ( ρ ) (cid:19) f ( ρ ) = E U f ( ρ ) , (A10)with the lowest eigenvalue E U , E U turns out to be an upper bound of E exact : E exact ≤ E U . Differentiating P ( ρ ) = 0with respect to ρ leads to h ∂∂ρ Φ( ρ, Ω) | ∂∂ρ Φ( ρ, Ω) i Ω + Q ( ρ ) = 0 . (A11)4Equation (A7) for E exact is recast to E exact = Z ∞ f ∗ ( ρ ) h − ~ m d dρ + U ( ρ ) + ~ m h ∂∂ρ Φ( ρ, Ω) | ∂∂ρ Φ( ρ, Ω) i Ω i f ( ρ ) dρ. (A12)Since the last term in the square brackets is non-negative, we obtain E exact ≥ Z ∞ f ∗ ( ρ ) (cid:18) − ~ m d dρ + U ( ρ ) (cid:19) f ( ρ ) dρ. (A13)By using the inequality U ( ρ ) ≥ U ( ρ ) and choosing f ( ρ ) to be the solution of the equation (the Born-Oppenheimerapproximation), (cid:18) − ~ m d dρ + U ( ρ ) (cid:19) f ( ρ ) = E L f ( ρ ) , (A14)with the lowest eigenvalue E L , we obtain a lower bound of E exact as E exact ≥ E L .If we calculate the expectation value of H for the wave function Ψ = ρ − ( d − / f ( ρ )Φ ν ( ρ, Ω) with the ν th channelwave function, we confirm Eq. (20) using the same argument as above. Appendix B: Matrix elements with correlated Gaussians
In this appendix we calculate F ρ ( ω ), Eq. (39), using as φ i ( x ) the generating function g ( s ; A, x ) [17, 18, 38, 39] ofthe CG: g ( s ; A, x ) = exp( − ˜ x A x + ˜ sx ) , (B1)where A is an ( N − × ( N −
1) symmetric, positive-definite matrix and s = ( s i ) is an ( N − A − + ˜ s A − s ≈ ρ . (B2)Note that for the special case that A is diagonal, A = ( a i δ i,j ), g ( s ; A, x ) reduces to a product of Gaussian wavepackets: g ( s ; A, x ) = Q N − i =1 exp[ − a i ( x i − s i ) + a i s i ]. We present formulas for F ρ ( ω ) calculated between g ( s ; A, x ) and g ( s ′ ; A ′ , x ). See Ref. [18] for details. The case with s = s ′ = 0 is given in Ref. [40]. Overlap
The function F ρ ( ω ) for O ( x ) = 1 is given by F ρ ( ω ) = (cid:18) (2 π ) N − det B (cid:19) / e − ˜ v B − v (B3)with B = ρ ( A + A ′ ) + 2 iωI, v = ρ ( s + s ′ ) . (B4)Here I is the ( N − × ( N −
1) identity matrix. Since A + A ′ can be diagonalized by an orthogonal matrix, thematrix B can be diagonalized as well. Kinetic energy
To calculate the matrix element for O ( x ) = T Ω = T − T ρ , we use the following relation [40]: ∂∂ρ g ( s ; A, x ) = 1 ρ ( − ˜ x A x + ˜ sx ) g ( s ; A, x ) ,∂ ∂ρ g ( s ; A, x ) = 1 ρ h ( ˜ x A x ) − ˜ x A x + (˜ sx ) −
2( ˜ x A x )˜ sx i g ( s ; A, x ) ,T g ( s ; A, x ) = − ~ m h − A + ˜ ss − s A x + ˜ x A x i g ( s ; A, x ) . (B5)5Combining these results, we obtain T Ω g ( s ; A, x ) = ~ mρ h ρ Tr A − ρ ˜ ss + 2 ρ ˜ s A x + ( d − sx + (˜ sx ) − ρ ˜ x A x − d ˜ x A x −
2( ˜ x A x )˜ sx + ( ˜ x A x ) i g ( s ; A, x ) . (B6)Here d is defined in Eq. (7). The change of variables, x → ρ ξ , yields F ρ ( ω ) as F ρ ( ω ) = ~ mρ Z e − ˜ ξ B ξ +˜ vξ h ρ Tr A − ρ ˜ ss + 2 ρ ˜ s A ξ + ( d − ρ ˜ sξ + ρ (˜ sξ ) − ρ ˜ ξ A ξ − dρ ˜ ξ A ξ − ρ ( ˜ ξ A ξ )˜ sξ + ρ ( ˜ ξ A ξ ) i d ξ . (B7)This integral can be performed analytically. In the case of s = s ′ = 0, we obtain F ρ ( ω ) = ~ mρ (cid:18) (2 π ) N − det B (cid:19) / h ρ A − B − ( ρ A ) − d Tr B − ρ A + 15(Tr B − ρ A ) − M ( B, ρ A ) i . (B8)Here use is made of the formula Z e − ˜ x B x ( ˜ x P x ) d x = (cid:18) (2 π ) N − det B (cid:19) / h B − P ) − M ( B, P ) i (B9)for an ( N − × ( N −
1) symmetric matrix P , where M ( B, P ) is defined by M ( B, P ) = 1det B N − X j>i =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) B B . . . B N − ... ... ... ... P i P i . . . P i N − ... ... ... ... P j P j . . . P j N − ... ... ... ... B N − B N − . . . B N − N − (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) . (B10) Potential energy
The matrix element for O ( x ) = V is conveniently calculated by expressing the distance vector of two nucleons as acombination of Jacobi coordinates r i − r j = ˜ ζ x , (B11)where ζ is an ( N − i and j . For a Gauss potential, O ( x ) = e − τ ( r i − r j ) = e − τ ˜ x ζ ˜ ζ x , F ρ ( ω ) reduces to that of the overlap. For s = s ′ = 0, we obtain F ρ ( ω ) = (2 π ) N − det( B + 2 τ ρ ζ ˜ ζ ) ! / = (2 π ) N − (1 + 2 τ ρ ˜ ζB − ζ ) det B ! / . (B12)In the last step Sherman-Morrison formula, det( B + cζ ˜ ζ ) = (1 + c ˜ ζB − ζ ) det B , is used, where c is a constant.By comparing this result with Eq. (B3) and by noting that B − is constrained by the condition (B2), we expectthat the contribution of the potential energy V to the adiabatic potential behaves as ( ρ/τ ) − for large ρ values. The ρ − dependence was found for a three-body system [25], but our result suggests that it is valid for many-body systemsas well.As an important application of Eq. (B12), we calculate the matrix element for the Coulomb potential, O ( x ) =1 / | r i − r j | . Using 1 | r i − r j | = 2 √ π Z ∞ e − τ ( r i − r j ) dτ, (B13)6and the integral Z ∞ (1 + aτ ) − / dτ = 1 √ a , (B14)we obtain the matrix element for the Coulomb potential as F ρ ( ω ) = r π ρ (cid:18) (2 π ) N − det B (cid:19) / (˜ ζB − ζ ) − / . (B15)As expected, the inverse ρ -dependence appears naturally. [1] See for example, K. Matsuyanagi, M. Matsuo, T. Nakatsukasa, N. Hinohara, and K. Sato, J. Phys. G: Nucl. Part. Phys. , 064018 (2010).[2] H. J. Krappe and K. Pomorski, Theory of Nuclear Fission , Lecture Notes in Physics, , (Springer, Berlin, Heidelberg,2012).[3] K. Hagino and N. Takigawa, Prog. Theor. Phys. , 1061 (2012).[4] P. Ring and P. Schuck,
The Nuclear Many-Body Problem , Texts and Monographs in Physics (Springer, New York, 1980).[5] M. V. Zhukov, B. V. Danilin, D. V. Fedorov, J. M. Bang, I. J. Thompson, and J. S. Vaagen, Phys. Rep. , 151 (1993).[6] C. D. Lin, Phys. Rep. , 1 (1995).[7] R. Krivec, Few-Body Syst. , 199 (1998).[8] E. Nielsen, D. V. Fedorov, A. S. Jensen, and E. Garrido, Phys. Rep. , 373 (2001).[9] J. Macek, J. Phys. B , 831 (1968).[10] N. Barnea, W. Leidemann, and G. Orlandini, Phys. Rev. C , 054001 (2000); Nucl. Phys. A 693 , 565 (2001).[11] N. K. Timofeyuk, Phys. Rev. C , 064306 (2002); ibid. , 054314 (2008); Phys. Rev. A , 032507 (2012).[12] N. P. Mehta, S. T. Rittenhouse, J. P. D’Incao, J. von Stecher, and C. H. Greene, Phys. Rev. Lett. , 153201 (2009).[13] D. Rakshit and D. Blume, Phys. Rev. A , 062513 (2012).[14] K. M. Daily and C. H. Greene, Phys. Rev. A , 012503 (2014).[15] K. Varga, Y. Suzuki, and R. G. Lovas, Nucl. Phys. A 571 , 447 (1994); K. Varga, Y. Suzuki, and I. Tanihata, Phys. Rev.C , 3013 (1995); K. Varga, Y. Suzuki, and R. G. Lovas, Phys. Rev. C , 041302(R) (2002).[16] K. Arai, Y. Ogawa, Y. Suzuki, and K. Varga, Prog. Theor. Phys. Suppl. , 97 (2001).[17] K. Varga and Y. Suzuki, Phys. Rev. C , 2885 (1995).[18] Y. Suzuki and K. Varga, Stochastic Variational Approach to Quantum-Mechanical Few-Body Problems , Lecture Notes inPhysics Vol. m54 (Springer, Berlin, 1998).[19] M. Cavagnero, Phys. Rev. A , 2877 (1986).[20] N. Barnea and A. Novoselsky, Ann. Phys. , 192 (1997).[21] P. Descouvemont, J. Phys. G: Nucl. Part. Phys. , 064010 (2010).[22] U. Fano, Phys. Rev. A , 2402 (1981).[23] A. F. Starace and G. L. Webster, Phys. Rev. A , 1629 (1979).[24] H. T. Coelho and J. E. Hornos, Phys. Rev. A , 6379 (1991).[25] I. J. Thompson, B. V. Danilin, V. D. Efos, J. S. Vaagen, J. M. Bang, and M. V. Zhukov, Phys. Rev. C , 024318 (2000).[26] H. Suno, Y. Suzuki, and P. Descouvemont, Phys. Rev. C , 014004 (2015).[27] O. I. Tolstikhin, S. Watanabe, and M. Matsuzawa, J. Phys. B: At. Mol. Opt. Phys. , L389 (1996).[28] H. Suno, J. Chem. Phys. , 064318 (2011); ibid. , 134312 (2011).[29] M. Brack, J. Damgaard, A. S. Jensen, H. C. Pauli, V. M. Strutinsky, and C. Y. Wong, Rev. Mod. Phys. , 320 (1972).[30] E. Hiyama, Y. Kino, and M. Kamimura, Prog. Part. Nucl. Phys. , 223 (2003).[31] A. I. Baz’ and M. V. Zhukov, Sov. J. Nucl. Phys. , 435 (1970).[32] H. Horiuchi, Prog. Theor. Phys. Suppl. , 90 (1977).[33] Y. Suzuki, R. G. Lovas, K. Yabana, and K. Varga, Structure and Reactions of Light Exotic Nuclei (Taylor & Francos,London, 2003).[34] J. McNamee, F. Stenger, and E. L. Whitney, Math. Comput. , 141 (1971).[35] R. B. Wiringa, S. C. Pieper, J. Carlson, and V. R. Pandharipande, Phys. Rev. C , 014001 (2000).[36] M. Girod and P. Schuck, Phys. Rev. Lett. , 132503 (2013).[37] M. Ito, Phys. Lett. B 636 , 293 (2006).[38] Y. Suzuki, W. Horiuchi, M. Orabi, and K. Arai, Few Body Syst , 33 (2008).[39] S. Aoyama, K. Arai, Y. Suzuki, P. Descouvemont, and D. Baye, Few Body Syst , 97 (2012).[40] J. von Stecher and C. H. Greene, Phys. Rev. A80