Selective state spectroscopy and multifractality in disordered Bose-Einstein condensates: a numerical study

We propose to apply a modified version of the excitation scheme introduced by Volchkov et al. on bosons experiencing hyperfine state dependent disorder to address the critical state at the mobility edge of the Anderson localization transition, and to observe its intriguing multifractal structure. An optimally designed, spatially focused external radio frequency pulse can be applied to generate transitions to eigenstates in a narrow energy window close to the mobility edge, where critical scaling and multifractality emerge. Alternatively, two-photon laser scanning microscopy is proposed to address individual localized states even close to the transition. The projected image of the cloud is shown to inherit multifractality and to display universal density correlations. Interactions – unavoidably present – are taken into account by solving the Gross-Pitaevskii equations, and their destructive effect on the spectral resolution and the multifractal spectrum is analyzed. Time of flight images of the excited states are predicted to show interference fringes in the localized phase, while they allow one to map equal energy surfaces deep in the metallic phase.


Introduction
Anderson localization is one of the most fundamental quantum interference phenomena in disordered quantum systems.As first pointed out in the seminal work of Anderson 1 , all eigenstates of a particle on a disordered lattice localize in space in sufficiently strong disorder.Though this statement is independent of dimensionality, the existence of delocalized states and the nature of eigenstates at various energies still depend on the spatial dimension of the system considered, and may also be affected by the specific structure of the disorder [2][3][4] .While in d = 1 dimension all eigenstates are always localized for uncorrelated diagonal disorder, in dimensions d > 2 a critical disorder strength exists 5,6 .Below this critical disorder, mobility edges E mob separating localized and delocalized states emerge.At the critical energies, E = E mob , a peculiar quantum phase transition takes place.Approaching E mob through localized states, the typical spatial extension ξ of the wave functions is found to diverge as a power law, ξ (E) ∼ |E − E mob | −ν , with ν a universal exponent that depends only on dimensionality and the symmetry of the underlying Hamiltonian 7 .In case of d = 3 spatial dimensions and in the absence of external magnetic fields and spin-orbit interaction, studied here, the exponent ν = ν orth = 1.58 . . .has been determined with great accuracy by transfer matrix methods, with the label 'orth' referring to orthogonal universality class 8 .At the critical energies, E = E mob , eigenstates of universal properties emerge; the absence of length scale, ξ → ∞ implies a self-similar character and entails a multifractal structure 9,10 of the wave function [11][12][13][14] .
In particular, at criticality, the probability distribution associated with the critical wave function ψ r is arXiv:1709.08993v1[cond-mat.dis-nn]26 Sep 2017 predicted to scale with system size as with f (α) the universal multifractal spectrum of the critical wave function and L the system's linear extension 15,16 .This multifractal structure has a deep origin: it is interpreted as a signature of the presence of infinitely many relevant operators at the Anderson transition 15 .Remarkably, multifractality is also suggested to emerge in the context of many-body localization: recent numerical studies of the Anderson problem on a random regular graph (RRG) reported the existence of the non-ergodic phase with eigenfunctions exhibiting multifractal structures with disorder dependent fractal dimensions [17][18][19] .Though Anderson's localization transition has been experimentally studied in a range of systems including doped semiconductors 20 , granular optical media [21][22][23][24] , disordered microwave waveguides 25 , elastic networks 26 , photonic lattices 27,28 , and cold atomic systems [29][30][31][32][33] , in spite of the considerable effort [34][35][36] , a convincing measurement of the critical state's predicted multifractal spectrum remains missing.
In conventional solid state systems such as doped semiconductors, e.g., Scanning Tunneling Microscopy (STM) provides a way to gain information on the detailed structure of the critical wave function 35,36 .Indeed, in doped ferromagnetic semiconductors, emergent short distance power law signatures have been demonstrated close to criticality at the Fermi energy 35 , but the predicted anomalous exponent has not been observed and the computed multifractal spectra remained far from the theoretically predicted universal form.Rather, measured amplitude distributions exhibited close to lognormal distributions, characteristic of disordered metals.Moreover, in these systems, strong electron-electron interactions have a large impact on the metal-insulator transition, and Anderson's non-interacting theory looses its validity 37,38 .STM signatures of multifractal properties have also been reported on iron-doped InAs surfaces in the context of quantum-Hall transition 36 , but the energy resolution was insufficient to observe the mulifractal spectrum of the critical quantum-Hall state 39 .Indications of multicriticality have also been reported in ultrasound experiments 34 , but to our knowledge the critical mulifractal spectrum has never been observed in these experiments, either.
Experiments with ultracold atoms provide powerful and versatile means to realize and study the Anderson transition.In these experiments laser speckles 31,32 , quasiperiodic optical lattices 30 , or holographic methods [40][41][42] can be used to produce disorder, and the interaction of the atoms can be controlled by changing the depth of an applied optical lattice or the strength of the confinement potential or using Feschbach resonances.Localization transition has been observed in one-and three-dimensional systems [29][30][31][32][33] , and weak localization and coherent back scattering in two dimensions has also been investigated experimentally 43,44 .None of these experiments could, however, detect the critical state so far.In fact, detecting the critical state turns out to be a very delicate task: importantly, the multifractal behavior is the property of a single eigenstate and, as we shall see, admixture of just a few eigenstates by interactions or inelastic processes even in a narrow energy window makes it hard to observe.
Here we propose to image the critical state and to detect its multifractal properties in peculiar cold atomic setups, sketched in Fig. 1.Following the method proposed and demonstrated by Volchkov et al 45 , we suggest using bosons with two hyperfine components, σ =↑ and ↓, with only the component ↑ being subject to a random potential V ↑ .Placing all bosons initially in the weakly interacting (and homogeneous) state ↓ and exciting from this spin ↓ condensate spin ↑ states with the external frequency ω, the energy E ↑ = E ↓ + hω of the final state can be selected and tuned through the mobility edge of the spin ↑ component either by changing ω or by varying the disorder amplitude.As we demonstrate through detailed large scale simulations, it is possible to single out only a few final states in a very narrow energy window by (1) carefully optimizing the shape of the excitation pulse, (2) switching off interactions in the final state, and (3) by using carefully designed excitation geometries.
In particular, we propose to use two different arrangements.In the first geometry, we propose to use a thin vertical laser beam with two frequencies to produce two-photon Raman transitions to localized states in a narrow spatial range and at a well-defined energy such that that the excited wave functions do not overlap.Projection imaging along the horizontal direction (see Fig. 1c) yields images of non-overlapping localized wave functions close to the mobility edge E mob .Alternatively, instead of one thin vertical beam with two frequencies, one could use crossed laser beams to address individual localized states via two-photon Raman transitions in just the small crossing region (see Fig. 1d).This method, developed in biological microscopy and known as 'Two Photon Laser Scanning Microscopy' (TPLSM) 46 , is similar in spirit to the one proposed by Kollath et al. 47 , but here no probe atoms are needed.Rather, transitions are generated from a Bose condensate of atoms, and the local character of the transitions is ensured through crossing the two laser beams.
As we show below through detailed simulations, multifractal properties of the critical state are inherited by the projected pictures of localized states with sufficiently large localization lengths.The projection of the squared amplitude of wave function, i.e., the density exhibits a non-trivial multifractal spectrum which we determine through large scale simulations, and the projected density-density correlation functions are predicted and demonstrated to display universal scaling collapse, thereby providing evidence for multifractal behavior.

