Nature Communications

Dicke model and chaos

The Hamiltonian of the Dicke model is written as

The eigenvalues j(j + 1) of the squared total spin operator $$ {\widehat{\mathbf{J}}}^2={\widehat{J}}_x^2+{\widehat{J}}_y^2+{\widehat{J}}_z^2$$ specify the different invariant subspaces of the model. We use the symmetric atomic subspace defined by the maximum pseudo-spin $$ j=\mathcal{N}/2$$, which includes the ground state. When the Dicke model reaches the critical value $$ {\gamma }_c=\sqrt{\omega {\omega }_0}/2$$, it goes from a normal phase (γ < γ_c) to a superradiant phase (γ > γ_c). Our studies are done in the superradiant phase, γ = 2γ_c, and we choose ω = ω₀ = 1. The rescaled energies are denoted by ϵ = E/j. For the selected parameters, ϵ_GS = −2.125 is the ground-state energy.

The classical Hamiltonian, h_cl(x) in the coordinates x = (q, p; Q, P), is obtained by calculating the expectation value of the quantum Hamiltonian under the product of bosonic Glauber and pseudo-spin Bloch coherent states $$ ∣\mathbf{x}⟩=∣q,p⟩\otimes ∣Q,P⟩$$ (see Methods) and dividing it by j. The effective Planck constant ℏ_eff = 1/j⁴² determines the resolution of the quantum states on the four dimensional phase space $$ \mathcal{M}$$. We are able to work with large system sizes (j ~ 100 and Hilbert space dimensions D ~ 6 × 10⁴), due to the use of an efficient basis that guarantees the convergence of a broad range of eigenvalues and eigenstates^43,44 (see Methods). The Dicke model displays regular and chaotic behavior. For the Hamiltonian parameters selected in this work, the system is in the strong-coupling hard-chaos regime for ϵ > −0.8 (see Supplementary Note 1).

Quantum scarring

The (unnormalized) Husimi function of a state $$ \widehat{\rho }$$ is defined as $$ {\mathcal{Q}}_{\widehat{\rho }}\left(\mathbf{x}\right)=⟨\mathbf{x}∣\widehat{\rho }∣\mathbf{x}⟩$$, which is the expectation value of the density matrix over the coherent state $$ ∣\mathbf{x}⟩$$ centered at x. This function is used to visualize how the state $$ \widehat{\rho }$$ is distributed in the phase space. Quantum scars are localized around the classical periodic orbits in an energy shell of the phase space. To visualize the scars, we consider the Husimi projection over the classical energy shell at ϵ,

By integrating out the bosonic variables (q, p), the remaining function can be compared with the projection of the classical periodic orbits on the plane of atomic variables (Q, P).

Identifying all periodic orbits that generate the scars of a quantum system is extremely challenging. We were able to identify two families of periodic orbits for the Dicke model, which we denote by family $$ \mathcal{A}$$ and family $$ \mathcal{B}$$⁴⁵. By calculating the overlap of the eigenstates with tubular phase-space distributions located around these orbits⁸, we selected twelve eigenstates $$ {\widehat{\rho }}_k=∣E_k⟩⟨E_k∣$$ scarred by those two different families. In Fig. 1 we plot their Husimi projections $$ {\tilde{\mathcal{Q}}}_k={\tilde{\mathcal{Q}}}_{{ϵ}_k,{\widehat{\rho }}_k}$$ at ϵ_k = E_k/j along with the corresponding periodic orbit of each family. Family $$ \mathcal{A}$$ (solid blue line in Fig. 1a1–a6) contains the periodic orbits of lowest period of the Dicke model, which emanate from one of the two normal modes around a stable stationary point at the ground-state energy. Family $$ \mathcal{B}$$ (solid red line in Fig. 1b1–b6) arises from the other normal mode around the same point. Scarring is clearly visible in all panels of Fig. 1. The quantum states are highly concentrated around the classical periodic orbits. This happens even in the chaotic region of high excitation energy, where the classical dynamics is ergodic, as seen in Fig. 1a5, a6, b5, b6. Notice that the eigenstates may be scarred by more than one periodic orbit. In fact, as we showed quantitatively in ref. ⁴⁵ after introducing a measure of scarring, an eigenstate may even be scarred by periodic orbits of different families. This means that different eigenstates exhibit different degrees of scarring.

