Semiclassical trace formula for the two-dimensional radial power-law potentials
aa r X i v : . [ n li n . S I] J un SEMICLASSICAL TRACE FORMULA FOR THE TWO-DIMENSIONAL RADIALPOWER-LAW POTENTIALS
A. G. Magner
Institute for Nuclear Research, 03680 Kiev, Ukraine
A. A. Vlasenko
Institute for Nuclear Research, 03680 Kiev, Ukraine andInstitute of Physics and Technology, NTUU “KPI”, 03056 Kyiv, Ukraine
K. Arita
Department of Physics, Nagoya Institute of Technology, Nagoya 466-8555, Japan (Dated: October 1, 2018)The trace formula for the density of single-particle levels in the two-dimensional radial power-lawpotentials, which nicely approximate up to a constant shift the radial dependence of the Woods-Saxon potential and its quantum spectra in a bound region, was derived by the improved stationaryphase method. The specific analytical results are obtained for the powers α = 4 and 6. Theenhancement of periodic-orbit contribution to the level density near the bifurcations are found tobe significant for the description of the fine shell structure. The semiclassical trace formulas forthe shell corrections to the level density and the energy of many-fermion systems reproduce thequantum results with good accuracy through all the bifurcation (symmetry breaking) catastrophepoints, where the standard stationary-phase method breaks down. Various limits (including theharmonic oscillator and the spherical billiard) are obtained from the same analytical trace formula. PACS numbers: 05.45.-a,05.45.Mt,21.60.Cs
I. INTRODUCTION
According to the shell-correction method (SCM) [1, 2],the oscillating part of the total energy of a finite fermionsystem, the so-called shell-correction energy δU , is associ-ated with an inhomogeneity of the single-particle energylevel distributions near the Fermi surface. Depending onthe level density at the Fermi energy – and thus the shell-correction energy δU – being a maximum or a minimum,the many-fermion system is particularly unstable or sta-ble, respectively. Therefore, the stability of this systemvaries strongly with particle numbers and parameters ofthe mean-field potential and external force.A semiclassical periodic orbit theory (POT) of shelleffects[3–6] was used for a deeper understanding, basedon classical pictures, of the origin of nuclear shell struc-ture and its relation to a possible chaotic nature of thedynamics of nucleons. This theory provides us with a nicetool for answering, sometimes even analytically, the fun-damental questions concerning the exotic physical phe-nomena in many-fermion systems; for instance, the originof the double-humped fission barrier and, in particular, ofthe creation of the isomer minimum in the potential en-ergy surface [7–11]. Some applications of the POT to nu-clear deformation energies were presented and discussedfor the infinitely deep potential wells with sharp edges inrelation to the bifurcations of periodic orbits (POs) withthe pronounced shell effects.In the way to more realistic semiclassical calculations,it is important to account for a diffuseness of the nuclearedge. It is known that the central part of the realisticeffective mean-filed potential for nuclei or metallic clus- ters are described by the Woods-Saxon (WS) potential V WS ( r ) [12]. The idea of Refs. [13, 14] is that the WS po-tential is nicely approximated (up to a constant shift) bymuch a simpler power-law potential which is proportionalto a power of the radial coordinate r α . The approximateequality V WS ( r ) ≈ V WS (0) + W r α (1.1)holds up to around the Fermi energy with a suitablechoice of the parameters W and α . In the case of thespatial dimension D = 2, one can use Eq. (1.1) for arealistic potential of electrons in a circular quantum dot[10, 11, 15]. We shall derive first the generic trace formulafor this radial power-law (RPL) potential in the case oftwo dimensions, and then discuss its well known limitsto the harmonic oscillator and cavity (billiard) poten-tials [11]. The main focus will be aimed to the non-lineardynamics depending on the power parameter α to showthe symmetry-breaking (bifurcation) phenomena. Theylead to the remarkable enhancement of PO amplitudes ofthe level density and energy shell corrections which wasfound within the improved stationary phase approxima-tion (improved SPM, or simply ISPM) [9, 10, 16, 17]. TheISPM means more exact evaluation of the trace formulaintegrals with the finite limits over a classically acces-sible phase-space volume and with higher-order (if nec-essary) expansions of the action phase of the exponentand pre-exponent factors up to the first non-zero termswith respect to the standard SPM (SSPM) [3–6]. In thisway, one may remove the SSPM discontinuities and di-vergences.The manuscript is organized as follows. In Sec. II theclassical dynamics is specified for the RPL potentials.The trace formulas for the RPL potentials in two dimen-sions are derived in Sec. III. Section IV is devoted tothe comparison of the semiclassical calculations for theoscillating level density and shell-correction energy withquantum results. The paper is summarized in Sec. V.Some details of our POT calculations, in particular fullanalytical derivations at the powers α = 4 (see alsoRef. [18]) and 6 for all POs and those at arbitrary α for the diameter and circle orbits, are given in the Ap-pendixes A–E. II. CLASSICAL DYNAMICS ANDBIFURCATIONS
The radial power-law (RPL) potential model is de-scribed by the Hamiltonian H = p m + E (cid:18) rR (cid:19) α , (2.1)where m is the mass of the particle; R and E are intro-duced as constants having the dimension of length andenergy, respectively, and are related with W in Eq. (1.1)by W = E /R α . (In practice, we fix E and adjust theWS potential by varying R and α .) This Hamiltonianincludes the limits of the harmonic oscillator ( α = 2)and the cavity ( α → ∞ ); realistic nuclear potentials withsteep but smooth surfaces correspond to values in therange 2 < α < ∞ . The advantage of this potential isthat it is a homogeneous function of the coordinates, sothat the classical equations of motion are invariant underthe scale transformations: r → s /α r , p → s / p , t → s /α − / t with E → sE . (2.2)Therefore, one only has to solve the classical dynamicsonce at a fixed energy, e.g., E = E ( s = 1); the resultsfor all other energies E are then simply given by thescale transformations (2.2) with s = E/E by definitionin the last equation of Eq. (2.2). This highly simplifiesthe POT analysis [13, 19]. Note that the definition (2.1)can also be generalized to include deformations (see, e.g.,Ref. [14, 19]).As we consider the spherical RPL Hamiltonian (2.1),it can be written explicitly in the two-dimensional (2D)spherical canonical phase-space variables { r, ϕ ; p r , p ϕ } ,where ϕ is the azimuthal angle (a cyclic variable), p ϕ = L is the angular momentum, and the radial momentum p r is given by p r ( r, L ) = r p ( r ) − L r ,p ( r ) = s m (cid:20) E − E (cid:18) rR (cid:19) α (cid:21) . (2.3) The classical trajectory (CT) r ( t ) can be easily found byintegrating the radial equation of motion ˙ r = p r /m withEq. (2.3). Transforming the spherical canonical variablesinto the action-angle ones, for the actions I r , I ϕ one has I r = 1 π Z r max r min p r d r ≡ I r ( E, L ) , (2.4) I ϕ = 12 π Z π p ϕ d ϕ ≡ L , (2.5)where r min and r max are the turning points which are thetwo real (positive) solutions of the equation p r ( r, L ) = 0.The definition (2.1) can be used in arbitrary spatialdimensions, as long as r is the corresponding radial vari-able. In practice, we are interested only in the 2D and3D cases. The spherical 3D and the circular 2D poten-tial models have common PO sets, see Fig. 1. For α > K = 1 (3) in the 2D(3D) cases] are specified by three integers and labeled as M ( n r , n ϕ ), where n r and n ϕ are mutually commensu-rable numbers of oscillations in the radial direction, andof rotations around the origin, each for a primitive or-bit, respectively; and M is the repetition number. Forthe isotropic harmonic oscillator ( α = 2), all the classicalorbits are periodic ones with (degenerate) ellipse shapes.By slightly varying α away from 2, the specific diame-ter and circle orbits appear separately, and they remainas the shortest POs with the corresponding degeneracies K = 1 and 0. With increasing α , the circle orbit and itsrepetitions cause successive bifurcations generating vari-ous new periodic orbits { n r , n ϕ } , n r > n ϕ . Fig. 1 showssome of the shortest POs M ( n r , n ϕ ). The shortest PO isthe diameter which has the degeneracy K = 1 in the 2Dproblem at α >
2. Other polygon-like orbits have K = 1at α > α bif , where α bif is a bifurcation value (see itsspecific expression below). The circle orbit having maxi-mum angular momentum is isolated ( K = 0) for the 2Dsystem (except for the bifurcation points).For the frequencies of the radial and angular motion ofparticle, one finds ω r = ∂H∂I r = (cid:18) ∂I r ∂E (cid:19) − L , ω ϕ = ∂H∂L = − ( ∂I r /∂L ) E ( ∂I r /∂E ) L , (2.6)where I r = I r ( E, L ) [Eq. (2.4)] is identical to the energysurface H ( I r , L ) = E . Thus, the PO condition is writtenas f ( L ) ≡ ω ϕ ω r = n ϕ n r , (2.7)where f ( L ) = − (cid:18) ∂I r ( E, L ) ∂L (cid:19) E = Lπ Z r max r min d rr p r ( r, L ) . (2.8)The energy surface I r = I r ( E, L ) is simply considered asa function of only one variable L [Eq. (2.4)]. The solu-tions to the PO equation [see Eq. (2.7)], L ∗ = L ∗ ( n r , n ϕ ), α τ (3,1)(4,1) (5,2) 2(3,1)(7,2)2(4,1) (7,3) (8,3) 3(3,1)(2,1) 2(2,1) 3(2,1)C 2C 3C circlediameterpolygon FIG. 1. Scaled periods τ PO of some short POs as functions of the power parameter α in dimensionless units m = R = E = 1(Appendix B). Thin solid curves are the circle orbits MC , dashed curves are the diameters M (2 , M ( n r , n ϕ ) ( n r > n ϕ ); their bifurcations from the MC are indicated by open circles. for the given co-primitive integers n ϕ and n r define theone-parametric families K = 1 of orbits M ( n r , n ϕ ) be-cause L is the single-valued integral of motion, whichis only one (besides the energy E ) in the 2D case[4, 7].The azimuthal angle ϕ can be taken, for instance, as aparameter of the orbit of such a family.According to the limit f ( L ) → / L →
0, one hasthe diameter orbits M (2 ,
1) as the specific one-parametric( K = 1) families related to the solution L = 0 of Eq. (2.7).The other specific solutions are the isolated ( K = 0) cir-cle orbits M C by which we represent the M -th repeti-tion of the primitive circle orbit C . The radius r C ofthe circle orbit is determined by the system of equations r min = r max ≡ r C , or equivalently by equations (A1) (seeAppendix A). Thus the angular momentum of the circleorbit is given by L C = r C p ( r C ). As seen obviously fromthe condition of the real radial momentum p r [Eq. (2.3)],this L C is the maximal value of the angular momentum L , i.e., 0 ≤ | L | ≤ L C .As shown in Appendix A, for the stability factor F MC of the circle orbit M C in the radial direction, defined inRefs. [3, 11] through the trace of the PO stability matrix,Tr( M MC ), one obtains F MC = 2 − Tr( M C ) M = 4 sin (cid:20) πM Ω C ω C (cid:21) = 4 sin (cid:2) πM √ α (cid:3) , (2.9)where Ω C [Eq. (A7)] and ω C [Eq. (A3)] are the radial andangular frequencies of the circle orbit. This factor F MC iszero at the bifurcation points α bif by the definition of thestability matrix, Tr( M C ) M = 2, for the POs ( Ω C /ω C ≡√ α = n r /n ϕ ), α bif = n r n ϕ − . (2.10) The PO family M ( n r , n ϕ ), which corresponds to the so-lutions L ∗ < L C of the PO equation (2.7), exists forall α > α bif . There is the specific bifurcation point α = 2 in the spherical harmonic oscillator (HO) limitwith the frequency ω ϕ = p E / ( mR ), where one hasthe two-parametric families at any L within a continuum0 ≤ L ≤ E/ω ϕ . In the HO limit, the above specifiedcircle and diameter orbits belong to these families. Inthe circular billiard limit α → ∞ , the isolated circle or-bit ( K = 0) is degenerating into the billiard boundary r C → R , L C → √ mE R [see the limit α → ∞ inEq. (A2) for r C ].Another key quantity in the POT is the curvature K of the energy surface I r = I r ( E, L ) given by K = ∂ I r ( E, L ) ∂L = − ∂f ( L ) ∂L , (2.11)where f ( L ) is the ratio of frequencies [Eq. (2.8)]. Asshown below, the curvature (2.11) and Gutzwiller factor(2.9) are the key quantities for calculations of the magni-tude of the PO contributions into the semiclassical leveldensity. III. TRACE FORMULAS
The level density g ( E ) for the Hamiltonian H ( r , p ) canbe obtained by using the phase-space trace formula (in D dimensions) [9, 16, 17, 20]: g scl ( E ) = 1(2 π ~ ) D Re X CT Z d r ′ Z d p ′′ δ ( E − H ( r ′′ , p ′′ )) × |J CT ( p ′⊥ , p ′′⊥ ) | / exp (cid:18) i ~ Φ CT − i π µ CT (cid:19) . (3.1)The sum is taken over all discrete CT manifolds for aparticle moving between the initial r ′ , p ′ ; and the final r ′′ , p ′′ points with a given energy E . Any CT can beuniquely specified by fixing, for instance, the initial con-dition r ′ , and the final momentum p ′′ for a given time t CT of the motion along the CT. For the action phase Φ CT inexponent of (3.1), one has Φ CT ≡ S CT ( p ′ , p ′′ , t CT ) + ( p ′′ − p ′ ) · r ′ = S CT ( r ′ , r ′′ , E ) − p ′′ · ( r ′′ − r ′ ) , (3.2)where S CT ( p ′ , p ′′ , t CT ) = − R p ′′ p ′ d p · r ( p ) and S CT ( r ′ , r ′′ , E ) = R r ′′ r ′ d r · p ( r ) are the actions in the mo-mentum and coordinate representations, respectively. InEq. (3.1), J CT ( p ′⊥ , p ′′⊥ ) is the Jacobian for the transfor-mation of the initial momentum p ′⊥ to the final one p ′′⊥ in the direction perpendicular to CT. µ CT is the Maslovphase related to the number of conjugate (turning andcaustics) points along the CT [21, 22].One of the terms in Eq. (3.1) is related to the localshort zero-action CT which is the well known Thomas-Fermi (TF) level density [10, 17]. For calculations of theother oscillating terms of the trace integral (3.1), onemay use the ISPM, expanding the action phase Φ CT andpre-exponent factor in both p ′′ and r ′ variables up to thefirst non-zero terms with the finite integration limits overthe classically accessible phase-space region [10, 17]. The stationary phase conditions are equivalent to the periodic-orbit equations, and therefore, the oscillating level den-sity can be presented as the sum over POs in a potentialwell [10, 11]. A. One-parametric orbit families ( K = 1) In order to obtain the contribution of the one-parametric families of the maximal degeneracy K =1 into the phase-space trace formula (3.1), it is use-ful to transform the usual Cartesian phase-space vari-ables { p ; r } to the other canonical action-angle ones { I ; Θ } , specified in the spherical action-angle variablesas Θ = { Θ r , Θ ϕ ≡ ϕ } ; I = { I r , I ϕ ≡ L } . The Hamil-tonian H , action phase Φ CT , and other related quanti-ties of the integrand in Eq. (3.1) [e.g. H = H ( I ) = H ( I r , I ϕ ) ≡ H ( I r , L )] are independent of the angle vari-ables Θ . Therefore, one can easily perform the integra-tion over these angle variables Θ , which gives the factor(2 π ) . Then, taking the integral over I r exactly by us-ing the energy conserving δ -function, for the oscillatingterms of the CT sum (3.1), one obtains δg scl ( E ) = 12 ~ Re X M,n r ,n ϕ Z d L ω r × exp (cid:26) πi ~ M [ n r I r ( E, L ) + n ϕ L ] − iπ µ M,nr,nϕ (cid:27) . (3.3) Here, the phase (3.2) is expressed in terms of the cor-responding action-angle variables through the actions inthe considered mixed representation, Φ CT = 2 πM [ n r I r ( E, L ) + n ϕ L ] , (3.4) n r and n ϕ are positive co-primitive integers, M is anonzero integer, ω r is the radial frequency in Eq. (2.6).We also omit the upper indexes in I (or { I r , L } ) variableswhich represent initial (prime) and final (double primes)values of Eq. (3.1), taking explicitly into account thatthese variables are constants of motion for the spher-ical integrable Hamiltonian. The integration limits inEq. (3.3) for L are − L C ≤ L ≤ L C , where L C is themaximum value corresponding to the circle orbit. Allquantities in the integrand are taken at the energy sur-face I r = I r ( E, L ) [Eq. (2.4)]. Thus, Eq. (3.3) is similarto the oscillating component of the semiclassical Pois-son summation trace formula which can be obtained di-rectly by using the EBK quantization rules [5, 11] forthe spherically symmetric Hamiltonian. Note that, be-fore taking the trace integral over the angular momen-tum L by the SPM in Eq. (3.3), one can formally con-sider positive and negative M , as those related to thetwo opposite directions of motion along a CT (with dif-ferent signs of the angular momentum). They give, ofcourse, equivalent contributions into the trace formula,due to a time-reversal symmetry of the Hamiltonian, andtherefore, one can write simply the additional factor 2 inEq. (3.3) but with a further summation over only posi-tive integers M . It is in contrast to the standard Poissonsummation trace formula [11] (except for its TF compo-nent) because there is no zero values of the integers inEq. (3.3), n ϕ /n r >
0. The essential point in the deriva-tions of Eq. (3.3) from Eq. (3.1) is that the generatingfunction Φ CT [Eq. (3.4)] is independent of the angle vari-ables for families of the maximal degeneracy K = 1 inthe integrable Hamiltonian. Notice that in these deriva-tions, the SPM conditions were satisfied simultaneouslywithin the continuum of the stationary points 0 ≤ ϕ , Θ r ≤ π , which form CTs, but they are not yet POsgenerally speaking for arbitrary angular momentum L .(Exceptions are the cases of the complete degeneracy asthe spherical HO; see below.) The integration range inEq. (3.3) taken from the minimum, L − = 0, to the max-imum, L + , value (for anticlockwise motion, for instance)covers the contributions of a whole manifold of closed andunclosed CTs of the tori in the phase space at the energysurface around the stationary point, L = L ∗ , which cor-responds to the PO [17]. We shall specify the integrationlimits L + for the contribution of the ( K = 1) diameterfamilies M ( n r = 2 , n ϕ = 1) into Eq. (3.3) in Appendix D.Then we apply the stationary phase condition withrespect to the variable L for the exponent phase Φ CT [Eq. (3.4)] in the integrand of Eq. (3.3),( ∂Φ CT /∂L ) ∗ = 0 , (3.5)which is equivalent to the resonance condition (2.7).This condition determines the stationary phase point, L = L ∗ = L PO , related to the families of the POs M ( n r , n ϕ ). All these roots of equation (2.7) for K = 1families M ( n r , n ϕ ) are in between the minimum value L = L ∗ = 0 for diameters, and a maximum one L = L C ,0 ≤ L PO ≤ L C (anticlockwise motion, for example). Ex-panding now the exponent phase Φ CT [Eq. (3.4)] in thevariable L up to the second order, and assuming thatthere is no singularities in the curvature (2.11) for thecontribution of all K = 1 families, one has Φ CT = S PO ( E ) + 12 J ( L )PO ( L − L ∗ ) + · · · , (3.6)where S PO ( E ) is the action along one of the isolated POfamilies determined by Eq. (2.7), S PO ( E ) = 2 πM [ n r I r ( E, L ∗ ) + n ϕ L ∗ ] . (3.7)In this equation, M is the number of repetitions of theprimitive ( M = 1) orbit, I r ( E, L ) is the energy surface[Eq. (2.4)], L = L ∗ ( n r , n ϕ ) is the solution of the POequations (2.7) or (3.5). The Jacobian J ( L )PO in Eq. (3.6)measures the stability of the PO with respect to the vari-ation of the angular momentum L at the energy surface, J ( L )PO = (cid:18) ∂ S CT ∂L (cid:19) L = L ∗ = 2 πM n r K PO , (3.8) K PO = (cid:18) ∂ I r ∂L (cid:19) L = L PO , (3.9)where K PO is the curvature (2.11), (B5) of the energysurface I r = I r ( E, L ) at L = L ∗ = L PO .For the sake of simplicity, we shall discuss the simplestleading ISPM taking up to the second order term in theexpansion over ( L − L ∗ ) for the action phase [Eq. (3.6)],and accounting for only the zeroth order component forthe pre-exponential factor in Eq. (3.3). Substituting nowthese expansions into Eq. (3.3), one can take the pre-exponential factor off the integral at L = L ∗ . Thus,applying Eq. (3.6), we are left with the integral over L of a Gaussian type integrand within the finite limitsmentioned above for contributions of the one-parametricpolygon-like and diameter families, including the contri-bution of boundaries for 0 < n ϕ /n r ≤ /
2. Taking thisintegral over L within the finite limits, one obtains theISPM trace formula, δg ( K ) ( E ), for contributions of theone-parametric ( K = 1) orbits, δg (1) ( E ) = Re X PO A (1)PO ( E ) × exp (cid:20) i ~ S PO ( E ) − i π σ PO − iφ d (cid:21) . (3.10)The sum is taken over the discrete families of the PO M ( n r , n ϕ ) with n r ≥ n ϕ , M ≥ S PO ( E ) is the action(3.7) along these POs. For the amplitudes A (1)PO , one finds A (1)PO = T PO π ~ / p M n r K PO erf( Z − PO , Z +PO ) , (3.11) just as for K = 1 families in the elliptic billiard [16],and the integrable H´enon-Heiles (IHH) potentials [17],with the period T PO = 2 πn r /ω r = 2 πn ϕ /ω ϕ along theprimitive ( n r , n ϕ ) PO. In the RPL Hamiltonian underconsideration, one has T PO = d S P O ( E )d E = π ( α + 2) αE [ n r I r ( E, L
P O ) + n ϕ L P O ](3.12)[see Eqs. (3.7) for the action S P O and (B2) with (B1) forthe scaling transformations]. In Eq. (3.11), K PO is thecurvature of the energy surface I r = I r ( E, L ) [ K PO > α > u, v ) = erf( v ) − erf( u ) with the standard errorfunctions, erf( z ), of the complex arguments z . Thesearguments are specified by Z ± PO = p − iπM n r K PO / ~ ( L ± − L PO ) , (3.13) L − = 0 and L + = L C for all K = 1 polygon-like PO fam-ilies (besides of the diameters, see below). For simplicity,the finite integration interval of the angular momentawas split into two parts, − L C ≤ L ≤ ≤ L ≤ L C ,where L C is the angular momentum of a circle orbit,as mentioned above. There are the symmetric station-ary points, ±| L ∗ | , related to the anticlockwise and clock-wise motions of the particle along the PO in two thesephase-space parts. As noted above, they give equiva-lent contributions to the amplitude, due to the indepen-dence of the Hamiltonian of time. Thus, we have reducedthe integration region to 0 ≤ L ≤ L C , accounting forthis time-reversibility symmetry simply by the factor 2in Eq. (3.11) (exceptions are the diameters, for which onehas the single stationary point L ∗ = 0, and therefore, thetime-reversibility degeneracy is one, as it is taken into ac-count automatically by the limits of the error functions).For all the polygon-like and diameter POs ( n r ≥ L − = 0 for the minimum value of the angularmomentum L .For the Maslov index of the considered K = 1 POfamilies and the constant phase φ d in Eq. (3.10), oneobtains σ (1)PO = 2 M n r , φ d = − π/ . (3.14)The Maslov index σ PO is determined in terms ofthe number of turning and caustic points by theMaslov&Fedoryuk catastrophe theory, see Refs. [17, 21,22]. Note that for the potentials with smooth edges, theexpression for the Maslov index σ PO differs from that forthe circular billiard [11, 23]. Note also that the totalMaslov phase, defined as a sum of the asymptotic part(3.14) and the argument of the complex density ampli-tude (3.11), depends on the energy E and parameter α ofthe RPL potential [Eq. (1.1); see Refs. [9, 16]]. This totalMaslov phase is changed through the bifurcation pointssmoothly, due to the phase of the complex error functionin the amplitude (3.11) in Eq. (3.10).For the stationary point L ∗ far from the ends of thephysical integration interval, one can extend the integra-tion range to the infinity from −∞ to ∞ (in the case ofdiameters from zero to ∞ ). We then arrive asymptoti-cally at the Berry&Tabor result [5] for the contributionof all K = 1 families (3.10) with the following amplitude: A (1)PO → d PO T PO π ~ / p M n r K PO , (3.15)where d PO accounts for the discrete degeneracy, d PO = 1for diameters M (2 ,
1) ( n r = 2 n ϕ ), and 2 for all other(polygon-like) POs ( n r > n ϕ ) [11]. In the circular bil-liard limit ( α → ∞ ), the action is given by S PO ( E ) → p L PO with the momentum p = √ mE , and the POlength L PO . For the curvature K PO [Eqs. (2.11) and(B5)], one can asymptotically ( α → ∞ ) obtain K PO → / [ πpR sin ( πn ϕ /n r )] . Substituting all these quantities, S PO , K PO , σ (1)PO [with accounting for the Maslov-phasecontribution of the turning points due to the pure re-flections from the infinite circle walls [23] as comparedto smooth potentials[17] in addition to Eq. (3.14)], andthe asymptotic amplitude (3.15) into Eq. (3.10), one ob-tains the well known trace formula for the circular billiard[11, 23]. Note that the amplitude (3.11) of the solution(3.10) is regular at the bifurcations which are the bound-ary points L = L ∗ = L C of the action ( L ) part of thetori as in the elliptic billiard [16].Our SSPM result (3.15) coincides with the Berry andTabor trace formula [5], as adopted to the 2D spherically-symmetric Hamiltonians by using the simplest expan-sions of the action phase and amplitude near the station-ary point (see above), instead of a more general but morecomplicated mapping procedure; see more comments inRef. [16]. The essential difference from the Berry&Tabortheory [5] is that Eq. (3.10) covers all the solutions ofthe symmetry-breaking problem for the highest degen-erate orbits, such as the one-parametric families in theIHH potential, or the elliptic and hyperbolic orbits in theelliptic billiard [16] (see also Refs. [10, 17]). Within theSPM of the extended Gutzwiller approach [4, 10, 17], wehave to derive separately the contributions of the otherorbits as the circle K = 0 POs in the RPL potentials be-yond the semiclassical Poisson summation-like trace for-mula (3.3) (with the restrictions to the range of the n r and n ϕ integer variables). We emphasize that the ISPMtrace formula (3.10) for the one-parametric families con-tains the end contributions related to the finite limits ofintegrations in the error functions. However, this traceformula can be only applied to the contribution of suchfamilies, as pointed out above in its derivation from thetrace formula (3.1). Therefore, there is no contributionsof the circle orbits in Eqs. (3.3) and (3.10). As shown be-low, these orbits correspond to the separate contributionof the isolated ( K = 0 ) stationary-phase point L ∗ = L C (as for the IHH potential [17], for example). B. Circle orbits ( K = 0 ) In contrast to the derivations of contributions of theorbits with the highest degeneracy K = 1, we now takeinto account the existence of the isolated stationary pointof the action phase Φ CT (3.2) in the radial sphericalphase-space variables r ′ ∗ = r ′′ ∗ = r C , p ′ ∗ r = p ′′ ∗ r =0 . After the transformation of the integration vari-ables in Eq. (3.1) to the spherical phase space coordi-nates { r ′ , ϕ ′ ; p ′′ r , L } , it is convenient first to perform theexact integrations over L by using the energy conserving δ -function, and over the cyclic azimuthal angle ϕ ′ leadingsimply to 2 π as above ( R d ϕ ′ /ω ϕ = T ϕ, CT is the primitiverotation period). Thus, one finds g scl ( E ) = 2(2 π ~ ) Re X CT Z d r ′ Z d p ′′ r T ϕ, CT × |J CT ( p ′ r , p ′′ r ) | / exp (cid:20) i ~ Φ CT − i π µ CT − iφ d (cid:21) . (3.16)The additional factor 2 accounts for the equivalent con-tributions of two CTs for the particle motion in the twoopposite directions (with the opposite signs of the angu-lar momentum as above). The stationary phase conditionfor the SPM integration over the radial momentum p ′′ r inEq. (3.16) is written as (cid:18) ∂Φ CT ∂p ′′ r (cid:19) ∗ ≡ ( r ′ − r ′′ ) ∗ = 0 . (3.17)The solution of this equation is the isolated stationarypoint p ′′ r = p ′′ ∗ r = p ∗ r = 0. The phase Φ CT [Eq. (3.2)] isexpanded in the momentum p ′′ r near this point p ′′ ∗ r = 0in power series, Φ CT = Φ ∗ CT + 12 J ( p )CT ( p ′′ r − p ∗ r ) + · · · , (3.18)where the Jacobian is given by J ( p )CT = (cid:18) ∂ Φ CT ∂p ′′ r (cid:19) ∗ = " πM n r K ( ∂p ′′ r /∂L ) ∗ . (3.19)The star implies again that the corresponding quantity istaken at the stationary point, p ′′ r = p ′′ ∗ r = 0. Using the2nd order expansion of the exponent phase (3.18) andtaking the pre-exponent amplitude factor off the integralat this stationary point, one gets the internal integralover p ′′ r in Eq. (3.16) in terms of the error function as inthe previous section. According to Eq. (3.2), with theradial-coordinate closing condition (3.17) for the CTs,the short phase Φ ∗ CT in Eq. (3.18) can be written in termsof the corresponding variables as Φ ∗ CT = R r ′′ r ′ p r d r . Tak-ing then into account the CT closing condition (3.17), r ′ = r ′′ = r , for the stationary phase equation in theintegration over the radial r coordinate perpendicular tothe circle orbit, one results in (cid:18) ∂Φ ∗ CT ∂r ′′ + ∂Φ ∗ CT ∂r ′ (cid:19) ∗ ≡ ( p ′′ r − p ′ r ) ∗ = 0 . (3.20)Therefore, together with Eq. (3.17), one has the PO con-ditions related to the circular orbit r = r ∗ = r C and L = L ∗ = L C (see Appendix A). As usually within theSPM, we expand now the phase Φ ∗ CT in the radial coor-dinate r near this r ∗ = r C , Φ ∗ CT = M S C + 12 J ( r ) MC ( r − r C ) + · · · , (3.21)where S C is the action along the primitive circle PO ( C ), J ( r ) MC = (cid:18) − ∂p ′ r ∂r ′ − ∂p ′ r ∂r ′′ + ∂p ′′ r ∂r ′′ (cid:19) ∗ MC . (3.22)Again, using the action phase expansion (3.21) at thesecond order as the simplest ISPM approximation, andtaking the pre-exponent amplitude factor at the isolatedstationary point r = r C off the integral, one finally ob-tains δg (0) { MC } ( E ) = Re ∞ X M =1 A (0) MC ( E ) × exp (cid:20) i ~ M S C ( E ) − i π σ (0) MC − iφ (0) d (cid:21) . (3.23)The sum runs all repetitions of the circle orbit M C with M = 1 , , · · · being positive integers. The time-reversalsymmetry of the Hamiltonian (equivalence of the contri-butions of both angular momenta and repetition numberswith opposite signs) was taken into account by the factor2 in Eq. (3.16). The action S C ( E ) along the primitive C orbit is given by S C ( E ) = I C p ϕ d ϕ = 2 π L C (3.24)with L C shown explicitly in Eq. (A2). In Eq. (3.23), σ (0) MC is the Maslov index determined by the number of causticand turning points along the circle orbit, according tothe Fedoryuk& Maslov catastrophe theory[17, 21, 22], σ (0) MC = 4 M, φ (0) d = 0 . (3.25)For the amplitudes A (0) MC ( E ) in Eq. (3.23), one finds A (0) MC = T C π ~ √ F MC erf( Z ( − ) p,MC , Z (+) p,MC ) erf( Z ( − ) r,MC , Z (+) r,MC ) , (3.26)where T C is the period of the primitive ( M = 1) orbit C , T C = d S C ( E )d E = πL C α + 2 αE ; (3.27)see Eqs. (3.24), (B2) and (B1). In Eq. (3.26), F MC is the Gutzwiller stability factor [3] of the circle orbits[Eq. (2.9)]. The arguments of the error functions in Eq. (3.26) can be transformed to the following invariantform (see Appendix E): Z ( ± ) p,MC = r − i ~ π M √ α + 2 K C ( L ± − L C ) ,L + = L C , L − = 0 , Z ( ± ) r,MC = s i F MC π M ~ ( α + 2) / K C Θ ( ± ) r ,Θ (+) r = 2 π, Θ ( − ) r = 0 . (3.28)Here, L ± are the maximum and minimum values of theangular-momentum integration variable for the contri-bution of the circle orbits, K C is their curvature (seeAppendix E), K C = ( α + 1)( α − √ α + 2) L C . (3.29)The simplest approximation in Eq. (3.28) is L + = L C , L − = 0; and Θ − r = 0, Θ (+) r = 2 π , which correspond tothe total physical phase space accessible for the classicalmotion. The factors √ α + 2 in front of the curvature K C appear because of the frequency ratio f ( L ) = ω ϕ /ω r forthe circle orbits for any parameter α ≥
2; see Eqs. (2.8),(A3) and (A7). For α = 4; the period T C [Eq. (3.27)],action S C [Eq. (3.24)], curvature K C [Eq. (3.29)], andstability factor F MC [Eq. (2.9)] for the circle orbits areidentical to those obtained in Ref. [18]. We used also theproperties of the Jacobians for transformations of thedifferent coordinates, in particular, given by Eq. (E2).Note that after applying the stationary phase conditions r ∗ = r C [Eq. 3.17] and p ∗ r = 0 [Eq. (3.20)], the angularmomentum L of the circular orbits as function of the r and p r becomes the isolated stationary point L ∗ = L C atthe boundary of the classically accessible phase space.Notice also that the asymptotic Maslov phase is de-fined traditionally in terms of the Maslov index σ (0) MC [Eq. (3.14)]. There is again the two components of theMaslov phase in the ISPM trace formula (3.23) for the M C orbits. One of them is the asymptotic constant part(3.14) independent of the energy. Another part is theargument of the complex amplitudes A (0) MC [Eq. (3.26)],that changes continuously through the bifurcation points.The total Maslov phase for the circle POs is given by thesum of these two contributions, which ensures a smoothtransition of the trace formula (3.23) for the contributionof the circle POs through the bifurcation points.In the asymptotic limit of the non-zero integrationboundaries, L − → −∞ and Θ + r → ∞ , i.e., far fromany bifurcations α bif [Eq. (2.10), including the HO sym-metry breaking at α = 2], the expression (3.26) tends(through the Fresnel functions of the corresponding realpositive arguments) to the amplitude of the Gutzwillertrace formula for isolated orbits [3, 11], A (0) MC ( E ) → π ~ T C √ F MC . (3.30)In this limit, the asymptotic Maslov index σ (0) MC and φ (0) d in Eq. (3.10) are given by Eq. (3.25). Notice thatthe number coefficient in Eq. (3.30) differs from theSSPM Gutzwiller’s expression (5.36) of Ref. [11] by fac-tor 1/4. The reason is that the two stationary-phasepoints r ′ ∗ = r C and p ′′ ∗ r = 0 belong to the boundaryof the physical { r ′ , p ′′ r } phase-space integration volumein Eq. (3.16); while in Ref. [3], all the stationary pointsare assumed to be internal ones which are far away fromthe integration boundary. Eq. (3.30) can be derived di-rectly from Eq. (3.16) by using the SSPM. To realizethis within the SSPM, one may extend in Eq. (3.16) the r ′ integration range from { r ′ = 0 , r C } to {−∞ , r C } , andsimilarly, the p ′′ r integration one to { , ∞} , assuming thatthe lower r ′ and upper p ′′ r integration limits are far awayfrom the corresponding other (stationary-point) integra-tion boundaries.For the opposite limit to the bifurcations ( F MC → α → α bif ), one finds that the both arguments ofthe second error function in Eq. (3.26) tend to zero as p | F MC | , see Eq. (3.28). The Gutzwiller stability factor F MC , going to zero, is exactly canceled by the same onein the denominator, and we arrive at A (0) MC ( E ) → T C ~ / p πM ( α + 2) / K C × erf (cid:16) Z ( − ) p, MC , Z (+) p, MC (cid:17) e iπ/ . (3.31)Thus, in contrast to the SSPM divergences, one obtainsthe finite results at the bifurcations within the ISPM.Notice that the enhancement in order of ~ − / withrespect to the Gutzwiller asymptotic amplitude (3.30)takes place locally near the bifurcation points. Note alsothat at the circular billiard limit, when K C → ∞ (sepa-ratrix), one finds a continuous limit which is zero in thecase of the RPL potential. C. Total trace formula for the oscillating leveldensity
The total semiclassical oscillating (shell) correction tothe level density (3.1) for the RPL potentials in two di-mensions is thus given by δg scl ( E ) = δg (1)scl ( E ) + δg (0)scl ( E ) , (3.32)where δg ( K )scl ( E ) = Re X PO A ( K )PO ( E ) × exp (cid:20) i ~ S PO ( E ) − i π σ ( K )PO − iφ ( K ) d (cid:21) . (3.33)The amplitudes A ( K )PO [see Eqs. (3.11) for K = 1 and (3.26)for K = 0], actions S PO , Maslov indexes σ ( K )PO , and con-stant phases φ ( K ) d [Eqs. (3.14) and (3.25)] were specifiedabove. Using the scale invariance (2.2), one may factorize theaction integral S PO ( E ) = (cid:18) EE (cid:19) + α I PO( E = E ) p · d r ≡ ετ PO . In the last equation, we define the scaled energy ε andscaled period τ PO by ε = (cid:18) EE (cid:19) + α , τ PO = I PO( E = E ) p · d r . (3.34)To realize the advantage of the scaling invariance (2.2),it is helpful to use the scaled energy (period) in place ofthe corresponding original variables. For the HO, one has α = 2, and the scaled energy and period are proportionalto the unscaled quantities. For the cavity potential ( α →∞ ), they are proportional to the momentum p and length L PO , respectively.Using the transformation of the energy E to the scaledenergy ε , one can introduce the dimensionless scaled-energy level density. The advantage of this transforma-tion is that a nice plateau condition is always found inthe Strutinsky SCM smoothing procedure by using thescaled spectrum ε i (see Refs. [9, 16] for the case of the bil-liard limit α → ∞ ). Then, one can use a simple relationbetween the original and scaled-energy level densities, G ( ε ) = X i δ ( ε − ε i ) = g ( E ) d E d ε . (3.35)For the semiclassical oscillating part of the level density(3.35), one finds δ G ( K )scl ( ε ) = d E d ε δg ( K ) ( E ) = X PO δ G ( K )PO ( ε )= Re X PO A ( K )PO ( ε ) exp (cid:20) i ~ ετ PO − iπ σ ( K )PO − iφ ( K ) d (cid:21) , A ( K )PO ( ε ) = d E d ε A ( K )PO ( E ) . (3.36)The simple form of the phase function (3.34) enables usalso to make easy use of the Fourier transformation tech-nique. The Fourier transform of the semiclassical scaled-energy level density with respect to the scaled period τ is given by F ( τ ) = Z d ε G ( ε ) e iετ/ ~ ≈ F ( τ ) + X PO e A PO δ ( τ − τ PO ) , (3.37)which exhibits peaks at periodic orbits τ = τ PO . F ( τ )represents the Fourier transform of the smooth Thomas-Fermi level density and has a peak at τ = 0 related tothe zero-action trajectory [10]. Thus, from the Fouriertransform of the scaled-energy quantum-mechanical leveldensity (3.35), F ( τ ) = X i e iε i τ/ ~ , ε i = (cid:18) E i E (cid:19) + α , (3.38)one can directly extract the information about classicalPO contributions. The trace formula (3.32) has the cor-rect asymptotic SSPM limits to the Berry&Tabor results(3.10), (3.15) for K = 1 polygon-like (including the diam-eters) and to the Gutzwiller trace formula (3.23), (3.30)for K = 0 circle POs. As shown in the sections III Aand III B, one obtains also the limit of the trace formula[Eqs. (3.32) and (3.33)] to that of the circular billiard α → ∞ [11, 23]. In this limit one has obviously zero forthe circle orbit contributions as for the potential barrierseparatrix in the IHH potential [17].For comparison with the quantum level densities ob-tained by the SCM, we need also to perform a local aver-aging of the trace formula (3.32) over the spectrum. Asthis trace formula is given through the sum of the in-dividual PO terms everywhere (including the bifurcationregions), one can approximately take the folding integralsover energies in terms of the Gaussian weight factors witha width parameter Γ ≪ E F . As the result, one obtainsthe Gaussian-averaged oscillating level density in the an-alytical form [4, 10, 11]: δg Γ ( E ) = X PO δg PO ( E ) exp h − ( t PO Γ/ ~ ) i . (3.39)Adding the TF smooth component g TF ( E ) [11] to thisoscillating component, one results in the total trace for-mula: g Γ ( E ) = g TF ( E ) + δg Γ ( E ) , (3.40)where g TF ( E ) = 1(2 π ~ ) Z d r Z d p δ (cid:18) E − p m − V ( r ) (cid:19) = mr ~ = 12 E (cid:18) EE (cid:19) /α . (3.41)Here, r max is the maximal turning point (one of solutionsof the equation V ( r ) = E ), which is given by r max = R ( E/E ) /α for the RPL Hamiltonian (2.1), and we put E = ~ /mR in the last expression of Eq. (3.41).Using the scaled-energy transformation (3.35) of theoscillating part (3.39) of the Gaussian-averaged level den-sity [Eq. (3.40)], one finally obtains the semiclassicalscaled-energy trace formula: δ G γ ( ε ) = X K =0 δ G ( K ) γ ( ε )= X K =0 X PO δ G ( K )PO ( ε ) exp (cid:20) − (cid:16) τ PO γ ~ (cid:17) (cid:21) , (3.42)Here, δ G ( K )PO ( ε ) is given by Eq. (3.36), γ is a dimensionlesswidth parameter used for the Gaussian averaging over thescaled spectrum ε i . For the scaled-energy Thomas-Fermidensity component, one finds G TF ( ε ) = g TF ( E ) d E d ε = α α ε . (3.43) D. The shell correction energies
The semiclassical PO shell correction energies δU scl isgiven by [4, 9–11, 16] δU scl = 2 X PO ~ t PO δg PO ( E F ) , (3.44)where t PO = M T PO ( E F ) is the period of particle motionalong the PO (taking into account its repetition number M ) at the Fermi energy E = E F . The Fermi energy E F as function of the particle number N is determined bythe particle number conservation, N = 2 X i n i = 2 Z E F d E g ( E ) , (3.45)where n i = θ ( E F − E i ) are the occupation numbers. Thefactors 2 in Eqs. (3.44) and (3.45) account for the spindegeneracy of Fermi particles with spin 1/2.Note that the shell correction energies δU which are theobserved physical quantities do not contain an arbitraryaveraging parameter Γ , in contrast to the level density g Γ ( E ). The convergence of the PO sum (3.44) to shorterPOs (if they occupy enough large phase-space volume) isensured by the additional factor in front of the oscillatingdensity components δg PO which is inversely proportionalto square of the PO period t PO .In the quantum SCM calculations, the shell correctionenergies are usually obtained by extracting the oscillat-ing part from a sum of the single-particle energies, Notethat the direct application of the SCM average procedureto the spectra E i of RPL potentials (except for the HOlimit) does not give any good plateau condition as for thelevel density g ( E ) in Eq. (3.35). However, one may findrather a good plateau in the SCM application to a sumof the single-particle scaled energies ε i , U = 2 P i n i ε i .Applying exactly the same derivations of Eq. (3.44) tothe semiclassical trace formula for the oscillating part of U , one gets δ U scl = 2 X PO ~ τ δ G PO ( ε F ) . (3.46)Here, the scaled Fermi energy ε F is determined by N = 2 Z ε F G ( ε )d ε . (3.47)Using now the obvious relations t PO = τ PO d ε/ d E and δg PO ( E ) = δ G PO ( ε ) d ε/ d E in Eq. (3.44), one obtains δU scl = (cid:18) d E d ε (cid:19) ε F δ U scl . (3.48)Thus, we arrive at the simple relation between the origi-nal shell-correction energy δU [Eq. (3.44)] and the scaled0one δ U , valid for both semiclassical and quantum (ne-glecting the second order terms in the shell fluctuationsof the Fermi energy) calculations: δU = (cid:18) d E d ε (cid:19) ε F δ U = E αα + 2 ε ( α − / ( α +2) F δ U . (3.49)This relation can be also directly obtained by using thestandard quantum SCM relations of the first-order shell-correction energy δU to the oscillating part of the leveldensity δg ( E ) up to the same second order terms in theFermi energy oscillations [2], and corresponding ones forthe scaled quantities, δ U = 2 X i δn i ε i = 2 Z ε F d ε ( ε − ε F ) δ G ( ε ) . (3.50)In these derivations, δn i = n i − e n i represents the oscillat-ing part of the occupation number defined by subtractingthe smooth part e n i from the exact one. We applied alsothe usual transformations from the Fermi energies to theparticle numbers by using Eqs. (3.45) and (3.47), as wellas the definitions of the averaged Fermi energy e E F , andthe scaled one e ε F , N = 2 Z e E F d E e g ( E ) = 2 Z e ε F d ε e G ( ε ) . (3.51) E. Harmonic oscillator limit
In the isotropic harmonic oscillator limit [ α → E = ω r I r + ω ϕ I ϕ = ω ϕ (2 I r + L ) . (3.52)Therefore, in this limit the curvature K PO for all POs[including the maximum value L = L C = E/ω ϕ for thecircle orbits, Eq. (3.29), and L = 0 for diameter ones,Eq. (D3)] and stability factor F MC [Eq. (2.9)] for the M C orbits turn into zero. However, there is no singularitiesin the ISPM trace formulas (3.10) for the contributionsof all K = 1 families and (3.23) for the circle orbits in thelimits K PO → F MC →
0. The arguments of botherror functions, ∝ √ K C and ∝ p F MC /K C in Eq. (3.26),for instance, approach zero and singularities are canceledwith the same ones in the denominators of the multipli-ers in front of them, and similarly, in Eq. (3.10) for oneerror function; see Eqs. (3.11), (3.13), (3.26) and (3.28)with the help of Eq. (3.31). Therefore, one has a con-tinuous limit of the total trace formula (3.32) for α → M C orbit contribution(3.23) and the M (2 ,
1) diameter one [Eq. (3.10)] up tothe relatively small higher-order corrections in ~ [see alsoEq. (3.68) of Sec. 3.2.4 in Ref. [11]], g (0) { MC } ( E ) → δg (2)HO ( E ) , g (1) { MD } ( E ) → δg (2) HO ( E ) . (3.53) Here, { M C } and { M D } represent sum of all repetitionsof circle and diameter orbits, M = 1 , , ... , respectively.Thus, the HO limit of the sum of the circle and diame-ter orbit contributions into the (averaged) level densityand the energy shell corrections is exactly analyticallygiven by the corresponding HO trace formulas. We pointout that for α →
2, the contributions of circle
M C anddiameter M (2 ,
1) orbits encounter local increases of thedegeneracies K by 2 and 1 units, respectively.As noted above, in the HO limit α →
2, only the diam-eter M (2 ,
1) and the circle
M C (both with repetitions)survive, and they form K = 2 families in the HO po-tential. Taking into account also that the angular mo-mentum for the diameters is always zero, L ∗ = 0, andfor the circle orbits L ∗ = L C , we shall assume that theintegration over L for the diameters is performed from L − = 0 to L + = L C / L − = 0 to L + = L C , such that they give naturallyequivalent contributions into the HO trace formula, asshown in Eq. (3.53), see also Ref. [17]. The difference is inthe integration limits for the circle orbits [Eq. (3.28)], incontrast to Eqs. (3.13) and (D2) for the diameter bound-aries. Notice that the contribution of the polygon-likeone-parametric orbits, δg (1) ( E ), disappears in the ghostHO limit. Thus, one obtains the continuous transitionof the oscillating part of the ISPM level density δg scl ( E )through all bifurcation points, including the HO symme-try breaking. IV. AMPLITUDE ENHANCEMENT ANDCOMPARISON WITH QUANTUM RESULTS
A remarkable enhancement of the ISPM amplitudes inPO sum for the oscillating level density (3.33) and shellcorrection energy (3.44) due to the bifurcation (symme-try breaking) is typically expected for some short periodicorbits. In Fig. 2, the scaled amplitudes |A PO | , dividedby ε / to normalize the energy dependence for K = 1orbits, are presented for several shortest POs as func-tions of the power parameter α in order to show the typ-ical bifurcation enhancement phenomena. In Fig. 2 (a) ,the enhancement of the primitive diameter (2 ,
1) ampli-tudes |A (2 , | [Eq. (3.11)], and those |A MC | [Eq. (3.26)]for the primitive circle orbit C are clearly seen in theHO limit α →
2; see also Eq. (3.36). Figure 2 (b) showsthe enhancement of the shortest orbit C around the bi-furcation point α = 7, and the birth of the triangle-likeorbit (3 ,
1) there. Note that the ISPM amplitude |A (3 , | [Eq. (3.11)] of the (3 ,
1) orbit keeps its magnitude up torather a large value of α above the bifurcation ( α > α bif ).The ISPM amplitude for the circle orbit C exhibits a re-markable enhancement at the bifurcation point α = 7.The divergence of its SSPM amplitude at the bifurca-tion point is successfully removed. As also seen fromFig. 2 (b) , the ISPM amplitude for the (3 ,
1) PO is con-tinuously changed through this bifurcation, in contrastto the discontinuity of the SSPM amplitude. This orbit1 | A ( ε ) | / ε / | A ( ε ) | / ε / C ISPMC ISPMC SSPM(3,1) SSPM(3,1) ISPM(2,1) ISPM (2,1) SSPM ε=40 αα (a)(b) FIG. 2. Moduli of the scaled ISPM (solid) and SSPM (dashedcurve) amplitudes |A PO ( ε ) | as functions of α for the primitive( M = 1) circle C [Eq. (3.26)], diameter (2 ,
1) and triangle-like(3 ,
1) [Eq. (3.11)] POs, in units of ε / at the scaled energy ε = 40. The panels show (a) the HO limits ( α →
2) for thecircle C (thin) and diameter (2 ,
1) (thick curve) orbits (thefilled circle denotes one half of the HO amplitude (3.53) at α = 2); and (b) the circle C (thin) and triangle-like (3 , exists, in fact, only at α ≥
7, and the amplitude in theregion α < |A (3 , | in the ghost region far from the bifurcationhas no physical significance, since it is washed out in theGaussian-averaged level density by an rapidly oscillat-ing phase of the complex amplitude A (3 , [9, 17]. Theseghost amplitude oscillations are suppressed even moreby using higher order expansions in the phase and am-plitudes in a more precise ISPM [9].Figures 3–7 show the oscillating part of the semiclassi-cal scaled-energy level density δ G γ ( ε ) [Eq. (3.42)] in unitsof ε / as functions of the scaled energy ε for several val-ues of the power parameter α and the Gaussian width γ . The ISPM semiclassical results show good agreementwith the quantum mechanical (QM) ones for a transitionfrom the gross to fine resolutions of the spectra. TheQM calculations are carried out by the use of the stan-dard Strutinsky averaging over the scaled energy ε , inwhich we find a good plateau around the Gaussian aver-aging width e γ = 2 − α = 4 . .
0, one finds a good agree-ment with the SSPM asymptotic behavior [Figs. 4 (a) and5 (a) ] because they are sufficiently far from the bifurcationpoints α = 4 .
25 and 7 . ,
2) and triangle-like (3 ,
1) POs (Figs. 6and 7). For the gross shell structure ( γ ≈ . α = 4 .
29 30 31 32 33 34 35 36 37-4-20246 29 30 31 32 33 34 35 36 37-4-20246 QM ISPM δ ( ε ) / ε G γ / G QM ISPM δ ( ε ) / ε / γ ε α=4.0, γ=0.1α=6.0, γ=0.1 (a)(b) FIG. 3. The oscillating part of the level density δ G γ ( ε ) in unitsof ε / vs the scaled energy ε for α = 4 . (a) , and α = 6 . (b) , at the (dimensionless) width parameter γ = 0 . give the leading contributions. (This is in contrast to the3D case where the circular orbits become also important[11, 18].) For instance, the gross shell structure in termsof the shortest POs for α = 6 . − . γ & .
3, unlike for the powers α = 4 . − .
25. Withdecreasing γ and increasing α , the POs for larger scaledperiods τ [or actions S , see Eq. (3.34)] become more sig-nificant [cf. Figs. 4 (b,d) and 7 (b,d) ]. In the case of thefine shell structure (e.g., γ ≈ .
03) the dominant contri-butions are due to the bifurcating K = 1 POs [polygon-like POs denoted by { M P } ; see Fig. 4 (b) ]. (This is simi-lar to the situations in the elliptic [16] and spheroidal [9]cavities, and in the IHH potential [17].) However, theinterference of these much longer one-parametric POs[such as M (7 ,
3) for α = 4 . M (5 ,
2) for α = 6 . M (2 ,
1) diameters explain some peaks,too. For smaller α = 4 . .
25, the circle orbit con-tributions are not shown because they are insignificantat these power parameters in the 2D case. (This is dif-ferent situation from the 3D case, see Ref. [18] for thetrace formulas based on the uniform approximation us-ing the classical perturbation approach [11, 24].) Thesecontributions into the trace formula (3.42) are increasingfunctions of α , and they become significant at α & (b,d) ]. An intermediate situationbetween the gross- and fine- shell structures where all ofPOs become significant are shown too at γ = 0 . γ = 0 . . (b,d) . Ourfull analytical expressions (accessible for any long peri-odic orbits) for the classical PO characteristics at α = 4and 6 are quite useful in the simple ISPM calculationsof the oscillating level density with a good accuracy upto the fine spectrum-structure resolutions by using, forinstance, γ ≈ .
03 and 0 .
1. Figures 6 and 7 show a nice2
27 28 29010 27 28 29010 α=4.0, γ=0.03 δ ( ε ) / ε G γ / QM ISPM SSPM δ ( ε ) / ε γ / G ε {MP}{MD} α=4.0, γ=0.03 {MP}{MD} (a)(b)
20 22 24 26 28 30 32 34 36 38 40 42-2-101220 22 24 26 28 30 32 34 36 38 40 42-2-1012 δ ( ε ) / ε G γ / α=4.0, γ=0.2 QM ISPM {MP} δ ( ε ) / ε γ G / {MD} α=4.0, γ=0.2 ε {MD} (c)(d) FIG. 4. The same as in Fig. 3 for α = 4 . γ = 0 . (a,b) and γ = 0 . (c,d) . Panels (a,c) :the solid, dashed and dotted lines are the QM, ISPM andSSPM [ the panel (a) ] results.
Panels (b,d) : { MD } (dashed)is the contribution of the diameters (including their repeti-tions) and { MP } (thin solid) for other K = 1 polygon-likePOs. δ ( ε ) / ε G γ / α=6.0 γ=0.03 ε QM ISPM ε SSPM δ ( ε ) / ε / G γ α=6.0 γ=0.2 QM ISPM (a)(b)
FIG. 5. The same as in Fig. 3 for α = 6 .
0, but with otherwidth parameters, γ = 0 . (a) , and γ = 0 . (b) .
30 31 32 33 34 35 36 37 38 39 40-4-2024620 22 24 26 28 30 32 34 36 38-2-1012 δ ( ε ) / ε / G γ δ ( ε ) / ε / γ G α=4.25 γ=0.1 εε ISPM QMQM α=4.25 γ=0.2 ISPM (a)(b)
FIG. 6. The same as in Fig. 5 for α = 4 .
25, but with thewidth parameters, γ = 0 . (a) , and γ = 0 . (b) .
21 22 23 24 25 26 27 28-4-2024 21 22 23 24 25 26 27 28-4-2024 α=7.0 γ=0.1 {MC} {MP} QM {MD} ISPM {MP} ε δ ( ε ) / ε G γ / δ ( ε ) / ε / G γ {MD} α=7.0 γ=0.1 {MC} (a)(b)
24 26 28 30 32-2-10123 24 26 28 30 32-2-10123 α=7.0 γ=0.2 δ ( ε ) / ε G γ / δ ( ε ) / ε γ G / ε QM{MD}{MP} {MC} α=7.0 γ=0.2 ISPM {MP} {MD} (c)(d)
FIG. 7. The same as in Fig. 4, but for α = 7 . γ = 0 . (a,b) and γ = 0 . (c,d) .
18 20 22 24 26 28 30 32-0.2-0.100.10.218 20 22 24 26 28 30 32-0.2-0.100.10.2 α=4.0 F δ / ε / N QM ISPM δ / ε F / {MD} {MP} α=4.0 {MD} (a)(b) UU
18 20 22 24 26 28 30-0.2-0.100.10.2 18 20 22 24 26 28 30-0.2-0.100.10.2 δ / ε F / QM α=6.0 δ / ε / F N {MD} α=6.0 ISPM{MP}{MD} (c)(d) UU FIG. 8. Scaled shell correction energies δ U , normalized by thefactor ε − / F , as functions of square root of the particle num-ber N / at the values of α , where the full analytical formulasare obtained for α = 4 . (a,b) and α = 6 . (c,d) . Panels (a,c) :QM (solid curve) represents the quantum-mechanical resultsusing the Strutinsky SCM, and ISPM (dashed curve) showsthe semiclassical result using the TF approximation in thecalculation of N ( ε F ) by Eqs. (3.47) and (4.1). Panels (b,d) :the contributions of several POs into the shell correction en-ergy δ U are shown. Other notations are the same as in Figs. 4and 7. agreement of the fine-resolved semiclassical and quantumlevel densities δ G γ ( ε ) as functions of the scaled energy ε at the critical bifurcation points α = 4 .
25 and 7 . ,
2) and triangle-like (3 ,
1) orbits,respectively.Figures 8 and 9 show the scaled shell correction ener-gies δ U [Eqs. (3.46) for the semiclassical and (3.50) forthe quantum results], normalized by the factor ε − / F ,as functions of the particle number variable N / . Agood plateau is realized for the QM calculations of thescaled shell-correction energies [see the first equation inEq. (3.50)] near the same averaging parameters e γ andcurvature corrections as mentioned above. In the semi-
16 18 20 22 24 26 28 30-0.2-0.100.10.216 18 20 22 24 26 28 30-0.2-0.100.10.2 α=4.25 α=4.25 δ / ε / F δ / ε F / N QM {MD}ISPM{MP}{MD} (a)(b) UU
17 18 19 20 21 22 23 24-0.2-0.100.10.2 17 18 19 20 21 22 23 24-0.2-0.100.10.2 δ / ε F / δ / ε U / α=7.0 α=7.0 N QM {MD}{MP} {MC}ISPM {MD} {MP} (c)(d) U F FIG. 9. The same as in Fig. 8; α = 4 .
25 for the (5,2) bifurca-tion, and α = 7 . (a,b) and (c,d) , respectively. classical calculations, the Fermi level ε F is determined bythe particle number conservation (3.47) with using thecoarse-grained scaled-energy POT level density, G γ, scl ( ε ) = G TF ( ε ) + X K =0 δ G ( K ) γ, scl ( ε ) . (4.1)The oscillating ISPM components δ G ( K ) γ, scl ( ε ) are givenby Eqs. (3.42) and (3.36). We evaluated the Fermi level ε F ( N ) by varying the averaging width γ and found thatthere is no essential sensitivity within the interval ofsmaller γ ( γ ≈ . − . G TF ( ε ) [Eq. (3.51] in Eq. (3.47) with G ( ε ) ≈ G TF ( ε ) pro-vides us a good value of ε F in the POT calculations of theshell correction energies (3.46). The PO sums at α = 7 . γ = 0 . M = 4) for the circle and diame-ter orbits, and a few first simplest other K = 1 ( P ) POs,4such as (3 , , ,
3) and (8 , . α .
6, one hasa similar PO convergence relation with the same γ ≈ . K = 1) PO families (P)are the (5 , , (7 ,
3) and (7 ,
3) POs at α = 4 . − . .
0, respectively [Figs. 8(a,b) and 4(c,d)]. As seen fromFigs. 6, 7 and 9, we obtain a nice agreement between thesemiclassical (ISPM, dashed) and quantum (QM, solidcurve) results exactly at the bifurcations α = 4 .
25 and7 .
0. Notice that the dominating contributions in thesesemiclassical results at the bifurcation point α = 7 . C and newborn (3 ,
1) orbits with the simplest diameters.As shown typically in Figs. 7 (d) and 9 (d) , one can seethat the circle C and triangle-like (3 ,
1) orbits are mainlyin phase, but the diameter (2 ,
1) is sometimes in phase tothem and sometimes out of phase. Thus, the occurrenceof a characteristic beating pattern in the level densityamplitude at α = 7 . C and (3 ,
1) with the shortest diam-eter (2 ,
1) having all the amplitude of the same orderin magnitude but different phases. The bifurcating cir-cle 2 C and star-like (5 ,
2) orbits [as expected from theenhancement of the amplitudes of the circular C andtriangular-like (3 ,
1) POs in Fig. 2] are more importantfor α = 4 .
25, though the primitive diameters become sig-nificant much compared to the bifurcation case α = 7 . ,
1) and (5 ,
2) yield more contributions neartheir bifurcation values of α , and even more on the right-hand side ( α & α bif ) in a wide region of α as mentionedabove. The bifurcation parent-daughter partner orbits { C , (3 , } and { C , (5 , } , taken together with the sim-ple diameter (2 , (b,d) and 9 (d) for the same α = 7 .
0. The diameter ISPM contributions are close tothe SSPM asymptotic ones near the bifurcation points α = 7 . .
25 (as for α = 4 . .
0) because theyare sufficiently far from their single symmetry-breakingpoint at the harmonic oscillator value α = 2.Figure 10 shows the Fourier transform of the quantum-mechanical scaled-energy level density [Eq. (3.38)]. Fora smaller α = 2 .
1, the diameter (2 ,
1) orbit gives thedominant contribution to the gross-shell structure as theshortest POs; see the peak at τ ∼ .
0. With increas-ing α , the amplitude of the circle orbit becomes againlarger due to a prominent enhancement around the bi-furcation point ( τ ∼ . α bif = 7 . , , ,
3) and (8 ,
3) give compa-rable contributions at α = 7 . ,
2) and (7 , α = 4 .
25] in nice agreement with thequantum Fourier spectra in Fig. 10. The contributions ofthe newborn triangle-like orbit family (3 ,
1) having rela-tively a smaller scaled period τ (3 , and higher degener-acy K = 1 become important and dominating for larger α & α bif = 7. The newborn (3 ,
1) peak cannot be distin-guished from the parent circle C orbit near the bifurca- | F ( τ ) | τα = . | F ( τ ) | α = . | F ( τ ) | α = . | F ( τ ) | α = . | F ( τ ) | α = . | F ( τ ) | α = . ( , ) C ( , ) ( , )( , ) C ( , ) ( , )( , ) C ( , ) C C ( , ) C ( , ) C ( , ) ( , ) C ( , ) ( , ) C ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) ( , ) C ( , ) ( , )( , ) ( , ) C ( , ) ( , ) ( , ) ( , ) ( , ) C C ( , ) C ( , ) C ( , ) C ( , ) C ( , ) C FIG. 10. Moduli of the Fourier transform | F ( τ ) | of the quan-tum scaled-energy level density (3.38) as functions of the di-mensionless variable τ are plotted for several values of α ; MC and M ( n r , n ϕ ) indicate the classical POs corresponding toeach peak (see Fig. 1). tion point α bif as well as the diameter and circle orbits at α close to the HO limit, see Ref. [10]. We emphasize thatthe shell correction energies δ U are similar to the oscil-lating parts of the level densities coarse-grained over thespectrum by using the Gaussian width γ = 0 . α = 4 . .
25, which have mainly the gross-shell structure dueto the shortest diameters. However, for this γ , the fine-resolved shell structures (due to their interference withthe other polygon-like and circular POs) are pronouncedat larger powers near α = 6 .
0, and especially, 7 . α → V. CONCLUSIONS
We presented a semiclassical theory of quantum oscil-lations of the level density and energy shell corrections for5a class of radial power-law potentials which turn out asgood approximations to the realistic Woods-Saxon poten-tial in the spatial region where the particles are bound.The advantage of the RPL potentials is that, in spiteof its diffuse surface, the classical dynamics scaling withsimple powers of the energy simplifies greatly the ana-lytical POT calculations. The quantum Fourier spectrayield directly the contributions of the leading classicalPOs with the specific periods and actions into the traceformulas.We described the main PO properties of the classicaldynamics in the RPL potentials as the key quantitiesof the POT. Taking the simplest two-dimensional RPLHamiltonian we developed the semiclassical trace formu-lae for any its power α , and studied various limits of α (the harmonic oscillator potential for α = 2 and thecavity potential for α → ∞ ). The completely analyt-ical results were obtained for the RPL powers α = 4and 6. This can be applied for both 2D and 3D casesand allow us to far-going fine-resolved shell structuresat γ = 0 . − .
1. This POT is based upon extendedGutzwiller’s trace formula, that connects the level den-sity of a quantum system to a sum over POs of the cor-responding classical system. It was applied to expressthe shell correction energy δU of a finite fermion systemin terms of POs. We obtained good agreement betweenthe ISPM semiclassical and quantum-mechanical resultsfor the level densities and energy shell corrections at sev-eral critical powers of the RPL potentials. For the pow-ers α = 4 and 6, we found also good agreement of theISPM trace formulas with the SSPM ones. The strongamplitude-enhancement phenomena at the bifurcationpoints α = 7 and 4 .
25 in the oscillating (shell) compo-nents of the level density and energy were observed inthe remarkable agreement with the peaks of the Fourierspectra. We found a significant influence of the PO bi-furcations on the main characteristics (oscillating compo-nents of the level densities and energy-shell corrections)of a fermionic quantum system. They leave signatures inits energy spectrum (visualized, e.g., by its Fourier trans-form), and hence, its shell structure. We have presenteda general method to incorporate bifurcations in the POT,employing the ISPM based on the catastrophe theory ofFedoryuk and Maslov, and hereby, overcoming the diver-gence of the semiclassical amplitudes of the Gutzwillertheory and their discontinuity in the Berry&Tabor ap-proach at bifurcations. The improved semiclassical am-plitudes typically exhibit a clear enhancement near a bi-furcation and on right side of it, where new orbits emerge,which is of the order ~ − / in the semiclassical parame-ter ~ . This, in turn, leads to the enhanced shell struc-ture effects. Bifurcations are treated, again, in the ISPMleading to the semiclassical enhancement of the orbit am-plitudes. The trace formulae are presented numericallyto show good agreement with the quantum-mechanicallevel density oscillations for the gross- (coarse-grainedwith larger averaging width γ and a few shortest POs),and the fine-resolved (with smaller γ and longer bifur- cating POs) shell structures. The PO structure of theshell-correction energies is similar to that of the coarse-grained densities for smaller powers α = 4 – 4 .
25, and ofthe fine-resolved densities for larger α & γ ≈ .
2. The fine-resolved and coarse-grained shell struc-tures were found at the same α in the corresponding av-eraged oscillating densities at smaller width parameters γ = 0 .
03 – 0 . γ & . .
3, respec-tively. The fine-resolved shell structure for larger powers, α &
6, occurs in a larger interval, γ = 0 .
003 – 0 .
2, includ-ing the essential contributions of the circle orbits alongwith the polygon-like and diameter orbits. Full explicitanalytical expressions for the diameters and circle orbitcontributions into the trace formula as functions of thediffuseness potential parameter α are specified too.For prospectives, we intend a further study of shellstructures in the 3D RPL potentials, within the ISPMand uniform approximations to treat the bifurcations, byvarying continuously the power parameter α from 2 (har-monic oscillator) to ∞ (spherical billiards). ACKNOWLEDGMENTS
Authors thank Profs. M. Brack and K. Matsuyanagifor many valuable discussions.
Appendix A: The stability factor, bifurcation powersand frequencies
Let us consider in more details the non-linear classicaldynamics in the RPL Hamiltonian (2.1) for any real α ≥
2. The critical values of the radial coordinate r = r C and angular momentum L = L C for the circle orbit ( C )are determined by the solutions of the system of the twoequations with respect to r and L : F ( r, L ) = 0 , ∂ F ∂r = 0 , where F ( r, L ) ≡ p r ( r, L ) , (A1)see Eq. (2.3). In the internal region where the stableorbits in the radial direction exist, one has a nonzero F ′′ C = ∂ F ( r C , L C ) /∂r <
0. First equation in Eq. (A1)means that there is no radial velocity, ˙ r = 0, and the nextequation is that the radial force is equilibrating by thecentrifugal force. For the Hamiltonian (2.1), the solutionsof the two these equations are the radius r C and angularmomentum L C [13], r C = R (cid:18) E (2 + α ) E (cid:19) /α , L C = p ( r C ) r C . (A2)Using Eq. (2.6) at L = L C for the rotational frequency, ω C = ω ϕ ( L = L C ) = L C / ( mr C ), and (A2) for r C and L C ,one finds [13] ω C = s αE mR (cid:18) E (2 + α ) E (cid:19) / − /α . (A3)6Applying now the second order expansion in r − r C toEq. (2.3), one gets the first-order ordinary differentialequation for the radial CT r ( t ) locally near the circle PO r = r C : ˙ r = ± r F ′′ C m ( r − r C ) . (A4)Integrating the dynamical equation in Eq. (A4), one ob-tains r ( t ) = r C + ( r ′ − r C ) exp ± r F ′′ C m t ! , (A5)where r ′ = r ( t = t ′ = 0). In the stable case, F ′′ C < r ( t ) locally near the circle orbit r = r C , one writes r ( t ) = r C + ( r ′ − r C ) exp ( ± iΩ C t ) , (A6)where Ω C is a positive radial frequency ω r at L = L C [Eq. (2.6)], Ω C = q |F ′′ C / (2 m ) | = ω r ( L = L C ) . (A7)For the Hamiltonian (2.1), this quantity is given by [13] Ω C = s αEmR (cid:20) (2 + α ) E E (cid:21) /α > . (A8)From Eq. (A6) after the period T C along the primitivecircle orbit, T C = t ′′ − t ′ = t ′′ = 2 πω C , (A9)one finds δr ′′ ≡ r ′′ − r C = δr ′ exp ( ± iΩ C T C ) , δr ′ = r ′ − r C . (A10)The eigenvalues of the stability matrix M C for M = 1in Eq. (2.9) are given by [11] (cid:18) ∂r ′′ ∂r ′ (cid:19) p ′ r = exp ( iΩ C T C ) , (cid:18) ∂p ′′ r ∂p ′ r (cid:19) r ′ = exp ( − iΩ C T C ) . (A11)These two eigenvalues of the stability matrix are complexconjugated in agreement with its general properties. As Ω C is real [ Ω C >
0, according to Eqs. (A7) and (A8)]the circle orbit is isolated stable PO. Substituting theexpressions (A11) into the first equation in Eq. (2.9) andusing Eqs. (A9) for the period T C , (A3) and (A8) for theC orbit frequencies ω C and Ω C , relatively, one obtains thelast equation in Eq. (2.9) for the stability factor F MC . Appendix B: Scaling properties
For convenience, let us consider the classical dynam-ics in terms of the variables in dimensionless units m = R = E = 1. Due to the scaling property (2.2) for theclassical dynamics in the Hamiltonian (2.1), the energydependence of the action I r ( ε ) [Eq. (2.4)], the angularmomentum L ( ε ), the frequency ω r ( ε ) [Eq. (2.6)] and thecurvature K ( ε ) [Eq. (2.11)] can be expressed in terms ofthe simple powers of the scaled energy ε , ε = E /α +1 / . (B1)In particular, one can express these classical quantitiesthrough their values at ε = 1 ( E = 1), I i = I i (1) ε, L = L (1) ε, ω − r = ω − r (1) ε (2 − α ) / (2+ α ) ,K = K (1) /ε . (B2)Therefore, due to the scaling properties (2.2) and (B2),we need to calculate these classical dynamical quantitiesonly at one value of the energy ε = 1. For simplicity ofnotations, we shall omit the argument ε = 1 everywhere,if it is not lead to misunderstandings.The radial action I r ( L, E ) [Eq. (2.4)] can be ex-pressed explicitly in terms of the frequencies ω ϕ and ω r [Eq. (2.6)], and their ratio f ( L ) [Eq. (2.8)], I r = 2 αα + 2 ω − r − L f ( L ) . (B3)To prove this identity, we express Eq. (2.6) for ω − r interms of the determinant, ω − r = ∂ ( I r , L ) ∂ ( E, L ) = ∂I r ∂E − ∂I r ∂L ∂L∂E . (B4)Calculating directly the derivatives in this equation byusing Eq. (B2), one obtains the expression for ω − r (1).Solving then this equation with respect to I r (1), one ar-rives at Eq. (B3). Differentiating the identity (B3) termby term over L and using the definition for the ratio offrequencies f ( L ) [Eq. (2.8)], for the curvature (2.11) onefinally obtains K = − α ( α + 2) L ∂ω − r ∂L = − απ ( α + 2) L ∂T r ∂L , T r = 2 πω r . (B5)According to Eqs. (2.6) and (2.8) with the help ofRef. [25], ω − r is obviously simpler quantity to differenti-ate over L than f ( L ), ω − r = 12 π √ Z x max x min d x p Q ( x, L, α ) , (B6) Q ( x, L, α ) = (cid:16) − x α/ (cid:17) x − L / , (B7)and x = r . The turning points x min ( L, α ) and x max ( L, α ) are determined by the equation: Q ( x, L, α ) = 0 . (B8)7Thus, we may calculate ω − r and f ( L ), and then, useEqs. (B3) and (B5) for the radial action I r , and curvature K at the scaled energy ε = 1. Then, one obtains theirenergy dependence through the scaling equations (B1)and (B2), respectively. Appendix C: Full analytical classical dynamics forpowers 4 and 6
For the powers α = 4 and 6, the roots of function (B7),in particular, the turning points x min and x max can beobtained explicitly analytically. Therefore, one can findthe explicit analytical expressions for the key quantitiesof the classical dynamics for the POT, namely, the radialfrequency ω r [Eq. (2.6)] (or the radial period T r ), and thefrequency ratio f ( L ) [Eq. (2.8)] in terms of the ellipticintegrals from Ref. [25] (all in dimensionless units).For α = 4, one has the cubic polynomial equation Q ( x, L, α ) ≡ x − x − L / x min , x max and x ; given by the Car-dano formulas explicitly as functions of L in the physicalregion L ≤ L C , r < ≤ r min ≤ r max ; x q = r q . Forthe radial period T r [Eqs. (B5) and (2.6)], one obtainsthe analytical expression through these roots in terms ofthe complete elliptic integral F( π/ , κ ) of the first kind[18, 25], T r = 2 πω r = √ √ x max − x F (cid:16) π , κ (cid:17) , (C1)where κ = [( x max − x min ) / ( x max − x )] / . For the ra-tio frequencies f ( L ) [Eq. (2.8)], one finds f ( L ) = Lπ √ x max √ x max − x Π (cid:18) r max − r min r max , κ (cid:19) , (C2)where Π( n, κ ) is the complete elliptic integral of the 3rdkind [25].For α = 6, one has the polynomial equation of the4th power, Q ( x, L, ≡ x − x − L / x + ix and x − ix ,and again, two real positive roots, x min and x max ; seeEqs. (B7) and (B8)]. The radial period T r is determinedthrough these roots by the expression [similar to Eq. (C1],see Refs. [18, 25], T r = √ p AB ( x max − x ) F (cid:16) π , κ (cid:17) , (C3)where κ = { [( x max − x min ) − ( A − B ) ] / (4 AB ) } / , A = [( x max − x ) + x ] / , B = [( x min − x ) + x ] / .(We reduced the 4-power polynomial equation to a cubicone and obtained its 4 analytically given roots, mentionedabove, in the explicit Cardano’s form as functions of L ). For f ( L ) [Eq. (2.8)] at α = 6, one obtains [25] f ( L ) = √ L ( A + B ) π √ AB ( Ax min − Bx max ) × (cid:20) β F (cid:16) π κ (cid:17) + β − β − β ) Π (cid:18) π, β − β , κ (cid:19)(cid:21) , (C4)where Π( ϕ, n, κ ) is incomplete elliptic integral of the3rd kind, β = ( Ax min − Bx max ) / ( Ax min + Bx max ) , β = ( A − B ) / ( A + B ) . The curvatures K [Eq. (B5)]for α = 4 and 6 are determined by taking analyticallythe derivative of the radial period T r [Eqs. (C1) and(C3)] over L through the derivatives of the roots x min ( L ), x max ( L ), x ( L ) and x ( L ) for the derivative of F( π/ , κ )over κ [25]. The expressions for the curvatures K at theboth powers α = 4 and 6 can be found in the closedanalytical form through a rather bulky formulas, whichcontain the complete elliptic integrals of the 1st and 2ndkind. Appendix D: Classical dynamics and boundaries forthe diameters
For the primitive diameter D = (2 , S D (all in this appendix in dimensionless units) is specifiedanalytically through the scaled period τ D and energy ε by S D = τ D ε , τ D = 4 √ πα + 2 Γ (cid:18) α (cid:19) Γ (cid:18)
12 + 1 α (cid:19) , (D1)where Γ( x ) is the Gamma function of a real positive ar-gument x . For the diameter PO boundaries, one can usethe same L − = 0, but L + = b D L C , where b D = 1 −
12 exp " − (cid:18) L HOD − L D ∆ D (cid:19) (D2)(see Ref. [17] and more details in relation to the HO limitin Sec. III E), L D = 0 is the stationary point, L HOD = L C / ε/ (2 √
2) is the upper angular momentum L + for the D orbits in the limit α →
2, in which b D → /
2. In the semiclassical limit ε ≫
1, one has b D → ∆ D = ( πM n r K D ) − / is the Gaussian width of thetransition region between these two asymptotic limits.The D curvature for α ≥ L = L D is given by K D = Γ (1 − /α ) ε √ π Γ (1 / − /α ) . (D3)This exact analytical expression for the curvature K D atany α was derived by using a power expansion in Eqs.(B8) and (B3) over the variable proportional to L near L = 0 up to the terms linear in L . The Maslov phasefor the diameter orbit was determined by Eq. (3.14) at n r = 2 and n ϕ = 1. Note that for the limit α → τ D and action S D [Eq. (D1)], the Maslov index σ D [Eq. (3.14)] with thesame asymptotic (SSPM) limit of the constant part of thephase φ ( D ) d = − π (4 b D − / → − π/
4, and the curvature K D [Eq. (D3)] for the diameter (2 ,
1) are identical tothose obtained in Ref. [18].
Appendix E: The boundaries and curvature forcircle orbits
For the arguments Z ( ± ) p and Z ( ± ) r of the error functionsin Eq. (3.26), one originally has Z ( ± ) p,MC = r − i ~ J ( p ) MC (cid:16) p ( ± ) r − p ∗ r (cid:17) , Z ( ± ) r,MC = r − i ~ J ( r ) MC ( r ( ± ) − r C ) , (E1)where p ∗ r = 0 is the stationary point, p ( ± ) r and r ( ± ) are maximal and minimal classically accessible values of p r and r as the finite integration limits for the corre-sponding variables. To express the integration bound-aries (E1) in an invariant form through the curvature K C (3.29), and stability factor F MC (2.9), one may usenow the simple standard Jacobian transformations, andthe definition of the angle variable Θ ′ r as canonicallyconjugated one with respect to the radial action vari-able I r by means of the corresponding generating func-tion. In these transformations, we apply simple lin-ear relations: p ′′ r − p ∗ r = ( ∂p ′′ r /∂L ) ∗ ( L − L ∗ ) , and r ′ − r ∗ = ( ∂r ′ /∂Θ ′ r ) ∗ ( Θ ′ r − Θ ∗ r ), where we immediatelyrecognize the Jacobian coefficients. Note that there isno crossing terms due to the isolated stationary point I ∗ r = 0 , Θ ∗ r = 0 and to equations for the canonical trans-formations. At the stationary point for the isolated cir-cle PO, one has f ( L C ) = − ( ∂I r /∂L ) L = L C = − / √ α + 2[Eqs. (2.8), (A3) and (A7)]. For the transformation ofthe derivative ∂r ′′ /∂Θ ′ r , one can apply the Liouville con-servation of the phase space volume for the canonicalvariables to arrive at ∂r ′′ /∂Θ ′ r = ( ∂I r /∂L ) / ( ∂p ′ r /∂L )and |J CT ( p ′ r , p ′′ r ) | = | ( ∂p ′′ r /∂L ) / ( ∂p ′ r /∂L ) | = 1 at the POconditions r ′ → r ′′ → r C , p ′ r → p ′′ r →
0. Using also theJacobian identity, F MC = −J ( p ) MC J ( r ) MC / J MC ( p ′ r , p ′′ r ) , (E2)one obtains Eq. (3.28) for the arguments of the errorfunctions in Eq. (3.26).The expression (3.29) for the C curvature K C (in di-mensionless units at ε = 1) was obtained from expan-sion of f ( L ) [Eq. (2.8)] as function of L in powers of L C − L = ǫ up to the 2nd order terms in ǫ . For thispurpose, by using standard perturbation theory, we haveto solve first Eq. (B8) for the turning points r min and r max , [the integration limits in Eq. (2.8)] in the following general form ( r is taken below in units of R ), r max = r C + c ǫ + c ǫ + c ǫ + c ǫ + · · · ,r min = r C − c ǫ + c ǫ − c ǫ + c ǫ + · · · . (E3)Existence of such form of the solutions follows from asymmetry of the equation (B8) with respect to the changeof the sign of ǫ . Substituting these solutions into Eq. (B8)for arbitrarily small ǫ , one gets the system of the recur-rent equations for the coefficients c n . The solutions ofthis system up to the 4th order in a perturbation param-eter ǫ is given by c = r L C α , c = − α + 16 c , c = ( α − α + 5)72 c ,c = − ( α + 1)(4 α + 8 α + 13)1080 c , (E4)and so on. We transform now the integration variable r in the integral of Eq. (2.8) for f ( L ) to y , r = r C (1 − y ),such that f ( L ) = − L − ǫ π Z y max y min d y (1 − y ) p Q ( y, L, α ) . (E5)Here, Q ( y, L, α ) is given by Eq. (B7), Q ( y, L, α )= 2 r C (cid:20) − α + 2 (1 − y ) α − L C + 2 L C ǫ − ǫ (cid:21) ≡ ( y max − y )( y − y min ) R ( y ) , (E6) y max = ¯ c ǫ − ¯ c ǫ + ¯ c ǫ − ¯ c ǫ ,y min = − ¯ c ǫ − ¯ c ǫ − ¯ c ǫ − ¯ c ǫ , (E7)where ¯ c n = c n /r C . We use the last representation inEq. (E6), introducing a new function R ( y ) of the newvariable y to separate the singularities of the integrandin Eq. (E5) due to the turning points. This integrandhas to be integrated exactly by using a smooth function R ( y ) of y , which can be expanded in y at y = 0 up tothe second order, R ( y ) = R (0) + R ′ (0) y + 12 R ′′ (0) y + · · · . (E8)In order to get analytically the final result, we note that y in this expansion is of the order of ǫ , according toEq. (E7). Substituting then these expansions (E7) and(E8) into very right of Eq. (E6), we expand their middlein y at y = 0 up to the 4th order. After the cancellationof ǫ from both sides, and simple algebraic transforma-tions, one has R (0) = 2 L C c (cid:2) ǫ c (1 − k ) (cid:3) + O ( ǫ ) , k = c c , R ′ (0) R (0) = 2 k + O ( ǫ ) , R ′′ (0) R (0) = 2 k (3 k + 1) + O ( ǫ ) . (E9)9For the calculation of the circle orbit curvature K C ,we obviously need only quadratic terms in ǫ [linear in( L C − L )]. Therefore, one may neglect the ǫ correc-tions in the second and third lines of Eq. (E9) becausethey are multiplied by y ∼ ǫ and y ∼ ǫ in the ex-pansion (E8), respectively. Substituting now expansions(E7) and (E8) into the integral over y in Eq. (E5), andtaking R (0) off the integral, one then expands to thesecond order all quantities of the integrand in y ∼ ǫ , ex-cept for ( y max − y )( y − y min ) under the square root (inthe denominator) which can be integrated exactly. Tak-ing remaining integrals as R dyy n / p ( y max − y )( y − y min ) from y min to y max [Eq. (E7)], and then, expanding finally f ( L ) [Eq. (E5)] in ǫ , we find that the linear terms ex-actly disappear. It must be the case because f ( L ) is aneven function of ǫ . Thus, the coefficient in front of ǫ with the expressions for c n ( n = 1 , ,
3) from Eq. (E4)is Eq. (3.29) for the curvature K C . We can also usethis perturbation method for calculations of the next or-der curvatures, for instance, ∂ I r /∂L , which appears inexpansion of the phase integral in the exponent up tothe third order terms near the stationary points within amore precise (3rd-order) ISPM [9]. [1] V. M. Strutinsky, Nucl. Phys. A , 420 (1967); A ,1 (1968).[2] M. Brack, J. Damg˚ard, A. S. Jensen, et al. , Rev. Mod.Phys. , 320 (1972).[3] M. C. Gutzwiller, J. Math. Phys. , 343 (1971); Chaos in Classical and Quantum Mechanics (Springer,New York, 1990).[4] V. M. Strutinsky, Nukleonika (Poland) , 679 (1975);V. M. Strutinsky and A. G. Magner, Sov. J. Part. Nucl. , 138 (1976).[5] M. V. Berry and M. Tabor, Proc. Roy. Soc. London, Ser. A349 , 101 (1976);J. Phys. A , 371 (1977).[6] S. C. Creagh and R. G. Littlejohn, Phys. Rev. A , 836(1991); J. Phys. A , 1643 (1992).[7] V. M. Strutinsky, A. G. Magner, S. R. Ofengenden, andT. Døssing, Z. Phys. A , 269 (1977).[8] M. Brack, S. M. Reimann and M. Sieber, Phys. Rev. Lett, , 1817 (1997).[9] A. G. Magner, S. N. Fedotkin, K. Arita, and K. Mat-suyanagi, Prog. Theor. Phys. , 853 (2002).[10] A. G. Magner, I. S. Yatsyshyn, K. Arita, and M. Brack,Phys. Atom. Nucl., , 1445 (2011).[11] M. Brack and R. K. Bhaduri, Semiclassical Physics (re-vised edition: Westview Press, Boulder, USA, 2003).[12] R. D. Woods and D. S. Saxon, Phys. Rev. , 577 (1954). [13] K. Arita, Int. J. Mod. Phys. E