B G Giraud

S Karataglidis

T Sami
email: sami@subatech.in2p3.fr

Fluctuations of Collective Coordinates and Convexity Theorems for Energy Surfaces

Keywords: Collective coordinates, convexity, energy landscapes, fluctuations

Constrained energy minimizations of a many-body Hamiltonian return energy landscapes e(b) where b ≡ B represents the average value(s) of one (or several) collective operator(s), B, in an "optimized" trial state Φ b , and e ≡ H is the average value of the Hamiltonian in this state Φ b . It is natural to consider the uncertainty, ∆e, given that Φ b usually belongs to a restricted set of trial states. However, we demonstrate that the uncertainty, ∆b, must also be considered, acknowledging corrections to theoretical models. We also find a link between fluctuations of collective coordinates and convexity properties of energy surfaces.

Introduction

Collective coordinates [START_REF] Bohr | Nuclear Structure[END_REF] have been of central importance in descriptions of structure and reactions in atomic, molecular, and nuclear physics. They generate models with far less degrees of freedom than the true number of coordinates, 3A, as needed for a microscopic description of a system of A particles. Often, the system's dynamics can be compressed into slow motions of a few collective degrees of freedom B, while the other, faster degrees can be averaged out. Also, for identical particles, such collective degrees can be one-body operators, B = spin, and if necessary isospin, respectively, of particle i. The summation over i provides for more inertia in B than in the individual degrees β i .

The concept of energy surfaces [START_REF] Nix | Studies in the liquid-drop theory of nuclear fusion[END_REF] has been as important. Given a "coordinatelike" collective operator B and its expectation value b ≡ B , most collective models use an energy function, e(b), and also a b-dependent inertia parameter, that drive the collective dynamics. Keywords such as "saddles", "barriers", etc., flourish [START_REF] Nix | Studies in the liquid-drop theory of nuclear fusion[END_REF].

Simultaneously, it is often assumed that the function, e(b), results from an energy minimization under constraint. Namely, while the system evolves through various values of b, it is believed to tune its energy to achieve a (local) minimum. This aspect of finding e(b) is central to many fields of physics. To illustrate, consider a Hamiltonian, H = i T i + i<j V ij , where T and V denote the usual kinetic and interaction operators. Given a trial set of density operators, D, in many-body space, normalized by TrD = 1, the energy function e(b) may be defined as, e(b) = inf

D⇒b

Tr {HD} ,

where Tr is a trace in the many-body space for the A particles. The constraint, D ⇒ b, enforces Tr {BD} = b.

There are theories which do not use, a priori, an axiom of energy minimization for the "fast" degrees of freedom. Time-dependent Hartree-Fock (HF) [START_REF] Bonche | Time dependent Hartree Fock calculations for heavy ion collisions[END_REF] trajectories, generalizations with pairing, adiabatic versions [START_REF] Villars | Adiabatic time-dependent Hartree-Fock theory in nuclear physics[END_REF], often show collective motions. Equations of motion [START_REF] Bulgac | Adiabatic time-dependent Hartree-Fock theory in the generalized valley approximation[END_REF] and/or a maximum decoupling [START_REF] Matsuyanagi | Open problems in the microscopic theory of large-amplitude collective motion[END_REF] of "longitudinal" from "transverse" degrees, have also shown significant successes in the search for collective degrees, at the cost, however, of imposing a one-body nature of both collective coordinates and momenta and accepting state-dependence of these operators. Such approaches define an energy surface once trajectories of wave functions have been calculated. But they are not the subject of the present analysis. Herein, we focus on fixed operators constraining strict energy minimizations within a fixed basis for single-particle and manybody states. The questions arise: are those constraints themselves subject to fluctuations, and what would be the effect of those fluctuations on the energy minimization? This is not a new question. The issue of constrained Hartree-Fock calculations was addressed in Ref. [START_REF] Giraud | Constrained Hartree-Fock calculations for a quadrupole generator coordinate[END_REF], which considered constrained Hartree-Fock calculations, and corrections to the energy surface. The issue has also been considered more recently in relation to high-energy (e, e ′ p) and (p, 2p) reactions [START_REF] Ryckebusch | Density dependence of quasifree single-nucleon knockout reactions[END_REF], where fluctuations in the position vectors of the target nucleons involved in those processes were considered as going beyond the mean-field approximation assumed for the structure of the target nucleus.

Ideally, to define mathematically a function e(b) of the collective coordinate, one should first diagonalize B within the space provided by the many-body states available for calculations [START_REF] Lesinski | Density functional theory with spatial-symmetry breaking and configuration mixing[END_REF]. The resulting spectrum of B should be continuous, or at least have a high density for that chosen trial space. Then, for each eigenvalue b, one should find the lowest eigenvalue, e(b), of the projection of H into that eigensubspace labeled by b.

In practice, however, one settles for a diagonalization of the constrained operator, H ≡ H -λB, where λ is a Lagrange multiplier, or at least for a minimization of H . Concomitantly, B is assumed to have both upper and lower bounds, or that the constrained Hamiltonian, H, always has a ground state. This returns the "free energy", ε(λ) ≡ H . The label b is no longer an eigenvalue but just an average value, b = B . A standard Legendre transform of ε(λ) then yields the "energy surface", e(b). This utilises the properties dε/dλ = -Tr{BD} = -b, and, de/db = λ. However we show in this work that constrained variation in a quantum system without additional precautions can raise at least two problems, namely: i) the parameter b may no longer be considered as a well-defined coordinate for a collective model due to non-negligible fluctuations; we report cases where the uncertainties, ∆b, can vitiate the meanings of both b and e(b); and ii) there is a link between strict minimization and the curvature properties of e(b) when fluctuations are at work.

Our argument is based fundamentally on a few theorems, presented in Section II, but also illustrated using a few explicit models, which may be solved analytically or numerically. Such solutions for those simple models are presented in Section III. A discussion and conclusion make Section IV.

Theorems linking strict minimization and convexity

Preliminaries

