{'fraction_non_alphanumeric': 0.037978573038081, 'fraction_numerical': 0.010508634287236003, 'mean_word_length': 4.5811135246370185, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Quantum and semiclassical phase functions for the quantization of symmetric oscillators

We investigate symmetric oscillators, and in particular their quantization, by employing semiclassical and quantum phase functions introduced in the context of Liouville-Green transformations of the Schr\"{o}dinger equation. For anharmonic oscillators, first order semiclassical quantization is seldom accurate and the higher order expansions eventually break down given the asymptotic nature of the series. A quantum phase that allows in principle to retrieve the exact quantum mechanical quantization condition and wavefunctions is given along with an iterative scheme to compute it. The arbitrariness surrounding quantum phase functions is lifted by supplementing the phase with boundary conditions involving high order semiclassical expansions. This allows to extend the definition of oscillation numbers, that determine the quantization of the harmonic oscillator, to the anharmonic case. Several illustrations involving homogeneous as well as coupling constant dependant anharmonic oscillators are given.


I. INTRODUCTION
The harmonic oscillator is quite often employed as a prototype for illustrating numerous phenomena in different areas of quantum mechanics. It is somewhat infortunate, since the harmonic oscillator is all but atypical, even within the family of symmetric oscillators. In particular the quantization of symmetric oscillators is radically different in the harmonic and anharmonic cases: as is very well-known, the first order semiclassical quantization condition, given by S(t 2 , E) − S(t 1 , E) = π (n + 1 2 ) (1) where S is the classical action, t 1 and t 2 are the turning points and n is the level integer for the energy E, is exact for the harmonic oscillator, but fails to capture even the most elementary aspects in the anharmonic case. For example, for homogeneous potentials x 2m (m = 1), Eq. (1) predicts that the ground state energy should decrease when m increases, whereas the exact quantum ground state levels behave the other way round (E increases with m). It is therefore no surprise that considerable theoretical and computational efforts have been made to investigate these systems. On the one hand, symmetric oscillators are among the simplest nonsolvable systems [1]. On the other hand they are employed in many branches of quantum physics, ranging from molecular vibrations to simple models of quantum field theories. The quantization of symmetric oscillators has largely focused on the accuracy of the standard semiclassical expansion. Eq. (1) can usually be improved by going to higher order in , as will be briefly recalled in Sec. 2. But the asymptotic series is factorially divergent, so at best the expansion must be truncated, although the accuracy rapidly increases with n [2]. The possibility of employing resummation techniques and resurgence analysis was also extensively investigated, especially on the quartic oscillator [3]. Proofs regarding the convergence of resurgence schemes do exist in some particular cases [4], though it seems hardly realistic to base practical calculations on such schemes. Alternatively, it was pointed out [5] that symmetric oscillators have turning points away from the real line and that including the subdominant contributions from these complex trajectories (where space-time coordinates become complex) accounts for most of the discrepancies between the best semiclassical estimate and the exact quantum result. It was even recently suggested that having particles moving along trajectories in the complex plane should be regarded as a real physical feature [6]. We also mention the importance of anharmonic oscillators of the form x 2 + λx 2m in the investigation of quantum perturbation theories [7].
Thus, even for simple systems such as symmetric oscillators with a single minimum, quantization is far from being understood. In the present work we analyze the quantization of symmetric oscillators (determination of bound states in potential wells with a single minimum) from the perspective of an exact quantum phase. The phase -a real function -is obtained from a Liouville-Green transformation of the Schrodinger equation. It achieves exact quantization and allows to retrieve the exact quantum mechanical wavefunction. It also extends to a general symmetric oscillator certain quantities that only exist for the harmonic oscillator, such as the oscillation number. Quantum phase functions suffer from ambiguities in their definitions; the idea to be implemented here is to obtain a phase function that will be close, in a sense to be precised below, to the predivergent semiclassical expansion. We will also give a simple and efficient numerical procedure to construct the phase function, achieving quantization with an a-priori unlimited precision. We will first recall the main semiclassical quantization schemes (Sec. 2). The standard semiclassical approximation is the most common one, but it is well-known (though too often overlooked) that expansions can be carried out by employing alternative schemes, based on different transformation functions. We will then introduce the quantum phase by making a specific transformation that will be justified, in particular by its relation to the semiclassical phase (Sec. 3). In Sec. 4 we will illustrate the roles of quantum and semiclassical phases in the harmonic oscillator case: this is a solvable problem that will further allow us to gain insight in the meaning of the quantum phase by working with the analytical solutions. We will in particular introduce the notion of optimal quantum and semiclassical phase functions by linking them to the exact and semiclassical oscillation number. After giving the numerical procedure to compute the quantum phase (Sec. 5), we will deal with the quantum phase investigation of anharmonic symmetric oscillators in Sec. 6. We will first insist on the pure quartic oscillator, which is the one that has been the most heavily studied. We will also illustrate some properties of sextic and octic oscillators as well as the behaviour as a function of λ for a harmonic oscillator perturbed by a λx 10 anharmonicity. Our closing comments and our conclusions will be given in Sec. 7.

A. Standard complex plane quantization
The standard semiclassical approach in one dimension -the WKB scheme -is based on a Riccati transform of the wavefunction followed by an expansion. The first step is to write the wavefunction ψ(x) as thereby transforming the Schrodinger equation for ψ, into a Riccati-type equation for the phase ζ. By the argument principle, the logarithmic derivative of ψ integrated along a contour encircling the n zeros of ψ yields the quantization condition (see eg [8]) Substituting the expansion into the Riccati equation and equating like powers of yields ζ 0 (x) = S(x) and the recurrence relation Eq. (4) was obtained by Dunham [9]. The solutions of Eq. (6) are then plugged into the asymptotic expansion (5) which is substitued in Eq. (4). Note that the contour in Eq. (4) must now enclose the turning points [10] since the solutions ζ k are singular at the turning points. The first odd term ∂ z ζ 1 (z) integrated along the contour yields iπ whereas all the other odd terms are total derivatives and therefore do not contribute [2]. Eq. (4) thus takes the form Except in specific cases (such as the quartic oscillator for which the integrals are known analytically), the integrals in Eq. (7) are taken on the real line and must be regularized since the functions ζ ′ 2k contain nonintegrable singularities at the turning points [11]. Finally inverting Eq. (7) for a finite value of k max yields the quantized eigenvalue E. As for any asymptotic series [12] the general trend as k max is increased is to obtain better approximations for the first few terms but quickly the series diverge (examples will be given below).

B. Alternative asymptotic quantization schemes
More general asymptotic expansions are readily obtained by employing alternatives to the transformation (2). Different asymptotic ansatzes can be found in Refs. [8,13]. The most useful form is based on the Liouville-Green transformation whereby the wavefunction is written as where u and w are two arbitrary (but sufficiently smooth) functions and ξ appears as a new dependent variable. This transformation can be restricted by requiring two linearly independent solutions of the Schrödinger equation ψ 1 and ψ 2 to have exactly the same form (8) with two different functions w 1 and w 2 . Recalling that the Wronskian W[ψ 1 ,ψ 2 ] is a constant, we are led to the transformation Assume w(ξ) fulfills the equation where the choice of the unspecified but smooth function R(ξ) determines the choice of w. ξ then obeys The procedure is now to employ an asymptotic expansion for ξ, which is substituted into Eq. (11) to obtain a recurrence system with the first term being found from Note that the odd terms in the expansion (12) are redundant. The advantage of the present scheme relative to the standard semiclassical treatment recalled above is twofold. First the choice of R, which determines both w(x) and ξ(x) allows a greater flexibility. In particular it allows to construct wavefunctions that are well defined on the entire real line, so that although the preceding analysis is valid for both complex and real variables, it is possible to work with real quantities only. This is of course not necessary, since the so-called 'phase integral method' [15] basically amounts to choosing the exponential function for w but keeping complex values for ξ, therefore differing from the standard method by requiring Eq. (9) to be enforced at each order. But as we shall argue below, working with real functions brings in a second advantage, which is of conceptual nature, namely that it allows to better understand the classical limit, or conversely to better follow what the classical quantities becomes in the quantum domain.

III. QUANTUM AND SEMICLASSICAL PHASE FUNCTIONS
A. Choice of the phase: carrier functions It is apparent that putting the wavefunction in the form given by Eq. (9) is tantamount to undertaking an amplitude-phase decomposition, where ξ(x) appears as a phase and appears as an amplitude function. Indeed Eq. (14) simply represents the one-dimensional continuity equation. w represents here a 'carrier' function, as it carries the phase. It is clear that defining ξ(x) as a quantum phase is ambiguous: first, the carrier function needs to be specified (thereby deducing R(ξ)) and second the boundary conditions of the third order nonlinear equation (11) need to be given. From a purely internal quantum mechanical viewpoint, it would appear at first sight that the choice of ξ(x) is meaningless, insofar as the quantum mechanical quantities are insensitive to a specific manner of cutting the wavefunction. From a semiclassical point of view things are different, since different choices of ξ(x) yield a different first order semiclassical phase ξ 0 (x). From Eq. (13) it is straightforward to see that the only choice leading to ξ 0 (x) = |S(x)| corresponds to |R(ξ 0 )| = 1, corresponding to exponential or circular carrier functions. Note that in that case, the terms ξ 2k (x) are identical in absolute value to the even terms of the expansion (5), thereby leading to the same quantization condition; the odd terms are different, although in both cases they are singular at the turning points. Of course this is not the case if Eq. (11) is to be solved exactly, since the exact solutions are well defined on the entire real line As an example of a semiclassical phase different from the classical action, consider taking R(ξ) = ξ, corresponding to Airy carrier functions. This case has a practical interest because it yields an asymptotic expansion not diverging at the turning points. Given the symmetry of the oscillators, it is convenient to restrict the analysis to the half line, say ] − ∞, 0]. The first order solution can be written in the handy form where the signs change when x crosses t 1 due to the Stokes phenomenon. Eq. (15) appears in the context of comparison equations when p(x) admits one single turning point; see eg Ref. [8], in particular Sec. II.3, where it is shown that ξ 0 (x) as well as the higher-order solutions ξ 2k (x) are well-behaved in the vicinity of t 1 . From a purely semiclassical perspective it is often advantageous to employ in the same problem different comparison equations to expand the solution in the neighbourhood where the relevant approximation holds best [14], leading to the matching of different phase functions. Whether ξ(x) can be called a 'phase' or not is a question of terminology. But it is clear that the choice of R(ξ) gives different carrier functions and different asymptotic expansions for ξ. What is interesting in the present context is that the phase corresponding to the standard semiclassical choice R(ξ) = 1 can be obtained in terms of other phase functions, opening the possibility of performing any asymptotic expansion for the standard semiclassical phase. This is done by recasting the phase equation in terms of Ermakov systems.