Fig. 1

Classical periodic orbits and scars in the Husimi projection of eigenstates.

a1–a6, b1–b6 Projected Husimi distribution $$ {\tilde{\mathcal{Q}}}_k\left(Q,P\right)$$ superposed by periodic orbits from the family $$ \mathcal{A}$$ (blue lines) [family $$ \mathcal{B}$$ (red lines)]. Dashed lines mark the mirror image (in Q and q) of the periodic orbits drawn with solid lines. The mirror images are also periodic orbits due to the parity symmetry of the Hamiltonian. In the projected Husimi distributions, lighter colors indicate higher concentrations, while black corresponds to zero. The values of the energy ϵ_k and localization measure $$ {\mathfrak{L}}_k$$ of each eigenstate k are indicated in the panels. Energies larger than −0.8 are in the chaotic region. The system size is j = 30.

It is evident from Fig. 1 that the degree of delocalization of the eigenstates in phase space also varies. The Husimi distribution of the eigenstates in Fig. 1a5, a6, for instance, is not entirely confined to the two periodic orbits drawn in blue. This contrasts with the high density concentration that the eigenstate in Fig. 1a4 shows around the plotted unstable periodic orbits. To quantify these differences, we introduce a measure of the degree of localization of a quantum state in the classical energy shell.

Scarring vs. phase-space localization

To measure the localization of a state in a Hilbert-space basis indexed by some letter n, the most commonly used quantity is the participation ratio P_R^46–48. Its inverse is given by $$ P_R^{-1}={\sum }_n{\mathcal{P}}_n^2$$, where $$ {\mathcal{P}}_n$$ is the probability of finding the state in the n’th basis vector. By analogy, we introduce a measure of localization in phase space that employs as basis the overcomplete set of coherent states within a single energy shell $$ {\mathcal{M}}_{ϵ}=\left\{\mathbf{x}=\left(q,p;Q,P\right)\mid h_{\mathrm{cl}}\left(\mathbf{x}\right)=ϵ\right\}$$, so that we replace the sum ∑_k with a three-dimensional surface integral $$ {\int }_{{\mathcal{M}}_{ϵ}}\mathrm{d}\mathbf{s}$$ over $$ {\mathcal{M}}_{ϵ}$$. For a given $$ \mathbf{x}\in {\mathcal{M}}_{ϵ}$$, the probability $$ {\mathcal{P}}_{\mathbf{x}}$$ of finding the state $$ \widehat{\rho }$$ in the coherent state $$ ∣\mathbf{x}⟩$$ is given by the Husimi function $$ {\mathcal{Q}}_{\widehat{\rho }}\left(\mathbf{x}\right)$$. With these replacements, we finally obtain a measure of phase-space localization $$ \mathfrak{L}\left(ϵ,\widehat{\rho }\right)$$ given by

The measure $$ \mathfrak{L}\left(ϵ,\widehat{\rho }\right)$$ is an energy-restricted second moment of the Husimi function⁴⁹. It is related to the second-order Rényi-Wehrl entropy⁴⁷, which, in turn, was shown for the Dicke model⁵⁰ to be linearly related to the first-order Rényi-Wehrl entropy⁵¹.