Before proving such theorems, we must recall that, with Hartree-Fock (HF) and Hartree-Bogoliubov (HB) approximations, both convex and concave branches are obtained for e(b) by the addition, in H with the constraint term -λ B , a square term, 1 2 µ B 2 , µ > 0 [START_REF] Giraud | Constrained Hartree-Fock calculations for a quadrupole generator coordinate[END_REF][START_REF] Flocard | Nuclear deformation energy curves with the constrained Hartree-Fock method[END_REF][START_REF] Flocard | Selfconsistent calculation of the fission barrier of 240 Pu[END_REF][START_REF] Staszczak | Augmented lagrangian method for constrained nuclear density functional theory[END_REF], allowing for adjustable values of λ and sometimes µ. It must be stressed that B 2 differs from B 2 . It must also be stressed that such mean-field results are approximations, and that all such methods using a quadratic term in the original function, while stabilizing the numerical procedure, amount to using an effective Lagrange multiplier. Typically, when the "completed" functional, H + 1 2 µ B 2 -λ B , is minimized by a trial function φ in a mean-field approximation, any bra variation, δφ|, induces the condition,

δφ |H| φ + (µ B -λ) δφ |B| φ = 0, (2)
and a similar result holds for any ket variation. Clearly, this means that the combination, H -Λ B , with the effective Lagrange multiplier, Λ = λ-µ B , has been made stationary (but not necessarily minimal). We, therefore, consider the generic form, H = H -λB, in the following discussions.

As far as we are aware, none of the literature on energy surfaces obtained with such quadratic cost functions is concerned with uncertainties in the collective labels. However, any e(b) may find an increase in its uncertainty if b itself, the expectation value of a quantum operator, is subject to large fluctuations. A surprising observation of the present work is that the convexity properties of energy surfaces and uncertainties in the collective coordinates are related. This work clarifies the situation.

Convexity from a convex domain of trial states

Assume that minimization under constraint is performed in a convex domain of many-body density matrices D. Namely, if D 1 and D 2 are trial states, any combination, D m = νD 1 + (1 -ν)D 2 , with 0 < ν < 1, is also a trial state. The constraint operator can be one bounded operator B, or several such operators, B 1 , . . . , B N , or an infinite number of them, such as, for density functionals, local field operators, ψ † r ψ r , labeled by a continuous position r. For notational simplicity, we denote the set of constraint operators by one symbol only, B, and the corresponding expectation values as b ≡ Tr(BD). We call "slice b" the set of all those D returning a given b. We also assume, for simplicity, that, inside "slice b", the minimum e(b) of Tr(HD) is reached for a state D mb , which is non-degenerate.

Then, given two constraint values b 1 and b 2 , with the associated minimalizing states D m1 and D m2 , it is clear that the interpolating state, D mix = νD m1 +(1ν)D m2 , returns both constraint and energy values, b mix = νb 1 + (1 -ν)b 2 and e mix = νe(b 1 )+(1-ν)e(b 2 ). However, there is no reason why D mix should be the energy minimizing state inside that slice labeled by b mix . Hence, by necessity, the energy minimum inside the slice obeys the inequality, e(b mix ) ≤ e mix . In short, e is necessarily a convex function (or functional).

Experimental tables of ground state energies E A are incompatible (see Ref. [START_REF] Barrett | Concavity for nuclear bindings, thermodynamical functions and density functionals[END_REF]) with such a convexity if a density functional is requested to be universal in terms of a particle number A. Indeed, there are many cases where E A ≥ (E A-1 +E A+1 )/2. In such a situation, a convex functional, with a mixed density, (ρ A-1 + ρ A+1 )/2, hence a fluctuation ∆A = 1, will return a lower energy than E A . See Ref. [START_REF] Barrett | Concavity for nuclear bindings, thermodynamical functions and density functionals[END_REF] for ways to "make Nature convex".

Many practical calculations of spectra or energy surfaces use at first mean field methods, where trial states do not belong to a convex domain. For instance, a weighted sum of two Slater determinants usually does not make a Slater determinant. A similar statement holds for HB states. The question of convexity [START_REF] Neergaard | Shell energies of rapidly rotating nuclei[END_REF], however, and its relation to fluctuations, often remains.

A general theorem

Given a constrained Hamiltonian, H -λB, consider a solution branch D(λ), expanding up to second order, assuming that the manifold of solutions is suitably analytic,

D(λ + dλ) = D(λ) + dλ(dD/dλ) + (dλ 2 /2)(d 2 D/dλ 2 ). (3)
In

The zeroth order of this, Eq. ( 5), is, Tr {HdD/dλ}. It vanishes, because of the first of Eqs. [START_REF] Villars | Adiabatic time-dependent Hartree-Fock theory in nuclear physics[END_REF]. The first order, once divided by dλ, gives,

-Tr {BdD/dλ} = -Tr Hd 2 D/dλ 2 . ( 6
)
The left-hand side of Eq. ( 6) is nothing but the second derivative, d 2 ε/dλ 2 . The right-hand side is semi-negative-definite, because of the second of Eqs. (4). Hence, the plot of ε(λ) is a concave curve and the plot of its Legendre transform, e(b), is convex. (Other papers [START_REF] Barrett | Concavity for nuclear bindings, thermodynamical functions and density functionals[END_REF] have the opposite sign convention of the second derivative to define convexity.) With the present sign convention, strict minimization of H necessarily induces convexity of e(b). Any concave branch for e(b) demands an explanation. Note that this proof does not assume any specification of D(λ), whether it is constructed either from exact or approximate eigenstates of H. Therefore strict minimization of (H -λB) can only return convex functions e(b). Maxima are impossible. In the generalization where several collective operators B 1 , . . . , B N , are involved, convexity still holds, so saddles are also excluded. Hence, only an absolute minimum is possible. (However, we shall recall below how to overcome the paradox: by keeping small enough [START_REF] Giraud | Symmetries and fluctuations in the nuclear density functional[END_REF] the fluctuations of the collective coordinate(s), one can deviate from convexity, and more important, validate the quality of the representation provided by branches D(λ).)