B. Choice of the phase: boundary conditions and Ermakov systems
Uncoupled Ermakov systems relate solutions of the linear differential Eq. (3) and of a nonlinear equation for α(x) similar to Eq. (11) when |R(ξ)| = 1, Eq. (16) was employed very early in a quantum mechanical context [16], but the connection with Ermakov systems, based on the nonlinear superposition principle [17], is quite recent (see [18,19,20] and Refs. therein). The only results we shall need are the following. Assume R(ξ) = 1 (ie, sin or cos carrier functions) and denote the phase in this case by σ(x) ≡ ξ(x). Let ψ 1 (x) be a solution of Eq. (3) regular at −∞ and ψ 2 (x) be a solution regular at +∞ (given the symmetry, we can take here ψ 2 (x) = ψ 1 (−x)) and W = W[ψ 1 , ψ 2 ] their Wronskian. Eq. (9) can be written as where I is the Ermakov invariant. Using σ(−x) = σ(∞) − σ(x) and σ(0) = σ(∞)/2 we also have the identity A solution of Eq. (3) independent from ψ 1 (x) and lagging π/2 out of phase is given by where c is an arbitrary constant. ψ 3 is in general irregular at ±∞. From these considerations, it follows that We have emphasized the dependence of σ(x) on I and c, the form taken in the present context by the two boundary conditions left from Eq. (11) once σ(−∞) = 0 is imposed as is done in Eq. (20). Actually once normalization is imposed (improper normalization except at the energy eigenvalues), there is only one free parameter left. Finally the quantization condition is obtained by noting that α(x → ±∞) → ∞ (because from Eqs. (20) and (14) we see that α 2 (x) is proportional to ψ 2 1 + ψ 2 3 ). It then follows from Eq. (17) that we must have σ(+∞) = (n + 1)π if ψ 1 is to be an eigenfunction and given the symmetry of the potential, the quantization condition thus reads where n is the integer counting the zeros of ψ 1 on the real line. Note that when E is not an eigenvalue, Eq. (21) is replaced by With regard to the semiclassical limit, we can draw two consequences from Eq. (20). First, by employing directly the asymptotic expansions for ψ 1 and ψ 3 , it can be seen that the semiclassical phase σ sc (x) is different from the classical action, except for a single value of c (see Appendix A). Therefore the classical action S(x), that appears as the first order solution of Eq. (11) for R(ξ) = 1 is indeed a solution of the semiclassical limit of the quantum phase, but one among others. Second, by employing Eq. (9) in Eq. (20), we obtain where w 1 and w 3 are independent solutions of Eq. (9) lagging π/2 out of phase. Eq. (23) expresses the phase obtained from Eq. (11) with R(ξ) = 1 in terms of solutions ξ(x) obeying the phase equation with a different function R(ξ). Now by replacing ξ(x) by its asymptotic expansion, we obtain an asymptotic expansion for σ(x). For example taking R(ξ) = ξ, a solution regular at −∞ is the Airy function Ai(x) and an irregular solution lagging π/2 out of phase is Bi(x), so that to first order we have where ξ 0 (x) is given by Eq. (15). Eq. (24) gives a uniform semiclassical phase σ sc (ie, with R(ξ) = 1) free of singularities at the turning points, with σ sc (−∞) = 0, σ sc (t 1 ) = π/6 and the value σ sc (0) defines the quantization condition when where the use of the relevant branch of the arctan is implicitly understood. Note that by using the expansions of the Airy functions for large negative values of the argument, for the regular solution and the relevant expansion for Bi, Eq. (24) takes the familiar approximate form Thus for large values of the action (so in practice for sufficiently large values of E) the semiclassical phase obtained from appropriate Airy carrier functions is approximately the same, in the classically allowed region, than the one obtained from the standard quantization scheme recalled above. However at low energies the two different expansions of σ sc lead to markedly different phase functions.