The value of $$ \mathfrak{L}\left(ϵ,\widehat{\rho }\right)$$ indicates the fraction of the classical energy shell at ϵ that is covered by the state $$ \widehat{\rho }$$. It varies from its minimum value $$ \mathfrak{L}\left(ϵ,\widehat{\rho }\right)~{\left(2\pi {\hslash }_{\mathrm{eff}}\right)}^2/\mathcal{V}\left(ϵ\right)$$, which indicates maximum localization, to $$ \mathfrak{L}\left(ϵ,\widehat{\rho }\right)=1$$, which implies complete delocalization over the energy shell. The former occurs for coherent states, and the latter happens if $$ {\mathcal{Q}}_{\widehat{\rho }}\left(\mathbf{x}\right)$$ is a constant for all $$ \mathbf{x}\in {\mathcal{M}}_{ϵ}$$, in which case the projection $$ {\tilde{\mathcal{Q}}}_{ϵ,\widehat{\rho }}\left(Q,P\right)$$ is also constant for the allowed values of Q and P.

All eigenstates in Fig. 1 have values of $$ {\mathfrak{L}}_k=\mathfrak{L}\left({ϵ}_k,{\widehat{\rho }}_k\right)$$ below 1/2. For the eigenstates in Fig. 1a1, a2, a4, b3, these values are very small, since the eigenstates are almost entirely localized around the plotted periodic orbits. The value of $$ {\mathfrak{L}}_k$$ is larger in Fig. 1a3, because at the center of the diagram, there is an unstable stationary point², which produces a one-point scar in addition to the scar associated to the blue orbit. The localization measure is larger in Fig. 1b1, b2 simply because in these cases the phase space is very small. It is also larger for the states in the high energy region in Fig. 1a5, a6, b4, b5, b6, because they spread beyond the marked periodic orbits. As the energy increases and one approaches the chaotic region, more unstable periodic orbits emerge in the classical limit, enhancing the likelihood that a single quantum eigenstate gets scarred by different periodic orbits. We stress, however, that even for those high-energy states with larger values of $$ {\mathfrak{L}}_k$$, the drawn periodic orbits cast a bright green shadow that is clearly visible in the Husimi projections.

It is important to make it clear that scarring and localization, despite related, are not synonyms. Of course, there is no eigenstate with a large value of $$ {\mathfrak{L}}_k$$ that would at the same time have a high value of our scarring measure, but there is more to the relationship between these two concepts. A highly localized eigenstate is scarred by few periodic orbits of a particular family and therefore has a high value of the scarring measure for that family, although it has low values for the other families of periodic orbits⁴⁵. Furthermore, even the most delocalized eigenstates are still scarred, but now by periodic orbits from different families²².

In Fig. 2, we take a step further in the analysis of localization and scarring. In the large panel in Fig. 2a, we show $$ {\mathfrak{L}}_k$$ against energy for all eigenstates between ϵ_GS and ϵ = 0.06. This plot is equivalent to a Peres lattice⁵² for expectation values of observables, as used in studies of chaos and thermalization. In the low-energy regular regime, $$ {\mathfrak{L}}_k$$ is organized along lines that can be classified with quasi-integrals of motion linked with classical periodic orbits⁵³. Conversely, as the system enters the chaotic region at higher energies, the distribution of $$ {\mathfrak{L}}_k$$ becomes dense and looses any order. Notice, however, that all eigenstates in the chaotic region have values of $$ {\mathfrak{L}}_k$$ much lower than 1, mostly clustering below 1/2.

Fig. 2

Husimi projection and localization measure of eigenstates and random states.

a Localization measure $$ {\mathfrak{L}}_k$$ as a function of energy for the first 17,000 eigenstates with ϵ_k ∈ [ϵ_GS, 0.1], j = 100. (s1–s22) Projected Husimi distributions for 22 eigenstates selected every 400 values of k, starting at k = 7700 in (s1) up to k = 16500 in (s22), j = 100. Lighter colors indicate higher concentrations, while black corresponds to zero. (r1–r4) Projected Husimi distributions for random states of energy width σ = 0.3 centered at ϵ = −0.6 in (r1), ϵ = −0.4 in (r2), ϵ = −0.2 in (r3), and ϵ = −0.1 in (r4), j = 100. Lighter colors indicate higher concentrations, while black corresponds to zero. b Distribution of $$ \mathfrak{L}$$ for 20,000 random states of width σ = 0.3 centered at ϵ = −0.5. c Distribution of $$ {\mathfrak{L}}_k$$ for the eigenstates in the chaotic region with ϵ_k ∈ [ − 0.8, 0]. The system sizes for b and c are indicated in the panels.

