Diagonalization scheme for the many-body Schr¨odinger equation 7 1 0 Lorenzo Fortunato∗, Tomohiro Oishi† 2 Dipartimento diFisica eAstronomia “G.Galilei”, n a I.N.F.N., Sezione diPadova, J via Marzolo 8, I-35131 Padova, Italy 7 1 January 18, 2017 ] h t - l Abstract c u A new convenient method to diagonalize the non-relativistic many- n body Schr¨odinger equation with two-body central potentials is derived. [ It combines kinematic rotations (democracy transformations) and exact 1 calculation of overlap integrals between baseswith differentsets of mass- v scaled Jacobi coordinates, thereby allowing for a great simplification of 4 this formidable problem. We validate our method by obtaining a perfect 8 correspondence with the exactly solvable three-body (N = 3) Calogero 6 model in 1D. 4 0 . 1 Outline of the method 1 0 7 Anew methodtodiagonalizethe non-relativisticmany-bodySchr¨odingerequa- 1 tion (1) with two-body central potentials, V is described. The programdevel- ij : v ops through the following main steps: 1) use of mass-scaled Jacobi coordinates i to separate out the center of mass motion; 2) expansion of the wavefunction X in a totally decoupled basis (tensor product of known orthonormal complete r a bases, one for each internal coordinate) ; 3) calculation of matrix elements in this basis: those of T and V are trivially calculated, while, for the calculation 12 of matrix elements of V , we propose the following steps: ij i) kinematicrotations(democracytransformations)toalternativesetsofJa- cobi coordinates in which ij becomes the ’first’ coordinate ii) analytic calculation of displaced overlapintegrals (σ coefficients) ∗[email protected] †[email protected] 1 The last two points allow a great simplification in the treatment of this long- standingquantumproblem,reducingessentiallyallcomputationstoone-dimensional integrals. This method has the advantage of pushing the analytical treatment as far as possible, without introducing too many complications. Points 1) to 3) are part of standard procedures, but 3), that is the core of the new method, contains new ideas. In particular, when the calculations in ii) is obtained ana- lytically, there is an advantage in our method. After this program is achieved in a general theory for N particles in 3D, we will demonstrate its validity in the case of three bosons in 1D, by reproducing in a numerical code the well- knownanalytic results ofCalogero[1, 2]. This opens the way to a largeclass of problems that can be easily studied with our formalism. 1.1 Mass-scaled Jacobi coordinates The three dimensional Schr¨odinger equation for N interacting particles with masses m , i=1,...,N is: i N ~p 2 N i + V ( ~r ) E Ψ(N)(~r ,...,~r )=0, (1) ij ij T 1 N " 2mi | | − # i=1 i<j X X where ~r =~r ~r . The problem of separating the center of mass motion has ij j i − been recently re-analyzed in Ref. [3] for equal masses, while here we are inter- ested in the generalization to different masses. Based on Refs.[4], we introduce N 1 Jacobi coordinates, ξ~ , plus the position of the center of mass : k − ξ~ = µk ~r ~c k µ k+1− k (2) ( ξ~N = q ~c(cid:0)N R~ (cid:1) ≡ with k = 1,...,N 1 and where the vectors~c = Pkmi~ri have been used for − k Pkmi the sake of simplicity. The ξ~ coordinates correspond to the relative position k vectorsconnecting the center of mass ofthe k th cluster (i.e. the cluster made − up of the first k particles, in some, previously adopted, ordering sequence) and the (k+1) th particle. The quantities − ( km )m i k+1 µ = k =1,...,N 1 (3) k k+1m − P i are intermediate reduced-mPassesbetween the cluster of the first k particles and the single (k + 1)-th particle, while µ = N−1 Nm / Nm is a common i i reduced-mass. The momenta ~πk (k = 1,...,Nq−Q1) are Pcanonically conjugate to ξ~ , and ~π P~ is the total momentum, canonically conjugate to the center k N ≡ of mass position. They are related to the lab frame momenta through ∂ξ~ l ~p = ~π . (4) k l ∂~r k Xl (cid:16) (cid:17) 2 The transformationof coordinates (2) leaves invariantthe quadratic forms cor- responding to total moment of inertia and kinetic energy: N N−1 m ~r 2 = µ ξ~ 2+MR~2 (5) i i i i=1 i=1 X X N p~ 2 N−1~π2 P~2 i = i + , (6) 2m 2µ 2M i i=1 i=1 X X whereM = Nm isthetotalmassofthesystem. Theproperty(6)allowsthe j exactseparationofthecenterofmassmotion. ByfactorizingΨ(N)(~r ,...,~r )= 1 N P f(R~)Φ(N)(ξ~ ,...,ξ~ ) with f(R~) = eiK~·R~ and using E = E + E , with 1 N−1 T cm E =~2K~2/2M,one obtains the Schr¨odingerequationin the intrinsic coordi- cm nates: N−1~π2 N i + V (r ) E Φ(N)(ξ~ ,...,ξ~ )=0. (7) ij ij 1 N−1 " 2µ − # i=1 i<j X X Note that we have not changed the potential energy term and its argument as there will be no need for this complication. It is useful to note that ~r = 12 µ /µ ξ~ are collinear proportional vectors. 1 1 p 1.2 Trivial matrix elements Onehasfreedomtochooseamongthemostsuitablebasisforhisscopes: leaving aside the square well, there are essentially four types of potentials that one can choosefrom,eachofwhichcanbeseenasaparticularcaseofSU(1,1)dynamical symmetry. They are illustrated in Fig. 1 with all the dualities that connect them. Totally decoupled basis states can be expressed as a tensor product of the N 1 relative motion wavefunctions that amounts to: − Φ˜(N) = φ1 (ξ~ ) ... φN−1 (ξ~ ) | c i | {ν1} 1 × × {νN−1} N−1 i φ1 (ξ~ )...φN−1 (ξ~ ) , (8) →| {ν1} 1 {νN−1} N−1 i where each ν = n ,ℓ ,m indicates the set of all quantum numbers needed i i i i { } { } tospecifyeachofthesinglemotionwavefunctions(namelytheprincipalq.n.,the angularmomentumq.n. anditsprojectiononthequantizationaxis). Thetensor productwithgoodangularmomentumisexpressedthroughangularmomentum coupling rules into sums of simple products, that can more practically be used as basis elements and one can forget the Clebsch-Gordan coefficients because these are reabsorbed into the amplitudes in the process of diagonalization. For ease of notation, in the following we will indicate each basis element with a single index, c, that counts the basis states. Within these basis states the matrix elements of T are trivially calculated, i because~π2 actsonlyonthei thcoordinate,whilealltheremainingcoordinates i − 3 DOUBLEDUALITYSCHEME “Spherical” ∞ Harm.Osc. Coulomb 0 → → mp. SU(1,1) mp. Asy Davidson Kratzer Asy “Deformed” Figure 1: Dualities between analytically solvable harmonic-type potentials (Harmonic Oscillator and Davidson potentials) and coulomb-type potentials (Coulomb and Kratzer potentials). The SO(2,1) SU(1,1) group structure is ∼ common to all cases. The two right quadrants correspond to potentials that asymptotically go to zero as 1/r, while the left ones go to infinity as r2 (there- forethey do notadmitunboundeigenstates). The twolowerquadrantsidentify potentialswith the minimum insomenon-nullfinite point,while the twoupper ones identify “spherical” cases (minimum in r = 0 for the harmonic oscillator 0 caseandsingularityfortheCoulombcase). Simplearrowsindicatehomologyin the algebraic treatment (same commutation relations, different representations in terms of differential operators),while double arrowsindicate an extension of the commutation relations that allows a formal analogy between the deformed and spherical cases (namely Zˆ Zˆ +k/Zˆ , with notation of Refs. [5]). 1 1 2 → lead to Kronecker delta’s in the respective quantum numbers: hΦ˜(cN′ ) |Ti |Φ˜(cN)i=δ{ν1′}{ν1}...δ{νi′−1}{νi−1} hφiνi′(ξ~i)|Ti |φiνi(ξ~i)iδ{νi′+1}{νi+1}...δ{νN′ −1}{νN−1}, (9) and the same is true for the matrix elements of V because its argument is 12 already one of the Jacobi coordinates (the first): hΦ˜(cN′ ) |V12 |Φ˜(cN)i=hφ1ν1′(ξ~1)|V12( µ1/µξ~1)|φ1ν1(ξ~1)i δ{ν2′}{ν2p}...δ{νN′ −1}{νN−1} . (10) 1.3 Democracy transformation and non-trivial matrix el- ements Now we turn to the core part of the paper, the determination of the matrix elements of the two-body potentials between each pair of particles, V ( ~r ). ij ij | | Clearly the choice of Jacobi coordinates operated above has the virtue of sep- arating the center of mass motion, but admittedly complicates the calcula- tion of the matrix elements we are after, except V , because, by construction 12 4 ξ~ µ/µ ~r . Noticethatthe’firstcoordinatespair’,namely12,hasthisspe- 1 1 12 ≡ cial property, while all other pairs ij don’t, therefore the calculation of matrix p elements of V ( ~r ~r ) is not separable. This is a long-standing problem of ij i j | − | whichwe offer a solutionby merginggrouptheoreticaland analytictechniques. i) First, we have to recall the concept of kinematic rotations, also called democracy transformations [6, 7, 8]. These are linear transformations between different choices (i.e. orderings) of Jacobi coordinates, viz. N−1 ξ~ = D ξ~ , (11) i ij j j=1 X where underlined denotes transformed coordinates. The set of matrices D with matrix elements D form the democracy group that includes rotations and ij reflections. As these transformations are norm-conserving, the group can be identified with the orthogonalgroup O(N-1). In the Jacobi coordinates system defined above, the calculation of matrix elements of V ( ~r ) involves exactly ij ij | | the set of the first j coordinates, namely ξ~ ,...,ξ~ . We call these hot coordi- 1 j { } nates and the others cold coordinates. For example,in the case of four particles depicted in Fig. 2, the calculation of the matrix element of V involves the 13 6-dimensional integral φ(ξ~ )φ(ξ~ ) V ( ξ~ ξ~ /2 ) φ(ξ~ )φ(ξ~ ) , while the 1 2 13 2 1 1 2 h | | − | | i third coordinate can be separated out. In order to achieve a full separation, it is necessary to perform a kinematic rotationin the subspace of hot coordinates that brings one of them to ~r . This can be viewed as a block-matrix of the ij form: ξ~1 ξ~1 ... Dj×j O ... ξ~ξ~Njξ~+...j−11 = O 1N−1−j ξ~ξ~Njξ~+...j−11 (12) where and are the identity and null (sub-)matrices of appropriate dimen- 1 O sions andthe hotcoordinates havebeen highlighted. We canalwayschoose the new set of hot coordinates in such a way that one of them coincides with ~r ij andthe othersconnect(sub-)clustersofincreasingsizeasusual(since theorder is immaterialto the finalresult,wechoosethe firsthotcoordinate,i.e. the 1-th as~r coordinate). Thekinematicrotationproposedhereallowstofactorizethe ij matrix element of V into the matrix element acting only on bra’s and ket’s of ij a single relevant coordinate ξ~ and a product of (j 1) Kronecker delta’s, thus 1 − simplifying noticeably the problem. In Fig. 2 we give an example for the case of 4 particles. The initial co- ordinates ξ~, that are suitable for calculating the matrix elements of the V i 12 interaction, are rotated into ξ~, that are suitable for calculating the matrix ele- i mentsofV . Noticethatthethirdcoordinate(thatisacoldcoordinatehere)is 13 5 4 4 2 2 ξ~3 D11 D12 0 ξ~3 D21 D22 0 ξ~2 0 0 1 ξ~ 1 1 ξ~2 3 1 ξ~ 3 1 Figure2: Initialcoordinatesystem(left),suitableforcalculatingtheinteraction between 1 and 2, rotated coordinate system (right), suitable for calculating the interactionbetween1and3andthematrixofkinematicrotationthattransforms the coordinates ξ~ into ξ~. i i notrotated. IfwecallU theN-dimensionalnunitarytransformation(2)that ord brings the initial coordinates of N particles to a certain set of mass-scaled Ja- cobi coordinates with a given of the labels, then the democracy transformation Dbetweenanytwosuchorderingisgivenbythe(N 1)-dimensionalsubmatrix − of the product Uord(Uord′)−1. For the sake of simplicity, we will be concerned only with orderings that brings the relevant pair ij to the ’first’ position. ii) Now that we have specified how to change coordinates, before proceed- ing to the calculation of matrix elements, we have to find the relation between the original and transformed basis states. This crucial step is not easily ac- complishedingeneral,butfor harmonicoscillatorbasesone canderive analytic formulas. The idea is that we need to expand each basis substate in terms of the basis statesofthe transformedcoordinatesasfollows(keeping inmind that the coordinates j+1, ,N 1 are not rotated): ··· − φ1 (ξ~ ) ... φN−1 (ξ~ ) = | {ν1} 1 × × {νN−1} N−1 i = σc φ1 (ξ~ ) ... φN−1 (ξ~ ) , (13) c | {ν1} 1 × × {νN−1} N−1 i X{c} where σc are expansion coefficients labeled by the set of all quantum numbers c needed in the transformed basis, i.e. c = ν ,...,ν . From Eq. (13), one 1 N−1 { } gets σc = φ1 (ξ~ ) ... φN−1 (ξ~ ) φ1 (ξ~ ) ... φN−1 (ξ~ ) , (14) c h {ν1} 1 × × {νN−1} N−1 | {ν1} 1 × × {νN−1} N−1 i that, in principle, cannot factorize into a product of N 1 overlaps, because each rotated ξ~ depends on all the non-rotated ξ~ , accord−ing to Eq. (11). This k k means that the integral would be an almost inextricable 3j-dimensional tangle for large values of j. Fortunately, at least in the case of a set of (N 1) − decoupledharmonicoscillatorbases,wecansolvethisproblemanalytically. The calculation of the σ coefficients can be recast in terms of summations of certain coefficients, that come from the repeated application of the umbral identities 6 for Hermite polynomials and from the series definition of the latter, and of Gamma functions, that come from the analytic integration of the remaining rotatedpolynomials. Thelengthyderivationandgeneralizationwillbedescribed elsewhere. Of course the numerical evaluation of these analytic expressions is hundreds of times quicker than the brute-force calculation of multidimensional integrals and this is a considerable advantage. Using expansion(13), the matrix elements that we are seeking are therefore given by hΦ˜(cN′ ) |Vij(|~rij |)|Φ˜(cN)i=δ{νj′+1}{νj+1}...δ{νN′ −1}{νN−1} σcc′′∗ σcchφ1ν1′(ξ~1)|Vij(ξ1)|φν11(ξ~1)iδ{ν2′}{ν2}...δ{νj′}{νj} , (15) c′,c X where some care must be taken with the labeling of indexes of ν (unprimed meansket, primedmeans bra,underlinedmeanstransformed). The lastmatrix element can be reduced to the one-dimensional integral in the radial variable φi (ξ ) V φi (ξ ) forcentralpotentials. Notice,infact,thatthecoefficients h νi′ i | ij | νi i i σ ofEq.(15)differ onlyinthe lastsetofquantumnumbers(associatedwiththe 1 coordinate) that can be further reduced to just the two principal quantum numbers n′ and n , respectively. 1 1 Eqs. (9),(10)and(15)giveanewmethodtodiagonalizetheN-bodySchro¨din- ger equation with two-body interactions. This method proceeds through exact analytical steps, albeit in practice one must decide a truncation in the number of quanta for each φ in Eq.(8). Numerical approximations come into play only in the calculation of one dimensional integrals and in the matrix diagonaliza- tion. Several exact formulas exist for certain analytic functions of the distance (take for example integrals of powers in the h.o.), therefore the calculation of matrixelementscanoftenbedoneexactly. Inmanycases,afterdiagonalization apropersymmetrizationproceduremustbeadopted. Thismethodtransfersthe (often untreatable) complexity of the calculationof matrix elements of the mu- tualparticle-particleinteractiontoanumberofsimpleoverlapsandwell-known matrix elements ofa single coordinate,thereby allowingastraightforwardsolu- tion of the formidable many-body problem. 2 Validation in the 1D three-body case Before embarking on a longer campaign of theoretical studies involving more sophisticated calculations, we need demonstrate the validity of our method on analyticcases. Thereareonlyahandfulsolvedmodels[1,2,9],amongwhichwe choosethe well-knownonedimensionalCalogeromodelwithharmonicpairwise interactions of the type V = ~ω(x x )2. The exact spectrum for N=3 ij i j − particles with Bose statistics on a line is found in Ref. [2], namely E = Cal √3~ω(2n+l+1) with n,l non-negative integers, and reproducing this result it’s easy with our method. Essentially, we just use the matrix elements of x2 7 deg. E (1,1) (0,3) 2 (1,0) (0,2) 2 (0,1) 1 (0,0) (n,l) 1 Figure 3: Spectrum of the Calogero linear model (1D) with harmonic pairwise interactions. Theanalyticresultintermsofquantumnumbers(n,l)isobtained easily with our numerical method. Degeneration of energy levels is shown on the right. between harmonic oscillator basis states and use the exact formula for the σ coefficientsinthe calculationsofV andV . We givehereonlythe finalresult 13 23 for three particles in 1D, where the sets of quantum numbers simply reduces to the number of oscillator quanta: σc = 2−(n1+n2+n1+n2)/2 n1 n1 n2 n2 (n ,n ,k,j,D) (n ,n ,k,j,D) c π n !n !n !n ! k j I 1 1 I 2 2 1 2 1 2 k=0(cid:18) (cid:19)j=0(cid:18) (cid:19) X X (16) p where aresomeone-dimensionalintegralsthatcanbewrittenintermsofsum- mationIs of Gamma functions and D is the democracy transformation between the set ’starting’ with 12 and the one with either 13 or 23 as ’first coordinate’. In a basis truncated at N = 15 quanta, we get the linear spectrum of q Calogeroeffortlesslyinthe turnofafewsecondsonatable-topmachine,witha numerical precision that is sufficient to all practical purposes and with correct degeneration of energy levels. The lowest eigenvalues are shown in Fig. 3. In our opinion, the mathematical derivation of the new method and the proof-of-principle discussed in the present paper hold great promises and pave the way to studies of greater impact, because we have now a handy theoretical tool that allows to explore several physics systems: from academic cases such as the spectrum and wavefunctions of N particles in 1D, to realistic models of few-bodysystems (atoms,molecules,nuclei,BEC)in3Dandto moreadvanced ideas such as for instance Efimov states, condensation, exotic systems, phase transitions in stable and unstable systems, etc. Acknowledgements In:Theory, PRAT Project n. CPDA154713,Univ. Padova (Italy). L.F. acknowledgesfruitful discussions with A.Richter and V.Efros (in 2011) andwithA.Vitturi(2011-2016). Mostofthepresentmaterialhasbeenconceived while working at the ECT*(Trento) in 2011, but could only be finalized now. 8 The computing facilities offered by CloudVeneto (CSIA Padova and INFN) are acknowledged. References [1] F. Calogero,J. Math. Phys., vol. 10, pp. 2191-2196(1969). [2] F. Calogero,J. Math. Phys., vol. 12, pp. 419-436(1971). [3] L. Fortunato, J.Phys.A:Math. Gen. 43, 065301(2010). [4] L.M. Delves, Nucl.Phys. 20, 275-308 (1960); V. Aquilanti, A. Beddoni, A.LombardoandR.Littlejohn,Int.J.QuantumChem.,Vol89,277-291 (2002). [5] D.J. Rowe, Prog.Part.Nucl.Phys. 37, 265-348 (1996); D.J .Rowe and C. Bahri, J.Phys.G A10, 4947 (1998); L. Fortunato and A. Vitturi, J.Phys.G:Nucl.Part.Phys. 29 1341-1349(2003). [6] J.D. Louck and H.W. Galbraith, Rev.Mod.Phys. 44, 540-601 (1972). [7] N. Barnea and A. Novoselsky,Ann.Phys. 256, 192-225 (1997). [8] U.Fano,D.Green,J.L.BohnandT.A.Heim,J.Phys.B:At.Mol.Opt.Phys. 32, R1-R37, (1999). [9] B. Sutherland, Beautiful Models: 70 years of exactly solved quantum many-body problems, World Scientific Co. Pte. Ltd., Singapore (2004) 9