C. Conclusion
Let us recapitulate. First we have recalled that the classical action is not the only semiclassical phase, nor necessarily the most useful one in practical computations. Certainly, the choice R(ξ) = 1 is the most natural one, since it is the only one that allows to recover the classical action, but we have seen that even in that case it is possible and advantageous to express the relevant phase σ(x) in terms of alternative phase functions ξ(x), leading to different semiclassical expansions σ sc (x). Second, quantum phase functions suffer from ambiguities, because they are irrelevant as quantum mechanical objects. Even if R(ξ) = 1 is chosen, there is still a free parameter (in the form of a boundary condition) that can be varied at will leading to different behaviour of the phase function. The phase ambiguity persists in the semiclassical limit: to first order in hbar there is a single value of the boundary condition leading to purely classical quantities -precisely the value that eliminates remnants of oscillating quantum quantities from the semiclassical phase (Appendix A). Now the question can be reversed. Assume that we have chosen the 'right' semiclassical phase, that is the one that only involves classical quantities. Can the corresponding quantum phase be constructed? This would amount to implicitly sum the diverging asymptotic expansions (5)- (7). The answer is positive in the case of the harmonic oscillator, due to the existence of analytic solutions. For other symmetric oscillators, it does not seem possible to determine this optimal solution exactly. Rather, Eq. (11) can be solved (numerically) with a boundary condition allowing to approximate this optimal solution. Notwithstanding exact quantization can be achieved and the exact wavefunctions can be retrieved as well. These points are developped in the subsequent paragraphs of the paper.