The value $$ \mathfrak{L}~1/2$$ marks a limit on the spreading of any pure state in the high energy shells of the phase space. To show this, we build random states $$ ∣{\mathcal{R}}_{ϵ}⟩={\sum }_k{c}_k∣E_k⟩$$, where c_k are complex random numbers from a Gaussian energy profile centered at energy ϵ (see Methods). The values of the localization measure $$ \mathfrak{L}\left(ϵ,{\widehat{\rho }}_{ϵ}\right)$$ for four different random states $$ {\widehat{\rho }}_{ϵ}=∣{\mathcal{R}}_{ϵ}⟩⟨{\mathcal{R}}_{ϵ}∣$$ centered at increasing energies ϵ between −0.6 and −0.1 are given in Fig. 2r1-r4, and indeed $$ \mathfrak{L}~1/2$$. This upper bound on the phase-space delocalization is not due to quantum scarring, but to quantum interference effects.

The panels in Fig. 2r1–r4 display the projected Husimi distributions $$ {\tilde{\mathcal{Q}}}_{ϵ,{\widehat{\rho }}_{ϵ}}\left(Q,P\right)$$ of the four random states and those in Fig. 2s1-s22 show the distributions for 22 eigenstates taken at fixed steps of k with ϵ_k ∈ [−0.62, 0.06] (various other examples are provided in the Supplementary Note 2). The difference between random states and eigenstates is clear. The Husimi projections of the random states do not show structures that resemble closed periodic orbits, so they are not scarred, while the Husimi projections of all eigenstates in the chaotic region do show those structures. In contrast to Fig. 1 (see also ref. ⁴⁵), we have not identified the periodic orbits associated with Fig. 2s1-s22, but the visible circular patterns are clear evidence of periodic orbits. Their existence is supported by the shape of the Husimi projections and by knowing the generic direction of the classical Hamiltonian flow. The patterns display all the features of periodic orbits: they always cross the line P = 0 perpendicularly, they are symmetric along both the horizontal Q and P axes, and they visibly form closed loops. There is no quantum effect other than scarring that would produce such patterns. One therefore deduces that deep in the chaotic region of the Dicke model, all quantum eigenstates are scarred, although they have different degrees of scarring, as seen by the patterns, and they have different degrees of localization, ranging from strong localization ($$ \mathfrak{L}~0.1$$) to states that are nearly as much delocalized as random states.

Our results are sharpened as one approaches the semiclassical limit. The patterns indicating scarring do not fade away as the system size increases. Quite the opposite, as j = 1/ℏ_eff increases, the periodic orbits get better defined in the Husimi projections (cf. the figures for j = 30 and j = 100 in the Supplementary Note 3). To study the dependence of the level of phase-space localization on system size, we show the distributions of $$ \mathfrak{L}$$ for random states (Fig. 2b) and for eigenstates in the chaotic region (Fig. 2c) for different values of j. For the random states, $$ \mathfrak{L}$$ is concentrated around 1/2 and the width of the distribution decreases as j increases, corroborating that $$ \mathfrak{L}=1/2$$ is indeed the delocalization upper bound. In contrast, for the eigenstates in the chaotic regime, the distributions are skewed and broader. The tail at small values of $$ {\mathfrak{L}}_k$$ does not change as j increases, showing that the highly localized states persist, while the portion of the states with large $$ {\mathfrak{L}}_k$$ decreases, suggesting that for large system sizes, none of the eigenstates reach values of $$ {\mathfrak{L}}_k>1/2$$.