Results
In the absence of radiation, we can describe the two-component Bose mixture in the optical lattice by the Hamiltonian, with the prime indicating summation over neighboring sites only.The onsite energies ε rσ incorporate the chemical potential, and ε r↑ also contains the random component V ↑ (r), responsible for Anderson localization of the ↑ bosons.In practice, the random potential V ↑ (r) is realized by fine-grained speckle potentials 48 or holographic disorder 42 .Here, for simplicity, we shall replace it by a uniform and independent distribution of ε r↑ ∈ [−W /2,W /2] on each lattice site.A possibility to generate hyperfine ↑ ⇒ ↓ transitions in a sufficiently narrow spatial region is to use two-photon processes 46 .Here we consider a stimulated Raman process 45 , which generates an effective RF field with a frequency corresponding to the beating frequency of the two lasers, and can be described in the rotating frame approximation by the coupling term with Ω r (t) the Rabi frequency, and the RF frequency appearing simply as an energy shift ε ↑,↓ → ε ↑,↓ ± hω/2.For Ω r (t) we assume a Gaussian profile, Ω r ∼ Ω 0 (t)e −2(x 2 +y 2 )/w 0 2 with a narrow waist w 0 in the range of a few lattice constants, w 0 ∼ 2 µm The shape of the pulse Ω 0 (t) must be determined carefully to generate transitions between the two hyperfine components in an energy window as narrow as possible.
Since the coupling happens between a discrete state and a quasi-continuum, the width of the final energy window follows from Fermi's Golden rule 49 , and can be adjusted by the duration and shape of the pulse Ω(t).However, interactions in the final state turn out to be particularly destructive: they generate unwanted transitions and broaden the spectrum of the final state.Therefore, in the following, we shall assume that interactions in the final state have been switched off by applying an appropriate external magnetic field, U ↑↑=0 .We assume further that the two hyperfine states are decoupled for t < 0, Ω 0 (t < 0) = 0, and the bosons form a condensate in the state ↓, described by a collective wave function a r,σ → ψ rσ at time t = 0. We describe the time dependence of this collective wave function in terms of the Gross-Pitaevskii equation (see Methods).