IV. A SPECIAL CASE: THE HARMONIC OSCILLATOR
The interest of the harmonic oscillator in the present context arises from its solvability: analytic solutions can be explicitly written down and the energy obtained as a function of the oscillation number in closed form. Though the approach given above is not necessarily useful for solving the harmonic oscillator problem as such, the harmonic oscillator represents a system for which the notions introduced in Secs. 2 and 3 can be illustrated on firm grounds. This study will in turn be valuable for understanding the approximate and numerical treatments that will be undertaken for general anharmonic oscillators.

A. Analytic solutions and the oscillation number N (E)
Let us solve the Schrödinger Eq. (3) with where we have written so that Eq. (3) takes the form of the Weber equation. A solution regular at −∞ is given by the parabolic cylinder function which can be written down in the integral representation as with the contour encircling the negative real axis and the factors ensure normalization for the eigenfunctions. A linearly independent solution regular at +∞ is D ν (−z) and their Wronskian is seen to be where we have used Γ(µ)Γ(1 − µ) = sin πµ/π. Note that the knowledge of Eq. (29) in conjunction with Eq. (31) suffices to determine the eigenvalues E, since the eigenfunctions are regular at both ±∞ and the Wronskian must therefore vanish. By doing so Eq. (29) is not interpreted as a simple change of variable, but as a functional relation by which the energy is given in terms of an oscillation number, which gives the number of oscillations of D v (z) between x = −∞ and x = ∞. We shall take for granted that ν indeed captures the entire oscillatory character of the parabolic cylinder functions; proofs may be obtained by employing an asymptotic expansion for D ν (z) in the interval [t 1 , t 2 ] [21] or by following Olver in introducing an auxiliary modulus function whose monotonicity on [t 1 , 0] may be proven [22]. Note that the first order semiclassical quantization condition (1) allows to define a first-order semiclassical oscillation number N sc (E) by the relation and the semiclassical quantization condition takes the form N sc (E) = n + 1. For the harmonic oscillator the action difference is Eπ, hence N sc (E) = E + 1/2, which by Eq. (29) yields the exact quantum relation (32).
B. Quantum and semiclassical phase functions

Quantum phases
We now examine, from the point of view of the formalism introduced in section 3, the choice of the quantum phase. We have already seen that working with circular carrier functions, R(ξ) = 1, is the most advantageous choice, so our problem is to find the most relevant phase function σ(x) in Eq. (17). Of course ψ 1 (by definition regular at −∞; recall we have imposed σ(−∞) = 0) is necessarily proportional to D ν , but there are infinite ways of decomposing ψ 1 as in Eq. (17). The main interest in using the amplitude-phase decomposition is probably that the oscillatory character of the wavefunction is entirely captured in sin σ(x), so that the amplitude does not oscillate. The unique quantum phase function displaying such a behaviour will be termed 'optimal'. In this case, we should have (or σ opt (0) = πN (E)/2 given the symmetry). We have just seen that the oscillation number is given in the harmonic oscillator case by N (E) = ν + 1. Therefore comparing with Eq. (22) gives c = − cot πν/2I and comparing Eq. (31) with Eq. (18) allows to set I = π −1 . These values of the parameters I and c ensure that the quantum phase function given by (20) is the optimal one. Note that we must allow for scale transformations ψ 1 → κψ 1 and ψ 2 → κψ 2 , implying W → κ 2 W, I → κ 2 I, c → c/κ 2 so that the parameters corresponding to the optimal phase actually belong to a class by which Ic and W/I are given in terms of the oscillation number,

Semiclassical phases
The semiclassical phase, defined as the semiclassical limit of σ(x) is given by Eq. (A3) in Appendix A. For arbitrary values of I and c, it is straightforward to see that ∂ x σ sc (x) will be a highly oscillating function of x. The optimal semiclassical phase is the classical action, giving a non-oscillating function ∂ x σ sc (x) = p(x), obtained for Note that here the set {Ic, W/I} refers to the semiclassical Ermakov system, and is generally different from the quantum set of parameters; this is why Eq. (38) is verified by construction. The very special property of the harmonic oscillator is that the right handside of Eqs. (37)-(38) is equal to the right handside of Eqs. (35)-(36), since to lowest order in we have φ = π/4 and S(t 2 ) is immediately evaluated as Eπ = π(ν + 1/2). This property explains why semiclassical quantization is exact, the mapping {Ic, W/I} sc → {Ic, W/I} quantum being the identity.