The ubiquitous scarring revealed by our studies motivates the question of whether scarring in other quantum models is also the rule. We have found hints in the literature suggesting that our findings may actually be quite general. For example, in ref. ²², the authors reconstruct significant portions of the spectrum of a quantum chaotic system using only periodic orbits. This means that all of these eigenstates are described by those periodic orbits and are therefore scarred. In ref. ²⁰, the authors claim that the great majority of the eigenstates of the hydrogen atom in a magnetic field may be related to periodic orbits, indicating that scars must be the rule. But to provide a definite answer, a phase-space analysis similar to the one presented here is needed.

Quantum ergodicity

We have so far discussed two concepts – quantum scarring and phase-space localization – that are related, but are not equal. How about their relationship with quantum ergodicity? In the classical limit, a system is ergodic if the trajectories cover the energy shell homogeneously. We then adopt the same definition for quantum ergodicity. To quantify how much of the energy shell is visited on average by the evolved state $$ \widehat{\rho }\left(t\right)=e^{-i{\widehat{H}}_{D}t}\widehat{\rho }{e}^{i{\widehat{H}}_{D}t}$$, we consider the infinite-time average^54,55,

How about non-stationary states, such as coherent states or random states, are they ergodic? We study the evolution of initial coherent states $$ ∣\mathrm{\Psi }\left(0\right)⟩=∣{\mathbf{x}}_0⟩={\sum }_k{c}_k∣E_k⟩$$ with mean energies ϵ = −0.5, that are in the chaotic region (see Methods). We select both coherent states that are highly localized and delocalized in the energy eigenbasis, with the degree of delocalization measured by the participation ratio $$ P_R={\left({\sum }_k\mid c_k{\mid }^4\right)}^{-1}$$. Their energy distributions are shown in Fig. 3a1–g1. The components of the states with low P_R are bunched around specific energy levels (Fig. 3a1, b1), exhibiting the comb-like structure typical of scarred states⁶. As P_R increases, the coherent states become more spread in the energy eigenbasis, looking more similar to the random state $$ ∣\mathrm{\Psi }\left(0\right)⟩=∣{\mathcal{R}}_{ϵ}⟩$$ shown in Fig. 3h1, whose mean energy is also ϵ = − 0.5.

Fig. 3

Dynamical behavior of coherent and random states.

a1–g1 Energy distribution for coherent states with different values of P_R and centered at ϵ = −0.5, j = 100. h1 Energy distribution for a random state $$ ∣{\mathcal{R}}_{ϵ}⟩$$ with energy width σ = 0.3 centered at ϵ = −0.5, j = 100. a2–h2 Quantum survival probability (gray solid line), its running average (blue solid line) and its asymptotic value (black dashed line) for the corresponding states of panels a1–h1. a3–h3 Projected Husimi distributions of the time-averaged ensemble, $$ {\tilde{\mathcal{Q}}}_{ϵ=-0.5,\overline{\rho }}$$ for the corresponding states of panels a1–h1 with the values of $$ \overline{\mathfrak{L}}\left(ϵ=-0.5,\widehat{\rho }\right)$$ indicated. i Distribution of $$ \overline{\mathfrak{L}}\left(ϵ,\widehat{\rho }\right)$$ for a set of 1551 coherent states that are evenly distributed along the atomic variables (Q, P) at ϵ = −0.5 and p = 0. This is done for system sizes j = 20, 30, 40, 50, 100, as indicated. j Distribution of $$ \overline{\mathfrak{L}}\left(ϵ,\widehat{\rho }\right)$$ for 500 random states with ϵ = −0.5 and σ = 0.3 for j = 20, 30, 40, 50, 100.