Same theorem, for exact solutions

Let ψ(λ) be the ground state of H. (For the sake of simplicity, we assume that there is no degeneracy.) The corresponding eigenvalue, ε(λ), is stationary with respect to variations of ψ, among which is the "on line" variation, dλ (dψ/dλ), leading to the well-known first derivative, dε/dλ = -b ≡ -ψ |B| ψ . Consider the projector, Q = 1-|ψ ψ| . Brillouin-Wigner theory yields the first derivative of ψ, viz.

d |ψ dλ = - Q ε -QHQ B |ψ . (7)
This provides the second derivative of ε,

- db dλ ≡ - d dλ ψ |B| ψ = 2 ψ B Q ε -QHQ B ψ . (8)
Since the operator (ε -QHQ) is clearly negative-definite, the eigenvalue, ε, is a concave function of λ. It is trivial to prove that the same concavity holds for the ground state eigenvalue ε(λ 1 , . . . , λ N ) if several constraints, B 1 , . . . , B N , are used. If, moreover, a temperature T is introduced, the thermal state, D = exp [-H/T ] /Tr exp [-H/T ] , replaces the ground state projector and the free energy also contains the entropy contribution, -T S, where S = -Tr {D ln D}.

A proof of the concavity of the exact ε(λ 1 , . . . , λ N ; T ) is also easy [START_REF] Balian | From Microphysics to Macrophysics: Methods and Applications of Statistical Physics[END_REF]. At T = 0, the usual Legendre transform expresses the energy, e ≡ ψ |H| ψ , in terms of the constraint values, b 1 , . . . , b N , rather than the Lagrange multipliers. For simplicity, consider one constraint only; the generalization to N > 1 is easy. Since e ≡ ε + λb, then de/db = λ, a familiar result for conjugate variables. Furthermore, the second derivative, d 2 e/db 2 , reads, dλ/db = 1/(db/dλ). From Eq. ( 8), the derivative, db/dλ, is positive-definite. Accordingly, e is a convex function of b. Now, if T > 0, the Legendre transform instead generates a reduced free energy, η ≡ (e -T S), a convex function of the constraint value(s). An additional Legendre transform returns e alone, as a convex function of the constraint(s) and S.

Let b -and b + be the lowest and highest eigenvalues of B. When λ runs from -∞ to +∞, then b spans the interval, [b -, b + ]. There is no room for junctions of convex and concave branches under technical modifications as used by Refs. [START_REF] Giraud | Constrained Hartree-Fock calculations for a quadrupole generator coordinate[END_REF][START_REF] Flocard | Nuclear deformation energy curves with the constrained Hartree-Fock method[END_REF][START_REF] Flocard | Selfconsistent calculation of the fission barrier of 240 Pu[END_REF][START_REF] Staszczak | Augmented lagrangian method for constrained nuclear density functional theory[END_REF]. For every exact diagonalization of H, or exact partition function, convexity sets a one-to-one mapping between b in this interval and λ. More generally, with exact calculations, there is a one-to-one mapping between the set of Lagrange multipliers, {λ 1 , . . . , λ N } , and that of obtained values, {b 1 , . . . , b N }, of the constraints. Convexity, in the whole obtained domain of constraint values, imposes a poor landscape: there is one valley only.

We tested this surprising result with several dozens of numerical cases, where we used random matrices for H and B, with various dimensions. As an obvious precaution, we eliminated those very rare cases where both H and B turned out, by chance, to be block matrices with the same block structure; such cases give rise to level crossings and degeneracies. Then every remaining situation, without exception, confirmed the convexity of e.

= 0 -2 0 3 -2 2 -2 -2 0 -2 -4 -1 3 -2 -1 7
for instance. Such branches are easily derived algebraically [START_REF] Giraud | Algebraic density functionals[END_REF] from the polynomial, P ≡ det(H -λB -ε). The convexity of the ground state branch is transparent. Clear also are its infinite derivatives when b reaches b -= -4.83 and b + = 9.60, and the vanishing derivative at the point corresponding to the unconstrained ground state, where b = -2.69 and e = e -= -9.71.

Same theorem, for approximations via constrained HF calculations

Consider now energy surfaces obtained from approximations. Typically, one uses a HF or HB calculation, at zero or finite T . Trial states in such methods span a nonlinear manifold; indeed, a sum of two determinants is usually not a determinant. Let D(λ) denote one A-body density operator where, within such nonlinear approximations, a minimum, ε(λ), of Tr {HD} or of (Tr {HD} -T S), is reached. It may be degenerate, but, in any case, it is stationary for arbitrary variations δD within the set of trial states. Accordingly, the first derivative again reads, dε/dλ = -b ≡ -Tr {BD}. Then, if a Legendre transform holds, defining η ≡ ε + λb in terms of b, the same argument that was used for the exact case again yields, dη/db = λ. With N constraints, the gradient of η in the domain spanned by {b 1 , . . . , b N } is the vector {λ 1 , . . . , λ N }.

To discuss second derivatives, consider, for instance, HF calculations, where A-body density operators are dyadics of determinants, D = |φ φ|. Normconserving variations of an HF solution, φ, can be parametrized as, |δφ = exp(iXδα) |φ , with X an arbitrary particle-hole Hermitian operator, and δα an infinitesimal coefficient. Under such a variation in the neighborhood of a HF solution, the first and second order variations of the free energy, ε ≡

Tr {H exp(iXδα)D exp(-iXδα)}, read, δε = iδαTr {[H, X] D} = 0, (9)
and

δ 2 ε = -δα 2 /2 Tr {[[H, X] , X] D} ≥ 0, (10)
respectively. If D is a HF solution, the first order vanishes ∀X. Since only those solutions that give minima are retained, the second order variation of ε is semipositive-definite, ∀X again. Now, when H receives the variation, -Bdλ, there exists a particle-hole operator, Y , a special value of X, that, with a coefficient dλ, modifies the solution. This reads |Φ = exp(iY dλ) |φ . Simultaneously, those particle-hole operators that refer to this new Slater determinant Φ become X = exp(iY dλ)X exp(-iY dλ). The new energy is,