V. NUMERICAL DETERMINATION OF THE QUANTUM PHASE
We present in this section an efficient numerical method allowing to compute the desired quantum phase σ(x). This is tantamount to solving directly Eq. (11) with R(ξ) = 1. From a computational point of view, it is advantageous to employ a function akin to the Riccati transformed ζ in Eq. (2) but incorporating from the start the amplitude-phase decomposition. Indeed, defining Eq. (11) becomes equivalent to which is a complex but first order nonlinear differential equation. M is found by an iterative linearization procedure. We introduce the functional and linearize Eq. (40) by expanding F to first order in the vicinity of an initial trial function M 0 (x). The resulting first-order linear differential equation is solved for M 1 (x), and the process is iterated until convergence is achieved after q iterations (details are given in Appendix 2). The exact form of the initial trial function is unimportant (but see Appendix 2) provided it is smooth and behaves as the converged solution. M 0 (x) can be conveniently built from the uniform semiclassical approximation to σ sc (x) given by Eq. (24), since the converged function M q (x) will be close to the trial function. The delicate and important point is to determine the boundary condition on Eq. (40), which holds for all the functions M i . It is convenient to choose the boundary condition at x = 0: symmetry imposes ∂ x α(x = 0) = 0 so only the real boundary condition ∂ x σ(x = 0) needs to be set. Let us take again the harmonic oscillator as a model. The knowledge of the analytic solutions and of the oscillation number allow us to determine the value ∂ x σ(x = 0) for what we called above the optimal quantum phase. Taking ψ 1 (x) = D ν (x), we have from Eqs. (17), (22), (32) and (35) where [22] D ν (0) = 2 ν/2 π 1/4 Elementary manipulations on the Γ functions yield Employing the boundary condition (44) in the iterative scheme yields a converged solution M ; integrating the real part of this solution gives us the optimal quantum phase having the property σ(∞) = πN (E) and whose derivatives ∂ m x σ are non-oscillating functions. If the boundary condition is (slightly) different, the quantum phase will also (slightly) differ from the optimal one: the phase at infinity will (slightly) differ from the oscillation number and (slight) oscillations will appear in the first derivatives (the oscillations will only become prominent for high order derivatives). Recall however that the wavefunctions (17) and (19) are exact solutions of the Schrödinger equation irrespective of the oscillations of the phase-derivative. Note also that the present numerical scheme does not allow to construct highly oscillating functions, given that the trial function is itself not oscillating.
Quantized energies are found by computing σ(x = ∞, E) = 2σ(x = 0, E) as a function of E, which in practice means determining σ(x = ∞, E) on a rather loose energy grid and then interpolate. It is then possible to solve for the eigenenergies E n by employing the quantization condition σ(∞, E) = π(n + 1), which in principle holds irrespective of the boundary condition. In practice it is however important to employ adequate boundary conditions so as to keep the interpolation tractable, which is possible provided σ(∞, E) is well behaved. Note that the accuracy of the energy eigenvalues can be improved by tightening the energy grid in the vicinity of each E n ; we have generally obtained eigenvalues with 24 decimal digits without much numerical effort.

A. General setting
Anharmonic symmetric oscillators are in general not solvable, except in some exceptional cases (e.g. some states in shape-invariant potentials in supersymmetric quantum mechanics [23]). There are no closed-from solutions to the Schrödinger equation and no such thing as the oscillation number, or more generally no such thing as a functional relation E = f (n * ) by which the spectrum is obtained as a function of a smoothly varying good quantum number n * . Semiclassically, such a relation exists by construction, by inverting ζ sc (E), ξ sc (E) or σ sc (E) (the absence of the space variable means that the phase has been integrated on the relevant contour or real line) and quantization precisely occurs when n * is an integer. But we know that all these semiclassical phase functions are diverging asymptotic . The quantum phase σ(x) (solid black curve) is plotted at the exactly quantized ground state energy. The classical action S(x) + π/4 is shown between the turning points at the WKB quantized energy (solid thick blue curve) and at the exact eigen-energy (thick dash-dotted red curve). The π/4 term is well-known to be obtained through connection formulae or via an asymptotic (in x) expansion of the uniform semiclassical approximation for σ, see Eq. (27). The first order uniform semiclassical phase obtained with Airy carrier-functions is given by the black dashed curve, seen to follow closely the quantum phase.
series, and thus the reciprocal relations (ζ sc ) −1 (n * ) are only approximate. Employing a quantum phase represents a compromise. It achieves exact quantization and allows to define relations of the type σ(E) or σ −1 (n * ) that are exact but arbitrary (except when n is an integer). As developed in the preceding sections the treatment is exact and satisfactory from a quantum-mechanical viewpoint: eigenfunctions and eigenvalues are determined with a much lower computational cost than the standard matrix methods, that need to employ large basis. Oscillation number functions N (E) can be constructed provided that the phase functions and their derivatives are sufficiently smooth; these oscillation number functions are approximate and not unique -they can be made as nearly exact as desired by optimizing the boundary condition on M , but an exact boundary condition of the type given by (44) in the harmonic oscillator case does not exist. The quantum phase employed in the results given below is defined from what is probably the most intuitive way of setting up the boundary condition: we determine the value at x = 0 of the local asymptotic expansion of the phase derivative (i.e. the inverse of amplitudes), see Eq. (12) with R(ξ) = 1. Solving the corresponding recurrence system directly gives us the elements ∂ x σ 2k (x). k max is set in principle by going to the highest possible order before the asymptotic expansion starts to diverge (so k max varies with the energy) and we use Stieljes' simple trick of terminating the series by multiplying the last retained term by 1/2. Indeed, the spirit of the present approach is obtain a quantum phase that would be 'close' in the classically allowed region to the semiclassical series if the latter converged. Obviously at some point (for large k max ) the extra effort imposed by the computation of high order terms ∂ x σ 2k (x) is not worth what is gained by including this term. The same can be said about using super-asymptotic or hyper-asymptotic methods [24]: getting into involved calculations that would result, say in changing the 20th decimal number, while not solving in principle the problem of achieving an optimal boundary condition is probably not advisable (given, to repeat, that the quantum phase is exact in all cases).