The evolution of the survival probability, $$ S_P\left(t\right)=\mid ⟨\mathrm{\Psi }\left(0\right)∣e^{-i{\widehat{H}}_{D}t}∣\mathrm{\Psi }\left(0\right)⟩{\mid }^2$$, for the coherent states with low P_R leads to large revivals before the saturation of the dynamics⁵⁶, as seen in Fig. 3a2, b2, c2. This contrasts with the evolution of the coherent states with large P_R, such as those in Fig. 3d2–g2, and the evolution of the random state in Fig. 3h2. In these cases, the approach to the asymptotic value of S_P(t) is much smoother and exhibits the so-called correlation hole, which corresponds to the ramp towards saturation. The correlation hole reflects the presence of correlated eigenvalues and is a quantum signature of classical chaos^57–59.

The values of $$ \overline{\mathfrak{L}}\left(ϵ,\widehat{\rho }\right)$$ for the states in Fig. 3a1–h1 are indicated in Fig. 3a3–h3. The random state is indeed ergodic, reaching $$ \overline{\mathfrak{L}}~1$$. For the coherent states, $$ \overline{\mathfrak{L}}$$ increases as P_R does, but even for the states with the largest values of the participation ratio, $$ \overline{\mathfrak{L}}$$ is still slightly under 1. The analysis of the dependence of $$ \overline{\mathfrak{L}}$$ on system size is done in Fig. 3i, j. The distribution of $$ \overline{\mathfrak{L}}$$ for random states (Fig. 3j) gets narrower and better centered at $$ \overline{\mathfrak{L}}=1$$ as j increases, confirming that these states indeed behave ergodically. The distribution for the coherent states (Fig. 3h) is much broader. The center moves towards larger values as j increases, but it is not clear whether there will ever be a significant portion of the initial coherent states with $$ \overline{\mathfrak{L}}~1$$.

In Fig. 3a3–h3, we plot the projected time-averaged Husimi distributions $$ {\tilde{\mathcal{Q}}}_{ϵ,\overline{\rho }}\left(Q,P\right)$$ for the states in Fig. 3a1–h1. Remarkably, even for the coherent states with high P_R, which do not exhibit any comb-like structure in their energy distributions and do not show revivals in the evolution of their survival probability, we still see an enhancement around unstable periodic orbits, as revealed by a careful inspection of Fig. 3d3–g3 and in contrast with the absence of any pattern for the random state in Fig. 3h3. This unexpected manifestation of dynamical scarring can only be observed in phase space, having no identifiable signature in the energy distribution of the initial states or in the evolution of the survival probability. Thus, revivals in the long-time quantum dynamics are signs of a scarred initial state, but lack of revivals do not exclude the presence of scarring, it just indicates that if it exists it is at a low degree.

Classical Hamiltonian

To construct the classical Hamiltonian in a four-dimensional phase space $$ \mathcal{M}$$ with coordinates x = (q, p; Q, P), we use the Glauber-Bloch coherent states $$ ∣\mathbf{x}⟩=∣q,p⟩\otimes ∣Q,P⟩$$^{30,34,63–66}. They are tensor products of the bosonic Glauber coherent states $$ ∣q,p⟩=e^{-\left(j/4\right)\left(q^2+p^2\right)}{e}^{\left[\sqrt{j/2}\left(q+ip\right)\right]{\widehat{a}}^{†}}∣0⟩$$ and the pseudo-spin Bloch coherent states $$ ∣Q,P⟩={\left(1-\frac{Z^2}{4}\right)}^j{e}^{\left[\left(Q+iP\right)/\sqrt{4-Z^2}\right]{\widehat{J}}_{+}}∣j,-j⟩$$, where Z² = Q² + P², $$ ∣0⟩$$ denotes the photon vacuum, $$ ∣j,-j⟩$$ designates the state with all atoms in their ground state, and $$ {\widehat{J}}_{+}$$ is the raising atomic operator. The classical Hamiltonian is given by