Numerical simulations.
We determined the time dependence of the collective wave function by solving Eq. ( 10) numerically.To maintain numerical accuracy, we suppressed fast phase oscillations of the wave function by means of a self-consistently determined global gauge transformation (see Eq. ( 12) in Methods).This allowed us to perform Gross-Pitaevskii simulations on systems as large as 40 × 40 × 40, where length is measured in units of the lattice constant a.Since the system size is not larger than the Rayleigh lenght πw 2 0 /λ of the laser beam, we neglected the spatial spreading of the waist, and assumed a simple cylindrical beam.The RF field Ω r (t) was turned on smoothly, kept constant and then turned off with an optimized pulse of length T = 400 h/J and turn on time τ = 60 h/J.This resulted in an energy window ∆ε of the final states encompassing a limited number of states (see the discussion below).
In small disordered systems, disorder averaging is necessary to eliminate sample to sample fluctuations.Averaging over a huge ensemble of disorder realizations is experimentally demanding.Disorder averaging can, however, also be replaced in part by averaging over the RF frequency ω in a small window, since the localization length ξ (E) displays only a weak energy dependence near the band center 5 .
The spectral and spatial structure of the state reached after the excitation pulse is shown in Fig. 2 for some typical parameters.As shown in panel (a), states are excited in a very narrow energy window.To determine the spectral weights w ↑ (E) and the density of states, displayed in panel (a), we employed the so-called Chebyshev polynomial expansion (see Methods), which made possible to obtain high accuracy results and thus circumvent the impossible task of performing a complete diagonalization of the final Hamiltonian.Although states in a very narrow window of width ∆ε ∼ 0.05 J could be excited by means of optimizing the pulse shape and using a small concentration of atoms, still, around 100 eigenstates can be found in this tiny interval even in a relatively small 40 × 40 × 40 lattice, containing 64,000 lattice sites.Mixing 100 eigenstates would make the multifractal properties completely invisible.Indeed, on the metallic side, states are uniformly excited over the whole sample (see Fig. 2b), and exciting and imaging of individual states does not seem to be possible.Fortunately, as shown in Fig. 2c, this problem can be avoided on the localized side of the transition by applying narrow laser beams of waist w 0 < ξ , and thereby reducing the number of excited eigenstates by a factor ∼ (ξ /L) 2 , where L denotes the linear size of the system.These localized states are excited along the waist of the laser beam, but they do not spatially overlap under the condition with BW standing for the the bandwidth of spin ↑ bosons.To derive this inequality, we notice that a narrow laser beam excites states only in a column of volume ξ 2 L. Consequently the total number of excited states is N ∼ ξ 2 L ∆ε/BW, and the average spatial distance between these states is d ∼ L/N.The overlap of the distinct excited states in the image is thus small if d > ξ , which leads us to the condition Eq. ( 4).This condition -meaning that the energy resolution be better than the level spacing in a localization volume -can be satisfied even relatively close to the transition.The same condition is obtained for the TPLSM protocol (Fig. 1d).We thus conclude that imaging individual states is possible on the localized side of the transition with the suggested protocols as long as Eq. ( 4) remains satisfied.