The homogeneous quartic oscillator
The pure quartic oscillator, with the classical momentum function given by is the simplest nonsolvable potential, and as such it has been the object of a great number of works that would be impossible to cite or summarize. Most of the work however still focuses on the same topics that motivated the  early papers of Bender et al. and Voros, namely the determination of the numerical properties of the semiclassical expansion for large orders [2] and the development of resummation procedures for the divergent series, in particular the understanding of analycities that produce the phenomenon of resurgence [3]. Our goal here is to give a few numerical results so as to show the relevance of the quantum phase approach. Fig. 1 shows the quantum phase for the ground state energy E n=0 = 0.53018... We also show the WKB phase (the classical action) calculated at the exact energy and the same phase determined at the WKB quantized energy E sc n=0 = 0.434... (which is off by almost 20%) and finally the first order semiclassical phase obtained with Airy carrier functions given by Eq. (24) at the corresponding quantized energy [cf. Eq. (25)] E sc,R=ξ n=0 = 0.480... Fig. 2 shows the derivatives of these different phases. The following comments can be made. The main observation is that quantum effects are important: the quantum phase keeps accumulating well beyond the turning points, in the classically forbidden region. The same feature is visible for the phase derivative: whereas the classical momentum vanishes at the turning points, the quantum phase is far from being negligible even at twice the value of the turning point. As required, the WKB phase is closer to the quantum one at the WKB quantized energy, but the phase derivative follows more closely the quantum curve when taken at the exact quantization energy. Note that the semiclassical phase obtained with Airy carrier functions is well behaved through the entire real line and thus follows the quantum phase in the classically forbidden regions; the quantum effects in this case show up in the dip visible in the phase derivative around x = 0. Note also that this semiclassical quantization is twice as accurate for the ground than the first order standard (WKB) semiclassical quantization. The standard semiclassical quantization scheme can be taken to higher order: divergence occurs after the third term, so we can go up to k max = 3. In this case, using our simple rule stated above for the terminant, we find the energy for the ground state to order 8 to be E sc n=0 = 0.483..., though stopping the expansion after the second term gives the slightly better result E sc n=0 = 0.490.... Fig. 3 shows the situation at a higher energy, for the 8th excited state with the exact energy E n=8 found by solving the quantum phase quantization condition. We only show the plot for the phase derivatives because the different phases would barely be distinguishable on the scale of the plot. The first order standard semiclassical phase derivatives taken at the exact and at the WKB quantized energies can not be distinguished on the figure, whereas the quantum phase In the classically allowed region the different curves are very close one to the other and appear as superposed. This trend is of course expected, given that for symmetric oscillators as E increases the first-order semiclassical approximation improves (in relative terms). Fig. 4 (a) gives σ(x = ∞, E) = 2σ(x = 0, E) interpolated as a function of the energy. Fig. 4 (b) zooms in the lower energy region, the boxes represent the quantized energies. The interest of such a curve is two fold. First, as already mentioned, we use this interpolated function to quantize the system with a high numerical precision and a modest computational cost. Second, this curve defines an approximation to the oscillation number N (E), introduced above in the context of the harmonic oscillator. Indeed, N (E) is readily obtained in the case of solvable potentials, but it is not defined for nonsolvable problems. The present scheme thus allows to define an approximate oscillation number to be denotedÑ (E), counting the number of half-wavelengths of the wavefunction.