ε Y = Tr {exp(-iY dλ)(H -Bdλ) exp(iY dλ)D} . (11)
The stationarity condition, Eq. ( 9), becomes,

0 = Tr {exp(-iY dλ) [(H -Bdλ), X ] exp(iY dλ)D} = Tr {[exp(-iY dλ)(H -Bdλ) exp(iY dλ), X] D} . ( 12
)
The zeroth order in dλ of this, Eq. ( 12), reads, Tr {[H, X] D}, and vanishes ∀X, because of Eq. ( 9). Then the first order in dλ gives, again ∀X,

Tr {[B, X] D} = iTr {[[H, Y ] , X] D} . (13)
The second derivative is,

d 2 ε/dλ 2 = -(d/dλ)Tr {exp(-iY dλ)B exp(iY dλ)D} = -iTr {[B, Y ] D} . ( 14
)
Upon taking advantage of Eq. ( 13), for Y as a special case of X, this becomes,

d 2 ε/dλ 2 = Tr {[[H, Y ] , Y ] D} , (15)
the right-hand side of which is semi-negative-definite, see Eq. [START_REF] Flocard | Nuclear deformation energy curves with the constrained Hartree-Fock method[END_REF]. The solution branch obtained when λ runs is, therefore, concave. Its Legendre transform is convex.

About concave branches

Let us return to the "quadratic cost" function, that redefines the "free energy" to be minimized, E, as, E = H + 1 2 µ B 2 -λ B , or, in a shorter notation, E = e + 1 2 µb 2 -λb. Let the positive number µ be kept constant, allowing λ to vary. It is easy to show again that dE/dλ = -b. But now one finds that, de/db = Λ ≡ λ -µb, where the effective Lagrange multiplier Λ is that noticed at the stage of Eq. ( 2). Assume again that variational solutions are analytic with respect to λ. For the sake of simplicity, use wave functions φ rather than density operators, hence D = |φ φ|. The analog of Eq. ( 3) reads, φ(λ + dλ) = φ + dλ φ ′ + 1 2 (dλ) 2 φ ′′ , with short notations for derivatives of φ. The stationarity of E then reads,

0 = E ′ ≡ φ ′ |(H -ΛB)|φ + φ|(H -ΛB)|φ ′ . ( 16
)
If this "quadratic cost" E is a strict minimum with respect to variations of φ, one finds,

0 < E ′′ ≡ φ ′′ | (H -ΛB) |φ + 2 φ ′ | (H -ΛB) |φ ′ + φ| (H -ΛB) |φ ′′ + µ ( φ ′ |B|φ + φ|B|φ ′ ) 2 . (17)
When λ becomes λ+dλ, the derivative of the wave function becomes (φ ′ +dλ φ ′′ ) and the functional to be minimized receives an additional contribution, -dλ B .

Then the stationarity condition at φ(λ + dλ) becomes,

0 = (E ′ + dλ E ′′ ) -dλ ( φ ′ |B|φ + φ|B|φ ′ ) . (18)
This simplifies into, E ′′ = db/dλ, since E ′ = 0, see Eq. ( 16), and since, simultaneously, db/dλ = φ ′ |B|φ + φ|B|φ ′ . Furthermore, since E ′′ > 0, see Eq. ( 17), and since we know that dE/dλ = -b, we find that d 2 E/dλ 2 = -E ′′ , a negative quantity. The function, E(λ), is, therefore, necessarily concave if a strict minimization of E has been performed in the domain of trial states. As stated at the beginning of this Subsection, the "energy surface" function, e(b), now obeys the condition, de/db = λ-µ b, hence d 2 e/db 2 = dλ/db-µ. While we found that dλ/db = 1/(db/dλ) = 1/E ′′ is a positive quantity, the term, -µ, competes with dλ/db and may induce negative values of d 2 e/db 2 , hence concave regions in the plot of e(b), where, incidentally, the second derivative, d 2 e/db 2 , cannot be more negative than -µ.

Consider two solutions, φ 1 , φ 2 , generating two points, (b 1 , e 1 ), (b 2 , e 2 ), that are separated by an inflection point of e(b) and that show the same derivative, Λ = (de/db) 1 = (de/db) 2 , in other words, λ 1 -µb 1 = λ 2 -µb 2 . Denote φ 1 the lower energy solution, hence e 1 < e 2 . Then φ 2 , compared with φ 1 , is an excited solution of the variational stationarity problem for the same constrained Hamiltonian, H -ΛB, that drives φ 1 . For the sake of the argument, assume that b 1 < b 2 . Then along the convex branch of e(b), the value of Λ increases when b increases, and, when b decreases along the concave branch, the value of Λ also increases. When one reaches the inflection point, (b i , e i ), a maximum value, Λ i , is reached. One must conclude that the two, formerly distinct solutions, φ 1 (Λ) and φ 2 (Λ), do not cross and, rather, smoothly fuse into a unique solution φ(Λ i ). This is a very unlikely situation for exact solutions, namely eigenstates, as is well known. Only approximations can afford such an anomaly. (Notice also that, because H -ΛB is assumed bounded from below, ∀Λ, there will exist a stationary solution φ + for any multiplier Λ + > Λ i , even infinitesimally close to Λ i . Then one must accept that Φ + and Φ i strongly differ from each other.)

The zoo of stationary solutions of approximate methods such as mean field methods (non linear!) can be rich enough to accommodate such singularities. This makes a paradox: would non-linear approximations generate a more flexible, physical tool than exact solutions? "Phase transitions", a somewhat incorrect wording for a finite system, are sometimes advocated to accept continuing branches of energy minima into metastable branches. But this definitely claims some caution with the axiom of strict energy minimization to freeze fast degrees. This need for excited solutions in constrained mean field calculations, namely two solutions for the same value of a Lagrange multiplier to describe both sides of an inflection point of the barrier, is well known and used. See in particular [START_REF] Flocard | Nuclear deformation energy curves with the constrained Hartree-Fock method[END_REF][START_REF] Flocard | Selfconsistent calculation of the fission barrier of 240 Pu[END_REF], where a tangent parabola rather than a tangent straight line is used to explore an energy surface by mean field methods.