Multifractal properties
Statistical properties of the critical state reflect the structure of the localization quantum phase transition.Close to the transition, the disorder-averaged correlation functions The exponent y q here turns out to be intimately related to the Legendre transform τ q of the multifractal spectrum f (α) with α ≡ dτ q dq , f (α) = α q − τ q , and d = 3 15 .The previous considerations refer to three dimensional wave functions and densities, which are quite difficult to determine experimentally 50 .Fortunately however, as we now show, it is also possible to observe the fingerprints of multifractality in the experimentally measured projected densities of the gas, ρxσ ≡ dz|ψ rσ | 2 (7)   with x = (x, y) referring here to coordinates in the plane of projection.While the rare regions of small amplitudes (large α) appear to be washed out by the projection and multifractality is destroyed there, multifractal scaling and power law correlations are found to survive the projection in the large amplitude (small α) regions and are inherited by the projected wave function.For q = 1, in particular, it is easy to prove that the projected densities display power law correlations with an anomalous exponent providing a clear analytical evidence for the survival of multifractality.This relation is indeed verified by our numerics (see online Supplemental Material), which also confirms that the projected densities ρx↑ display multifractal scaling according to Eqs. ( 1) and ( 5) with, however, a modified multifractal spectrum, f (α), and modified exponents ŷq and correlation functions Ĉ(q) (x) for spin ↑ bosons The universal collapse (9) of the correlation functions, as extracted from the horizontally projected image after the excitation pulse is displayed in Fig. 3a Ĉ(q) (x), the collapse onto a universal crossover function and the divergence of the localization length upon approaching W c are both clear.
Density-density correlations are found to be sensitive to the system size.In contrast, the multifractal spectrum characterizing the probability distribution of the projected wave functions seems to be much more robust.We extracted this latter simply by measuring the distribution of ln(ρ x↑ ) after a pulse targeting the critical state.The final results are shown in Fig. 3b.In the absence of interactions, the spectrum of the state reached after the Rabi pulse follows very closely the spectrum of the true critical state itself.The maximum of f (α) characterizes the scaling of the typical wave function amplitudes at criticality, and the fact that α max > 2 is a clear evidence of the fractal nature of the projected wave function.Notice that for α α 3D max − 1 ≈ 3 we do not expect a true multifractal spectrum: large values of α correspond to tiny wave function amplitudes which would be mixed with larger, typical wave function amplitudes upon projection, leading to the destruction of multifractality in the large α region.Indeed, a finite size analysis seems to confirm that the f (α) function is well defined only below a critical value of α.
Having the spectra f (α) (and f (α)) at hand, we can also perform a Legendre transformation to determine the critical exponents τq and ŷq (τ q and y q ), and compare them to the exponents extracted directly from the correlation functions (see online Supplemental Material).In particular, for q = 1 we obtain y 1 ≈ 1.8 and ŷ1 ≈ 0.8, in good agreement with the exponents extracted from the scaling collapse of the density correlation functions, thereby verifying the exact relation (8).
We need to mention two experimental aspects that may reduce the visibility of the multifractality of the excited condensate.1][52] One typically cannot therefore access ρx↑ through a single shot measurement, rather, one has to prepare the same wave function multiple times and average the measured densities to obtain ρx↑ .High density regions, in particular, can be measured accurately, while low density regions are less visible.Fortunately, the multifractal spectrum depends only logarithmically on these statistical fluctuations.Therefore, as shown in the Supplemental Material online, already hundreds of measurements on the same wave function can be sufficient to extract ρx↑ and reveal f (α).Particle number fluctuations may also lead to trouble: such fluctuations influence the background potential and lead to a slight detuning of the final state in each measurement.As shown in the Supplemental Material online, below few percent particle number fluctuations are, however, not detrimental, and do not destroy the multifractal spectrum.We remark that post-selection is a natural possibility to reduce the error induced by particle number fluctuations.