Sextic and octic oscillators
Pure sextic and octic oscillators, ie V (x) = x 2m /2 with m = 3 and 4, have received little attention compared to the homogeneous quartic case. Indeed, except in the complex trajectory approach, where the number of turning points in the complex plane rapidly increases with m [6], all the homogeneous potentials of higher degree have the same basic properties that can be seen on the quartic oscillator. In particular the functions appearing in the standard semiclassical expansion (24) can be obtained in closed form. As m increases, the standard semiclassical quantization procedure becomes worse (except in the m → ∞ limit, where a solvable potential -the infinite square well is obtained but then (7) does not apply). For example for the ground state of the octic oscillator, standard semiclassical quantization breaks down after the crude WKB term and is 38% too low. The semiclassical quantization condition obtained with Airy carrier functions is 30% too low, and the dip around the origin is more severe than in the quartic case, as portrayed in Fig. 5. The oscillation numbersÑ (E) for the pure sextic and octic oscillators are displayed in Fig. 4. These approximate oscillation number functions were obtained by employing the boundary condition given above. As noted earlier, if a different boundary condition is set, the resulting phase function (and thusÑ (E)) will be different. This is illustrated in Fig. 6 for the amplitude function α = (∂ x σ) −1/2 of the sextic oscillator near the 4th excited level at E = 10.8571: we have plotted the amplitude functions α(x) and the derivatives ∂ 5 x α and ∂ 6 x α corresponding to two close but different boundary conditions. The left panel is obtained with the boundary condition (45) given above up to o(h 16 ) whereas the right panel corresponds to the 'WKB' boundary condition ∂ x σ(x)| x=0 = p(0). These two boundary conditions are very close, less than 3 parts in 10 −5 , and α(x) is indeed seen to be almost identical in both cases. However when higher derivatives of the amplitude are plotted, differences appear: the oscillatory structure is already visible in the fifth derivative of the p(0) amplitude function, leading to radically different sixth derivatives. This means that when the classical boundary condition is employed, the quantum amplitude retains part of the oscillatory character of the wavefunction.

Anharmonic perturbations
By far, the main interest in anharmonic oscillators has focused on solving problems in which anharmonic corrections must be added to the harmonic case. The simple Hamiltonian in which a quartic perturbation is added to the harmonic oscillator has given rise to an incredible number of works since the pioneering paper by Bender and Wu [7]. The two issues here are the aymptotic semiclassical series, as in the homogeneous case, but also the behaviour of the eigenvalues E(λ) as a function of the coupling parameter λ. The resulting perturbative series is divergent, hence the development of different resummation schemes based on scaling transformations or renormalization (eg [25]). The summation of the perturbative series for a decadic perturbation, i.e.
is considered to be particularly challenging [26]. Although the method presented in this paper could be advantageously put to work to study the quantum/classical correspondence of the perturbative series, and in particular the behaviour of the different semiclassical expansions introduced above as a function of λ, this topic is out of the scope of this work. We only wish to stress that the method based on the quantum phase is a powerful tool to investigate the behaviour of the energies and wavefunctions: rather than having recourse to different perturbative expansions depending on whether λ lies in the weak or strong coupling regimes, the quantum phase allows to compute essentially exact quantum results within a unified scheme. It behaves as a semiclassical phase function (in the sense that physical properties are extracted from the semiclassical and quantum phase functions in the same way) in situations in which the semiclassical approximation fails. This is illustrated in Fig. 7, where we display the oscillation numberÑ (E, λ) for the decadic perturbation (48) for different values of λ ranging from λ = 0.01 to λ = 1000. Fig. 8 zooms on the lowest states; the broken lines correspond to the semiclassical oscillation number N sc (E, λ) obtained from the classical action [cf Eq. (33)] for the same values of λ. The difference between the quantum and semiclassical curves gives a measure of the accuracy of the first-order semiclassical quantization, which as expected gets worse as λ increases and as E decreases.
In the strong coupling regime the semiclassical phase function can barely said to constitute an approximation: for example the predicted semiclassical energy for the ground state of the λ = 1000 oscillator is seen to lie between the exact ground state energies of the λ = 5 and λ = 50 oscillators; the slopes of the semiclassical and quantum phase functions, from which periodic time scales are determined, are markedly different. Finally Fig. 9 gives the behaviour of the oscillator's energy as a function of λ andÑ . By fixingÑ at an integer value n the figure gives the quantization energy of the nth level as a function of λ. Alternatively by fixing λ we can follow the energy of the system as the oscillation number increases.