Bimodal solutions in mean field approximations

Besides the caution about the "fast degree minimization hypothesis" one should consider whether such mean field solutions, in convex or concave branches, might be vitiated by large uncertainties for b. The following solvable models give a preliminary answer.

Consider N identical, 1-D fermions with Hamiltonian,

H = N i=1 p 2 i /(2m) + M ω 2 R 2 /2 + N i>j=1 v ij , (19)
where R = i r i /N is the center-of-mass (c.m.) position, p i , r i , m denote the single particle momentum, position and mass, respectively, of each fermion, and M = N m is the total mass. The c.m. momentum is, P = i p i . We use a system of units such that = m = ω = 1, where ω denotes the frequency of the c.m. harmonic trap. The interaction, v, is set as Galilean invariant and so is the sum, V = i>j v ij . In the following, v is taken as a spin and isospin independent and local force,

v ij = v(|r i -r j |).
The collective operator we choose to constrain H is a half sum of "inertia" (mass weighted square radii), B = N i=1 mr 2 i /2. The constrained Hamiltonian then reads,

H = i p 2 i /(2m) + M ω 2 R 2 /2 + i>j v ij -λ i mr 2 i /2. ( 20
)
Let

ξ 1 = r 2 -r 1 , ξ 2 = r 3 -(r 1 + r 2 )/2, ..., ξ N -1 = r N -(r 1 + r 2 + • • • + r N -1
)/(N -1) denote the usual Jacobi coordinates with Π α and µ α , α = 1, . . . , (N -1), the corresponding momenta and reduced masses. In this Jacobi representation the constraint becomes, B = M R 2 /2 + α µ α ξ 2 α /2. Accordingly, the constrained Hamiltonian decouples as a sum of a c.m. harmonic oscillator,

H c.m. = P 2 (2M ) + 1 2 M Ω 2 R 2 , Ω 2 = ω 2 -λ, (21)
provided λ < ω 2 , and an internal operator,

H int = α Π 2 α /(2µ α ) -λµ α ξ 2 α /2 + V. ( 22
)
With the present power of computers and present experience with Faddeev(-Yakubovsky) equations, this choice of H and B, with its ability to decouple, provides soluble models with, typically, N = 2, 3, 4. (Decoupling also occurs if B is a quadrupole operator.) Exact solutions can thus be compared with mean field approximations and validate, or invalidate, the latter. Here, however, we are not interested in the comparison, but just in properties of the mean field solutions as regards B.

For this, as long as N does not exceed 4, we assume that the Pauli principle is taken care of by spins and isospins, understood in the following, and that the space part of the mean field approximation is a product, φ = i ϕ(r i ), of identical, real and positive parity orbitals. The corresponding Hartree equation reads,

p 2 2m + Λr 2 + u(r) ϕ(r) = ε sp ϕ(r), (23)
with u(r) = (N -1)

∞ -∞ v(r -s)[ϕ(s)] 2 and Λ = M ω 2 /(2N 2 ) -λm/2.
The term, M ω 2 /(2N 2 ), clearly comes from the c.m. trap. The same trap induces a two-body operator, M ω 2 i =j r i r j /(2N 2 ), which cannot contribute to the Hartree potential, u, since every dipole moment, r j , identically vanishes here.

Once ϕ and ε sp are found, one obtains the free energy, We show now, among many cases we studied, Hartree results when

ε Hart = H = N ε sp - N ϕ |u| ϕ /2,
v ij = -2[exp(-2(r i -r j + 8) 2 ) + exp(-2(r i -r j -8) 2 ) + 2 exp(-2(r i -r j + 4) 2 ) + 2 exp(-2(r i -r j -4) 2 ) + exp(-2(r i -r j ) 2 )], see Fig. 2.
The results shown in Figs. 234correspond to N = 2. Similar results hold with N = 3, 4. Orbitals ϕ are expanded in the first 11 even states of the standard harmonic oscillator. The only difference between the two orbitals shown in Fig. 3 is the value of λ. Both orbitals, and many other ones, when λ runs, show a bimodal structure, their left-hand-side and right-hand-side bumps being equivalent with respect to the even observable, B ∝ r 2 . Whether N = 2, 3, or 4, we found many cases where the value of B has nothing to do with the positions of the peaks of ϕ. The bad quality induced by the corresponding uncertainties is illustrated by Fig. 4. Therefore, little trust is available for the energy curve, e Hart (b), that results from this model.

Models showing the existence and extent of quantum fluctuations of collective coordinates

We now turn our attention to illustrating the effects of fluctuations of collective coordinates, as discussed above, with various examples.

Fluctuation of the mean square radius in a harmonic oscillator shell model

Consider A fermions with Hamiltonian,

H = 1 2 A i=1 (p 2 xi + p 2 yi + p 2 zi + x 2 i + y 2 i + z 2 i ), or, in a shorter notation, H = i (p 2 i + r 2 i )/2.
Ignore spin (and isospin) and fill completely each shell up to that one with energy n + 3 2 . It is trivial that the sequence of individual shell populations reads, 1, 3, 6, . . . , (n + 1)(n + 2)/2, and that, when all shells up to and including "shell n" are filled, the particle number reads, A = (n + 1)(n + 2)(n + 3)/6. Denote φ that Slater determinant made of such (n + 1) filled shells and choose now the collective operator, B = A -1 A i=1 r 2 i . By definition, b ≡ B represents the square of an average radius, r = √ b. For the sake of mathematical rigour, it might be necessary to redefine B with a cut-off, such as,