Time of flight images
While time of flight images do not show any particular features at criticality, they do exhibit peculiar structures on both sides of the localization transition.For strongly localized final states, in particular, excited atoms happen to localize at only a few resonant localized states.Releasing the atoms from the trap then leads to interference fringes created by atoms released from these localized states, similar to the fringes observed in the case of split condensates 53 (see Fig. 4).In the delocalized phase, on the other hand, weak disorder leads to the mixing of states on a relatively narrow energy shell, and, correspondingly, "Fermi surface"-like features emerge in the time of flight image, tracing approximately the equienergetic surfaces in the former Brillouin zone.

Discussion
Here we proposed an experimental protocol that makes it possible to image and observe the so far elusive Anderson critical state using weakly interacting ultracold bosons with two hyperfine states.To reveal the critical correlations and amplitude fluctuations, one needs to excite atoms into a few eigenstates in an extremely narrow energy window, otherwise the fragile multifractal structure becomes completely invisible.We performed detailed simulations to verify the feasibility of the proposed experiment.Three important ingredients have been employed to reach the energy resolution needed: (1) interactions in the final hyperfine state were suppressed, (2) the shape of the excitation pulse has been optimized, and (3) excitation beams narrower than the localization length have been combined with a horizontal imaging to reduce the number of excited eigenstates on the localized side of the transition.We demonstrated that combining these three ingredients, one can reach energy resolutions sufficient to observe even single localized states, and to detect critical multifractal correlations as well as multifractal amplitude fluctuations, even after projecting the image of the excited cloud.
We have verified that the projected image inherits the multifractal structure.In particular, we have determined the multifractal spectrum of the image in the case of orthogonal symmetry, and predicted and have demonstrated by our simulations a universal power low scaling for the projected density-density correlations as a clear signature of multifractality.Time of flight images of the excited states have also been shown to display striking features in both phases: strong localization-induced interference fringes appear in the localized phase, while momentum space equal energy surfaces can be imaged on the delocalized side of the transition.
Cold atomic systems and the method proposed here thus provide a unique framework to reach the critical state of Anderson localization and to study its multifractal properties, observed solely numerically so far.The proposed method experiment is, of course, not restricted to Anderson's localization transition.It can also be used to shed light on the critical state of disordered quantum Hall systems 54 , disordered systems with spin-orbit coupling 55,56 , or even more exotic topological phase transitions 57,58 .The crossed beam two-photon laser scanning imaging method would allow one to create and study local excitations.In particular, in the case investigated here it allows to excite single localized eigenstates and to study their spatial structure in detail.

Gross-Pitaevski simulation
The time evolution of the collective wave function is described by the mean field Gross-Pitaevskii equation, Substituting the Hamiltonian (Eq.( 2)) in the above expression we arrive at the equations The equation of motion for Φ rσ (t) is then solved numerically by using a standard 4th-order Runge-Kutta scheme.