VII. CONCLUSION
In summary we have investigated symmetric oscillators, and in particular their quantization, with the help of semiclassical and quantum phase functions. Although different semiclassical phase functions can be defined, the most natural choice corresponds to the scheme that yields to first order the classical action S(x). When the classical action and its associated amplitude (∂ x S) −1/2 are plugged into the semiclassical wavefunction, the entire oscillatory character of the wavefunction is captured by the oscillations of the phase. We have introduced quantum phase functions that mimic this property: an arbitrary quantum phase function would not be useful, but when supplemented by semiclassical boundary conditions it behaves in a similar way as the semiclassical phase, allowing to define oscillation number functions and to retrieve the exact quantum mechanical wavefunctions and eigenvalues. However the construction in principle of an 'optimal' phase function, the one that would generate entirely the oscillations of the wavefunction remained out of reach: this was possible for the harmonic oscillator case, thanks to the existence of closed form solutions, but not for the more general nonsolvable oscillators. The retained solution was to match the quantum phase to high order pre-divergent semiclassical expansions in the classically allowed region, so as to obtain nearly optimal phase functions.
Indeed, our underlying working hypothesis has been that the optimal quantum phase represents implicitly the resumed divergent series constituting the semiclassical expansion. It would therefore be interesting to connect the present approach with the so-called exact WKB analysis [3,27]. The quantizing equation for oscillators in exact WKB analysis arises from a quantum relation, namely the Wronskian of the function called ψ 2 (x) in this paper (appropriately normalized however) and of ψ 2 (q) where q represents the rotation of the real line in the complex plane by a spectral symmetry angle [27]. The Wronskian is determined at x = ∞ from an exact WKB representation of ψ 2  and at x = 0 in terms of spectral determinants. In principle a system of equations that can be solved numerically can be extracted by equating the expressions for the Wronskian, but expect for a few particular cases, the exact WKB analysis as well as methods based on resurgent functions [4] seem to be better qualified to obtain general proofs rather than numerical results. The method developed in this paper is some extent the opposite: the QLM method employed to determine the quantum phase is numerically straightforward and transparent, but the relation to semiclassics is indirect. The WKB type form of the initial trial function ensures that the converged result will not display strong oscillations and the main semiclassical input enters through the boundary condition. The oscillation number we define is a rewriting of the Wronskian W of the real solutions ψ 2 and ψ 1 (resp. recessive and dominant at ∞) appropriately renormalised by a factor that takes here the form of the invariant I. The present impossibility of ascribing a unique value to W/I explains why the oscillation number, i.e. the total phase accumulated on the real line, retains some arbitrariness (except at the eigenvalues). Lifting this arbitrariness should precisely be equivalent to giving a value to the divergent semiclassical series for the phase. and using the symmetry of the potential with W[ψ 1 , ψ 2 ] = 2I sin (S(t 2 ) + 2φ) . These representations of the wavefunctions are valid far from the turning points; we have taken S(t 1 ) = 0 and the connection rules at the turning point t 1 impose φ = π/4 (although this is not important in this context). These quantities are substituted into Eq. (19) to obtain ψ sc 3 (x; I, c) from which the semiclassical phase σ sc (x) is found as cot σ sc (x; I, c) = cot (S(x) + φ) − [cot (S(t 2 ) + 2φ) + 2Ic] . (A3) The standard semiclassical result σ sc (x) = S(x) + φ is therefore obtained when the bracket in the Eq. above vanishes, namely for a single value of c, given as a function of I and W. The same reasoning can be made for higher order expansions [19]. Note that Eq. (A3) for the semiclassical phase involves in the bracket quantum quantities through the Wronkian and the invariant (these quantities control the quantization condition and the normalization of the wavefunction).

APPENDIX B
We detail the method employed in Sec. 5 to evaluate the quantum phase. The nonlinear Eq. (40), where F is given by Eq. (41) is linearized by expanding the functional to first order in the vicinity of a function M q (x) which is almost identical to M (x) within a given preset accuracy (i.e. it represents the converged answer). This results in the equation F (M q (x), x) is in turn obtained from M q−1 (x) by solving the linear differential equation and so on. Of course in practice the process works the other way round: an initial trial function M 0 (x) is chosen and fed into Eq. (B3) which is solved for M q=1 , which is a better approximation to M than M 0 . The process is repeated with q = 2 in Eq. (B3), which yields a better approximation M 2 ; the iteration stops when M q (x) − M q−1 (x) falls below the preset precision. This process where the solution of the nonlinear differential equation (B1) is replaced by iteratively solving Eq. (B3) is known as the quasilinearization method (QLM). Introduced three decades ago in the context of linear programming it was first employed in a quantum mechanical context as an analytical tool to quantize the Coulomb problem [28]. However the extension of QLM to handle functions with singularities on unbounded domains, as is the case in quantum mechanics, is much more recent: Mandelzweig and co-workers [29] proved that the important property of quadratic convergence which makes this method powerful was still verified, and developed QLM based numerical tools to solve the Schrödinger equation in several cases of potential scattering and bound state problems [30]. The determination of the trial function seems at first sight delicate. In Ref. [28] M 0 (x) = ip(x) was taken to be the most natural choice because it is based on the standard semiclassical scheme (2), and it was shown that in that case M 1 (x) can be obtained as a series of the form c k ζ ′ k (x) where the coefficients c k remind us that M 1 is only an approximation to the converged solution M (x). The important point is that taking M 0 (x) = ip(x) gives iterated solutions that are singular at the turning points. Although this property was employed in [28] to quantize a solvable case (the Coulomb problem) by using an appropriate contour in the complex plane, it is clearly of limited use in the more general case of nonsolvable potentials. This is why implementing the QLM within a real quantum phase approach is advantageous. M (x) has the form given by Eq. (39), and the initial trial function M 0 is sought for in the form A good choice for σ 0 is to take the uniform semiclassical phase σ sc given by Eq. (24). This gives a smooth function M 0 without singularities nor dicontinuities that is known to have the same behaviour as the converged solution. We first thought that σ 0 = σ sc would the only available choice to implement the QLM with a quantum phase approach, but the algorithm turned out to be much more flexible. Indeed, taking (or the restriction of this equation to the half-line which is sufficient for symmetric oscillators) works as well and is faster and easier to implement. Eq. (B5) gives the main components of M in the classically allowed and forbidden regions: indeed in the forbidden region, ∂ x σ is small and tends to 0 whereas ∂ x α/α ∼ ∂ x ψ j /ψ j ∼ ±i |p| where j = 2 (1) in the left (right) forbidden region, whereas ∂ x α/α is small and ∂ x σ ∼ p in the classically allowed region. The discontinuities of Eq. (B5) at the turning points are smoothed out by the finite step size employed in the numerical integration routine.