B c = r 2 i exp(-r 2 i /R 2 c
), with R c significantly larger than atomic, molecular or nuclear radii, to ensure that constrained Hamiltonians, H -λB c , do have a ground state if λ > 0. But this technicality can be neglected in practice, with calculations using a large but finite basis of finite range variational states.

Each shell with index m, m = 0, . . . , n, contributes A -1 (m+3/2)(m+1)(m+ 2)/2 to φ|B|φ , hence a total B = A -1 (n + 1)(n + 2) 2 (n + 3)/8. Consider now B 2 , to calculate the fluctuation. Let h be those states filled in φ and P the other eigenstates of H. The one-body operator B can excite at most one-particle-one-hole states when acting upon φ, hence,

B 2 = φ|B|φ φ|B|φ + P h φ|B|φ P h φ P h |B|φ . ( 24
)
Accordingly, (∆b) 2 reduces to a sum of particle-hole squared matrix elements, (∆b) 2 = A -2 P h ( P |r 2 |h ) 2 . Slightly tedious, but elementary manipulations, yield the result, (∆b)

2 = A -2 (n + 1)(n + 2) 2 (n + 3)/8. This gives ∆b/b = √ 8/[(n + 2) (n + 1)(n + 3)] ∝ A -2/3
. In systems with ≃ 10 3 particles a "horizontal" bar of order 1% might likely be neglected, but many systems studied in atomic, molecular or nuclear physics rather deal with a few scores of particles or hardly two or three hundred of them, and a horizontal uncertainty ranging between ≃ 8% and ≃ 2% might trigger some attention.

Fluctuation of the quadrupole moment in a deformed harmonic oscillator shell model

The Hamiltonian now reads,

H = i [p 2 i + x 2 i + y 2 i + z 2 i /γ 4 ]/2
, hence prolate shapes if γ > 1 and oblate ones if γ < 1. The collective operator is now chosen as, B = i β i with β i = 2z 2 i -x 2 i -y 2 i , without a precaution cut-off. For values of γ close enough to 1 to induce a weak splitting of levels we can study the same scheme of fully filled shells as for the model just above. Orbitals ϕ i in a "shell m" have wave functions χ m ′ (x) χ m ′′ (y) χ m ′′′ (z/γ)/ √ γ, where the χ's are the standard 1-D harmonic oscillator states and

m ′ + m ′′ + m ′′′ = m. The orbital energies are, m ′ + m ′′ + 1 + (m ′′′ + 1/2)/γ 2 , obviously.
When all the orbitals ϕ i of the shells with labels 0, 1, . . . , n are filled, the particle number for spinless and isopinless fermions is again, A = (n + 1)(n + 2)(n + 3)/6 and the value of the quadrupole is, obviously, b ≡ B = A i=1 ϕ i |β|ϕ i . Its quantum fluctuation is given by the same particle-hole summation, (∆b) 2 = P h ( P |β|h ) 2 . This also reads, (∆b) 2 = Tr β (1-σ) β σ, where we use the onebody density operator, σ = A i |ϕ i ϕ i |, namely the projector upon holes. Its complement, (1 -σ), is the projector upon the particle subspace. The trace, Tr, runs in one-body space. This yields, (∆b) 2 = Tr β 2 σ-Tr β σ β σ. Since σ is here real symmetric and β is a local operator, the first and second terms read, in coordinate representation, d r σ( r, r) [β( r)] 2 , and, -d r d r ′ β( r) [σ( r, r ′ )] 2 β( r ′ ), respectively. The symbol r is here a short notation for the three coordinates x, y, z.

Quadrupole values, b = (γ 2 -1)(n+1)(n+2) 2 (n+3)/12, are easily found, with n the label of the highest filled "shell", and, accordingly, the particle number, A = (n + 1)(n + 2)(n + 3)/6. Because of the factor, (γ 2 -1), such values of the quadrupole do vanish when γ = 1, as expected for spherical shells, and also show the correct signs for prolateness and oblateness.

Similar brute force considerations provide (∆b) 2 = (1 + 2γ 4 )(n + 1)(n + 2) 2 (n + 3)/12. It would not be significant to discuss the ratio, ∆b/|b|, when γ is close to 1, obviously, but for a value of γ such as 3/2, for instance, the "full shell" particle numbers considered above return the following sequence of values for this relative uncertainty, 2.7, 1.1, 0.60, 0.38, 0.26, 0.19, 0.15, 0.11, 0.09, 0.08, 0.06. If γ = 1/2, the sequence becomes, 1.4, 0.58, 0.32, 0.20, 0.14, 0.10, 0.08, 0.06, 0.05, 0.04, 0.03. Both sequences show that the relative uncertainty shrinks when the particle number increases, as should be expected, but, as long as the particle number does not exceed a few hundred, fluctuations of ≃ 5% cannot be neglected.

Actually, level crossing occurs when deformation sets in. Full fillings of previously spherical shells do not represent ground states of deformed Hamiltonians. Set γ = 4/3 for instance, and, given a particle number A, fill the lowest orbitals of the corresponding deformation scheme. This gives quadrupole, mean square deviations and relative uncertainty values listed in Table I. In If spin, for electrons, and both spin and isospin are reinstated, for nucleons, keeping the same occupied levels while particle number is multiplied by 2 and 4, respectively, it is obvious that ∆b/|b| is divided by √ 2 and by 2, respectively. Actually, for nucleons, the proton number is usually smaller than that of neutrons, and orbitals may somewhat differ, hence the exact reduction factor of ∆b/|b| will slightly differ from 2, but this changes nothing to the fact that quantum fluctuations exist and a relative uncertainty of at least a few percent cannot be avoided.

Fluctuation in a constrained diagonalization

We consider in this section a one-dimentsional Hamiltonian, H = -d 2 /dr 2 + v(r), with a double hump potential, v = (r -1/5) 6 /5000 e -r 2 /12 + (r/12) 6 /8, shown as a full line in Fig. 7.