Efficient basis and system sizes

The efficient basis is the Dicke Hamiltonian (1) eigenbasis in the limit ω₀ → 0, which can be analytically obtained by a displacement of the bosonic operator $$ \widehat{A}=\widehat{a}+\left(2\gamma /\left(\omega \sqrt{\mathcal{N}}\right)\right){\widehat{J}}_x$$ and a rotation of − π/2 around the y axis of the collective pseudo-spin operators $$ \left({\widehat{J}}_x,{\widehat{J}}_y,{\widehat{J}}_z\right)\to \left({\widehat{J}}_z^{\prime },{\widehat{J}}_y^{\prime },-{\widehat{J}}_x^{\prime }\right)$$,

This basis allows to work with larger values of j by reducing the value of $$ N_{\max }$$ required for convergence of the high-energy eigenstates. With j = 100 and $$ N_{\max }=300$$ (D = 60,501), we are able to get 30,825 converged eigenstates covering the whole energy spectrum up to ϵ = 0.853. Having converged eigenstates in such a high-energy regime would be infeasible with the usual Fock basis for j = 100^44,67–69.

Husimi projection and localization measure

To compute the Husimi projection given in Eq. (2) and the localization measure given by Eq. (3), one has to compute integrals of the form

By using properties of the δ function, those integrals are reduced to

It is worth noting that a quantum state with a Wigner distribution that is constant in an energy shell will lead to a Husimi function that needs not to be constant within the same energy shell. This is because the Husimi distribution is the convolution of the Wigner distribution with the Gaussian Wigner distributions of the coherent states, which have different energy widths due to the geometry of the energy shells in the phase space. This rather marginal effect may be seen in Fig. 3h3, where there is a barely visible weak concentration towards the center of the plot causing $$ \overline{\mathfrak{L}}$$ to be slightly under 1. We stress that this effect is not related to quantum scarring. It is just a manifestation of the phase-space geometry in the Husimi distributions.

Coherent states

The coherent states $$ ∣\mathbf{x}⟩$$ are described by four coordinates x = (q, p; Q, P). To select coherent states at a given energy ϵ, we solve the second-degree equation h_cl(q, p; Q, P) = ϵ for q, which yields two solutions $$ q_{+}\left(ϵ,p,Q,P\right)\ge q_{-}\left(ϵ,p,Q,P\right)$$⁵⁹. All of the initial coherent states shown in Fig. 3 have bosonic variables given by $$ q=q_{+}\left(ϵ=-0.5,p,Q,P\right)$$ and p = 0. The coherent states shown in Fig. 3a1-g1 have atomic coordinates given by (Q, P) = (1.75, 0) a1, (0.5, 0.75) b1, (0.75, 0.5) c1, (1.25, 0.25) d1, (−1.25, 1) e1, (−1.25, 0.75) f1, and (−0.75, 0.5) g1. The atomic coordinates of the 1551 coherent states whose distributions of $$ \overline{\mathfrak{L}}$$ are shown in Fig. 3i were selected by constructing a rectangular grid with step ΔQ = ΔP = 0.05 from $$ \left(Qi,Pi\right)=\left(2,2\right)\left(Qi,Pi\right)=\left(2,2\right)$$ to (Q_f, P_f) = (2, 0). Of the 3321 points inside of this grid, 1551 have allowed values of Q, P (i.e. they satisfy Q² + P² ≤ 4) and fall inside of the energy shell at ϵ = −0.5 (i.e. there exists $$ q_{+}\left(ϵ=-0.5,p=0,Q,P\right)$$). We use these 1551 coherent states to compute the distributions in Fig. 3i.

Random states

The state $$ ∣{\mathcal{R}}_{ϵ}⟩={\sum }_k{c}_k∣E_k⟩$$ is built by sampling random numbers r_k > 0 from an exponential distribution λe^−λx and random phases θ_k from a uniform distribution in [0, 2π). We use