Chebyshev method for spectral properties
To get spectral information about the final state Φ r↑ (t end ) the kernel polynomial method is used 59 .At the end of the excitation process the external RF field is turned off, Ω r (t end ) = 0, and the mean-field Hamiltonian has the diagonal form The single particle density of states for the ↑ bosons is approximated by the series where E α denotes the energy of the eigenstate |α of Ĥ↑ (ψ), T n (x) stands for the n'th Chebyshev polynomial and µ n denote the Chebyshev expansion coefficients, which can be efficiently determined using the recursion relation for the Chebyshev polynomials 59 .The spectrum is transformed to the [−1, 1] interval by shifting energies to the band center, BC and normalizing them by the bandwidth BW of the mean-field Hamiltonian.The factors g (n max ) n are the expansion coefficients of the so called kernel function 59 that regularizes the finite order approximation of the δ -functions.In our calculations we used the Lorentz kernel g (n max ) n = sinh(λ (1 − n/n max ))/ sinh(λ ) that guarantees a positive definite density of states 59 .The parameter of the kernel was set to λ = 0.1.
A similar calculation is performed to approximate the spectral distribution of the final state of the Gross-Pitaevskii simulation yielding with Φ ↑ denoting the wave function of the ↑ bosons in the final state.The cutoff was set to values as large as n max = 30, 000 to resolve the sharp spectral peak of the final state.

Figure 1 .
Figure 1.(a) Projected image of a critical multifractal eigenstate of size (L/a) 3 = 40 3 , exhibiting clustering.(b-c) Schematics of the proposed experiment: (b) The condensate is prepared in the hyperfine state |↓ where no disorder is present.The external field excites the atoms to the hyperfine state |↑ , interacting with the disorder potential V ↑ .The energy E ↑ of the final state is controlled by the frequency ω.(c) Excited atoms in the |↑ state are imaged by a horizontal laser beam.The overlap between different eigenstates can be reduced by using a narrow excitation beams of waist w 0 ξ , where ξ is the typical size of a localized state.(d) Even a single localized state can be excited and imaged by crossing two laser beams.

Figure 2 .
Figure 2. (a) Spectral weight w ↑ (E ↑ ) of the bosons as a function of the final state's energy E ↑ after the excitation pulse for a lattice of (L/a) 3 = 40 3 sites, W = 17J and excitation beam waist w 0 = 4a.We simulated N = 3200 atoms with typical interaction strengths U ↑↓ = U ↓↓ = 2J, excitation frequency ω ≈ 8J/h.The grey area indicates the density of states.Inset: same spectral peak on a smaller scale.(bc) Projected images of a metallic (W = 10J) and a localized (W = 28J) final state.The excitation beam has vertical direction in both cases.(d) Projected image of a single localized state, excited by two crossed laser beams.All figures show images of size 40a × 40a.

Figure 3 .
Figure 3. (a) Rescaled correlation function of the projected densities ρx↑ for different disorder strengths but fixed excitation frequency ω ≈ 8J/h in the localized regime.Correlations were measured along the y direction in Fig. 1.Other parameters were set as in Fig. 2. Values of ξ (W ) were determined by collapsing of the curves.Inset: divergence of ξ for W → W c ≈ 16.5J.(b) Multifractal spectrum extracted from the probability distribution of the projected density ρ↑ .The measured spectrum after a pulse for non-interacting bosons (red circles) follows closely the spectrum obtained from the critical state (thick green line).Interactions broaden the spectrum of the final state (see inset) and slightly deform the spectrum f (α) (blue triangles).The location of the maximum, α max ≈ 2.3, characterizes the anomalous scaling of the typical density, ρtyp x↑ ∼ L −α max .

/ 4 W=25JW=10JFigure 4 .
Figure 4. Time of flight image of excited atoms in the localized (left) and in the delocalized phase (right).In the localized phase interference fringes emerge, while in the delocalized phase a "Fermi surface" structure tracing smeared equal energy surfaces in momentum space appears.