We take the constrained Hamiltonian as H = H -λ r and diagonalize it in a subspace spanned by shifted Gaussians. While the eigenstate, ψ λ (r), shows a single wave-packet when the average value, r ≡ ψ λ |r| ψ λ , sits near a minimum of v, an expected tunnel effect occurs when r sits near a maximum of v. There ψ shows two connected packets, one at each side of the barrier, inducing a lowering of the energy. Such a bimodal (even multimodal in several extreme cases we tested) situation induces a very bad probing of the barrier. The plot of the energy, e = H , (multiplied ten times for graphical reasons) in terms of r does not reflect v(r) in any way, see the thin full line in Fig. 8.

When tunnel effects occur, fluctuations, ∆r = r 2 -r 2 , are dramatically larger than when a unique packet sits in a valley. The label, r in our case, is thus misleading. Although H was exactly diagonalized, constrained variation generated a bad quality representation of v(r). Fig. 8 illustrates how big the uncertainty on r can become. Moreover, this "energy surface" turns out to be convex.

The previous models in this Section show that somewhat large quantum fluctuations of collective coordinates do exist, but the present model, in this Subsection, illustrates a new fact, namely that such fluctuations may vary along an energy surface and may influence the surface itself. A trivial way to prevent fluctuations from arbitrarily varying is to introduce a double constraint via the square, B 2 , of the initial constraint operator [START_REF] Giraud | Symmetries and fluctuations in the nuclear density functional[END_REF]. One adjusts the second Lagrange multiplier so that the fluctuation, ∆b, remains small, and, for a stable quality of the representation, reasonably constant. (Alternately one can tune the second Lagrange multiplier to ensure more or less constant and/or small enough ratios ∆b/|b|.) Again with 1-D toy Hamiltonians of the form, H ′ = -d 2 /dr 2 + v(r)-λ 1 r + λ 2 r 2 , or, equivalently, H ′ = -d 2 /dr 2 + v(r)+ λ 2 (r -λ 1 ) 2 , we tuned λ 2 into a function λ 2 (λ 1 ) to enforce a unimodal situation with ∆r kept constant when r evolves. A sharp, stable probe of the barrier results. Convexity is then "defeated". The true shape of v is recovered [START_REF] Giraud | Symmetries and fluctuations in the nuclear density functional[END_REF].

Discussion and conclusion

The present work shows that expectation values of collective coordinate operators may lead to misleading coordinates for an energy surface. Convexity is a major property of any energy surface obtained by an exact minimization of the energy under constraint(s), and may occur through tunnel effects and collective coordinate fluctuations. If an energy landscape with "saddles" is needed, such deviations from convexity contradict the requirement that the energy must transit through exact minima. The success of collective models that use a nontrivial landscape via mean field approximations is too strong to be rejected as physically and/or mathematically unsound, but its validation likely relates to further methods such as, resonating group methods [START_REF] Wildermuth | A unified theory of the nucleus[END_REF], generator coordinate (GC) ones [START_REF] Hill | Nuclear constitution and the interpretation of fission phenomena[END_REF][START_REF] Griffin | Collective motions in nuclei by the method of generator coordinates[END_REF], Born-Oppenheimer approximations, influence functionals [START_REF] Feynman | The theory of a general quantum system interacting with a linear dissipative system[END_REF], deconvolutions of wave packets in collective coordinate spaces [START_REF] Giraud | Energy surfaces, collective potentials and collective inertia parameters[END_REF], etc. In particular, one can argue that, while individual exact states ψ λ or mean-field ones φ λ may carry large uncertainties for the label b, such states may still provide a good global set for a GC calculation. But then the physics lies as much in non-diagonal elements H(b, b ′ ) and N (b, b ′ ) of the GC energy and overlap kernels, respectively, than in the diagonal, e(b) ≡ H(b, b).

Even so, except for anharmonic vibrations, where one valley is sufficient, one has to justify why concave branches can be as significant as branches obtained from exact minimizations of the energy.

Because of kinetic terms, which enforce delocalizations, a full Hamiltonian is often not well suited to ensure a good localization of the operators, B, a necessary condition for the exploration of an energy surface parametrized by their expectation values, b. Recall that, in the Born-Oppenheimer treatment of the hydrogen molecule, the proton kinetic energy operator is initially removed from the Hamiltonian, allowing the collective coordinate, namely the interproton distance, to be frozen as a zero-width parameter. Most often in nuclear, atomic, and molecular physics, such a removal is not available. Hence, such an approach may lead to "dangerous consequences", as illustrated by Fig. 8, from the toy Hamiltonian introduced at the stage of Fig. 7. The growth of fluctuations due to tunnel effects is spectacular. Fortunately, if one forces the constrained eigenstate to retain a narrow width while the collective label runs, convexity is defeated at the profit of a reconstruction of the potential shape.

The technical devices used with success in Refs. [START_REF] Giraud | Constrained Hartree-Fock calculations for a quadrupole generator coordinate[END_REF][START_REF] Flocard | Nuclear deformation energy curves with the constrained Hartree-Fock method[END_REF][START_REF] Flocard | Selfconsistent calculation of the fission barrier of 240 Pu[END_REF][START_REF] Staszczak | Augmented lagrangian method for constrained nuclear density functional theory[END_REF] require a further investigation of the solutions they generated. Uncertainties in the collective coordinates must be acknowledged. "Mountains" can be underestimated, as is the case when tunnel effects deplete the energy probing wave function. (A hunt for multimodal solutions would therefore be useful.) Since fluctuations are important at "phase transitions", collective operators must be completed by their own squares, in combinations of the form, K = H -λB + µ(λ)B 2 , with µ(λ) adjusted to avoid wild increases of ∆b . Such operators K govern both a constraint and its fluctuation, but obviously differ from a double constraint form, H = H -λ 1 B + λ 2 B 2 , with two independent parameters, λ 1 , λ 2 . The second derivative, d 2 ε/dλ 2 contains additional terms due to dµ/dλ and d 2 µ/dλ 2 , hence a one-dimensional path with non convex structures can be induced by K inside that convex two-dimensional landscape due to H. We verified this "adjustable B 2 " method with many numerical tests.

In theories using, partly at least, liquid drop models, see Ref. [START_REF] Moller | Heavy-element fission barriers[END_REF] for instance, labels b are purely classical. Such theories are thus safe from the present concerns. But many other energy surfaces used in realistic situations come from mean-field constrained calculations. It remains to be tested whether their solutions, stable or metastable, carry a mechanism that diminishes the fluctuation of collective degrees of freedom. This mechanism, if it exists, deserves investigation. We conclude that a review of landscapes obtained by constrained HF or HB is in order, to analyze the role of collective coordinate fluctuations. It is clear that such surfaces deserve corrections because of likely variable widths ∆b of their collective observables, and also, obviously, because convolution effects and zero-point energies must be subtracted.

To summarize, we first found that collective operators carry significant fluctuations. We also found that fluctuations, and convexity situations, can enforce poor energy landscapes if constrained energy minimizations are used with fixed operators B and fixed trial spaces. Unacceptable uncertainties, ∆b, can vitiate the meaning of collective labels. We did even discover bi-or multimodality in mean field approximations. Fortunately, given the same fixed operators and trial spaces, a modest deviation from fixed constraints, namely adjustable combinations of B and B 2 , commuting operators indeed, allows an analysis "at fluctuations under control", with unimodal probes of landscapes and a controlled quality of the collective representation. A puzzling question remains: is the good quality of constrained mean field solutions in the literature [START_REF] Flocard | Nuclear deformation energy curves with the constrained Hartree-Fock method[END_REF][START_REF] Flocard | Selfconsistent calculation of the fission barrier of 240 Pu[END_REF] the result of a "self damping" of collective coordinate fluctuations? Are multimodal situations blocked when nuclear mass increases? SK acknowledges support from the National Research Foundation of South Africa. BG thanks N. Auerbach for stimulating discussions and the University of Johannesburg and the Sackler Visiting Chair at the University of Tel Aviv for their hospitality during part of this work.

Figure 1

 1 shows all branches of e(b) for H =

Figure 1 :

 1 Figure 1: Branches obtained for e(b) by algebraic elimination of the Lagrange multiplier λ between b(λ) and e(λ) in a toy case with 4 × 4 matrices. Convexity of the ground state branch.

  then the value of the constraint, B Hart = N m ϕ r 2 ϕ /2, and the square fluctuation, (∆b) 2 Hart = N m 2 ϕ r 4 ϕ -ϕ r 2 ϕ 2 /4. The physical energy, e Hart (b), in this Hartree approximation, clearly obtains by adding λ B Hart to ε Hart .

Figure 2 :

 2 Figure 2: An interaction giving multimodal Hartree solutions.

Figure 3 :

 3 Figure 3: Dashes, Hartree orbital ϕ(r) if N = 2, λ = .47 for the interaction shown in Fig. 2. Dots, the same for λ = .55.

Figure 4 :

 4 Figure 4: Same Hartree model. ∆b as a function of B .

  Figure 5 we show, in terms of A, quadrupole values b along a full line and values b ± ∆b along dashed lines. With γ=3/4, typical results are listed in

Figure 5 :

 5 Figure 5: Quadrupole values b (full line) and b ± ∆b (dashed lines) for γ = 4/3.

Figure 6 :

 6 Figure 6: Quadrupole values b (full line) and b ± ∆b (dashed lines) for γ = 3/4.

Figure 7 :Figure 8 :

 78 Figure 7: Toy model, 1-D potential, full curve. Tunnel effect of constrained eigenstate under left barrier, dashed curve. Same effect under right barrier, dotted curve.

  those cases where pure states, D = |φ φ|, are used, the wave function of interest is assumed to be analytical, and, clearly, dD/dλ = |φ dφ/dλ| + |dφ/dλ φ|, and also, d 2 D/dλ 2 = |φ d 2 φ/dλ 2 | + 2 |dφ/dλ dφ/dλ| + |d 2 φ/dλ 2 φ|.

	The free energy ε is also stationary for D(λ + dλ), but the Hamiltonian is now, H(λ) -B dλ, and the derivative of the state is, dD/dλ + dλ(d 2 D/dλ 2 ) + O(dλ 2 ), hence,
	Tr (H -Bdλ) dD/dλ + dλ(d 2 D/dλ 2 ) + O(dλ 2 ) = 0.
	The stationarity and minimality of Tr {HD} with respect to any variation of D, and in particular with respect to that variation, D(λ + dλ) -D(λ), induce,
	Tr {HdD/dλ} = 0, Tr Hd 2 D/dλ 2 ≥ 0.	(4)

Table 1 :

 1 Table II and shown in Fig. 6. Typical results for γ = 4/3.Relative fluctuations ∆b/|b| decrease when A increases, as should be expected, but an order of magnitude of several percent at least cannot be avoided.

	A	Quadrupole Square deviation Relative uncertainty
		27	575/3	0.51
		920/9	44200/81	0.23
		1450/9	72230/81	0.19
		2300/9	110500/81	0.14
		259	14069/9	0.15
		3340/9	173300/81	0.13
		4535/9	222241/81	0.10
		4777/9	252863/81	0.11
		6260/9	311020/81	0.089
	102	6886/9	353930/81	0.086
	109	2695/3	133841/27	0.078
	118	2776/3	142952/27	0.079
	130	9338/9	484054/81	0.075
	A	Quadrupole Square deviation Relative uncertainty
	9	-169/16	3647/128	0.51
		-389/8	6019/64	0.20
		-225/4	4095/32	0.20
		-173/2	3067/16	0.16
		-1185/8	17895/64	0.11
		-1399/8	23585/64	0.11
		-3685/16	56867/128	0.092
		-3815/16	64945/128	0.094
		-2205/8	36459/64	0.087
	101	-5531/16	90301/128	0.077
	107	-5837/16	97339/128	0.076
	116	-3281/8	54247/64	0.071
	129	-7973/16	125443/128	0.063

Table 2 :

 2 Typical results for γ = 3/4.