Feedback cooled Bose-Einstein condensation: near and far from equilibrium

Continuously measured interacting quantum systems almost invariably heat, causing loss of quantum coherence. Here, we study Bose-Einstein condensates (BECs) subject to repeated weak measurement of the atomic density and describe several protocols for generating a feedback signal designed to remove excitations created by measurement backaction. We use a stochastic Gross-Pitaevskii equation to model the system dynamics and find that a feedback protocol utilizing momentum dependant gain and filtering can effectively cool both 1D and 2D systems. The performance of these protocols is quantified in terms of the steady state energy, entropy, and condensed fraction. These are the first feedback cooling protocols demonstrated in 2D, and in 1D our optimal protocol reduces the equilibrium energy by more than a factor of 100 as compared with a previous cooling protocol developed using the same methodology. We also use this protocol to quench-cool 1D BECs from non-condensed highly excited states and find that they rapidly condense into a far from equilibrium state with energy orders of magnitude higher than the equilibrium ground state energy for that condensate fraction. We explain this in terms of the near-integrability of our 1D system, whereby efficiently cooled low momentum modes are effectively decoupled from the energetic `reservoir' of the higher momentum modes. We observe that the quench-cooled condensed states can have non-zero integer winding numbers described by quantized supercurrents.


I. INTRODUCTION
The repeated measurement of an interacting closed quantum system typically leads to runaway heating and the ultimate destruction of quantum correlations. In some cases, the interplay of unitary evolution, measurement, and feedback can lead to physics beyond what is possible in closed quantum systems. For example, recent developments in quantum circuits with randomized measurements combined with (sometimes random) entangling operations have shown remarkable robustness of entanglement [1][2][3][4][5][6][7][8][9][10]. Furthermore, many quantum error correction algorithms rely on measurement of ancilla qubits followed by feedback to protect the entanglement between the constituent (logical) qubits [11,12]. To clarify the connection between these ideas and their potential realization in steady-state many body systems, we focus on theoretically modeling Bose-Einstein condensates (BECs) in which feedback conditioned on the history of measurement outcomes provides a mechanism to retain the system in a state of low entropy, low energy, and with long-range coherence.
Many feedback schemes have been proposed to bolster quantum simulation capabilities. Ref. [13] showed that weak measurement combined with classical feedback can produce dynamics described by a non-linear Schrodinger equation. More recent developments have shown the utility of feedback control for creating new phase transitions, simulating spin models and preparing targeted states in optical lattices [14][15][16][17][18]. Furthermore, Refs. [19,20] demonstrated the utility of weak measurement and feedback for generating effective interactions in continuum 1D BECs, with sufficient cooling to achieve a quasi-steady state. The ability to maintain the system in a quasi-steady state with long-range coherence is paramount for experimental realizations of these proposed ideas. To this end, here we describe and optimize cooling protocols for continuously monitored BECs in one and two dimensions.
We build upon a previous framework of weak measurement and classical feedback used to control BECs developed in Refs. [19,20], where the system is modeled by a stochastic Gross-Pitaevskii equation (GPE). These feedback schemes work by measuring the position resolved atomic density and applying a feedback potential that has the observed density distribution as its ground state.
Because the measurement provides no information regarding the condensate phase, this procedure is imperfect. By ignoring the condensate phase, these methods neglected dynamics which are most simply described by phonons. The improved cooling scheme presented here (see Fig. 1) uses signal filtering techniques to create feedback potentials that incorporate these effects into the feedback protocol, thereby making it more effective.
Finally, we explore the utility of weak measurement and feedback cooling in analog quantum simulation by preparing far from equilibrium Bose-condensed states. Using optical measurements as a model for weak measurement implies that there is a resolution limit below which an observer obtains no information about the system, even in principle. We show that 1D systems with significant excitations below this scale can still be feedback-cooled to Bose condensation. In this case, the long-wavelength modes are rapidly cooled while the short wavelength modes are left relatively untouched. Because the 1D system is nearly integrable, there is negligible energy transfer between these groups of modes, yielding a low energy system weakly coupled to a highly excited reservoir. We observe that these far from equilibrium condensed states may have non-zero integer winding numbers, described here by quantized supercurrents.
This manuscript is organized as follows: In Sec. II we outline our simulation methods and define important metrics for exploring condensate properties at finite entropy. The rest of Sec. II outlines our measurement model for phase-contrast imaging, starting from the lightmatter interaction Hamiltonian, and briefly reviews the concepts of signal filtering and feedback cooling essential to the protocols presented here.
Our main results are presented in Sec. III, showing feedback cooling protocols that have been significantly improved by incorporating momentum-dependent feedback. We also present the parameter space for which these protocols are optimal according to the final energy. In Sec. IV, we demonstrate that in 1D our feedback cooling protocols can cool from highly-excited un-condensed states to high energy, far from equilibrium condensed states with non-zero integer winding numbers. In Sec. V we conclude and present some directions for future research. Additional details regarding the measurement model are left to the appendices.

II. STATE OF THE ART
A. GPE with non-zero entropy The usual GPE accurately describes the low temperature properties of atomic BECs, discussed in 1D and 2D in this section. Here, we use the GPE to model BECs coupled by measurement and feedback to the environment: an open quantum system. Inspired by classical field methods [21] to describe finite temperature systems, we begin by studying the generic equilibrium properties of a GPE when used to model a Bose system with nonzero entropy.
In the following sections, we model measurement backaction as a stochastic process introducing number fluctuations conditioned on measurement outcomes [19]. This is reminiscent of stochastic noise terms used in the stochastic-projective GPE, but which are only present while the measurement process is active. Prior to developing our feedback cooling methodology, we identify the equilibrium properties of the GPE by adding a deterministic energy in the form of white noise and allowing the system to equilibrate under energy-conserving closed system evolution. We numerically model this microcanonical ensemble by simulating many trajectories with different noise realizations. Each trajectory can be characterized by its per-particle energy and momentum ground state (k = 0) occupation probability. An ensemble of many trajectories can be assigned a condensed fraction and Von Neumann entropy.

The single component GPE
describes atomic BECs in the mean field limit. Here, Ψ(r, t) is a classical field describing N particles with density n(r, t) = |Ψ(r, t)| 2 at position r and time t; µ is the chemical potential; and m a is the atomic mass. Since the interaction strength u 0 depends on dimensionality, we will use the subscripts u 0,1D and u 0,2D to mark this difference where necessary. For each trajectory, we obtain the total energy an extensive thermodynamic quantity, where the dimension D = 1 or 2. An estimate of the condensed fraction can be obtained from the occupation probability of the k = 0 mode, denoted P (k = 0). We use the ensemble formed by N t trajectories to compute additional thermodynamic quantities. The ensemble is characterized by the one body density matrix where j labels different members of the ensemble. The Von Neumann entropy S = −Tr [ρ logρ] can be directly derived from ρ(r 1 , r 2 ), giving an extensive thermodynamic quantity, where we implicitly set the Boltzmann constant k B = 1. Numerically, the entropy computed by this method is bounded from above by S max = log(N S ), where N S is the smaller of N t or the number of points in our simulation grid, N x . The Penrose-Onsager [22] condensate number N c is the largest eigenvalue of ρ(r 1 , r 2 ), giving another definition of condensed fraction f c = N c /N t . A condensate is present when this eigenvalue is ∼ O(N t ) or f c ∼ O(1). This definition of condensed fraction can be directly extended to fragmented condensates where more than one eigenvalue of ρ(r 1 , r 2 ) is ∼ O(N t ).
Simulation parameters. We numerically simulate the the 1D and 2D GPE using split step Fourier methods [23]. We simulate atomic systems with density n = 10 3 particles per micrometer, a spacing δx ≈ 0.1 µm between grid points, a time step δt = 3.7 µs, and use periodic boundary conditions. We selected 1D and 2D interaction constants u 0,1D u 0,2D to give a chemical potential µ = h×165 Hz, small but reasonable for 87 Rb. The chemical potential µ determines the healing length, where ξ = 2 /2m a µ ≈ 0.6 µm for the system. Energy is presented in units of µ throughout the manuscript.

B. Continuously monitored GPE
Here we briefly review our general formalism describing cold atomic gases [18][19][20] continuously monitored with an off-resonant homodyne technique such as phase contrast imaging [24]. We show that the measurement backaction appears as a random variable in the atomic equation of motion, and we show that the dynamics of the BEC are described by stochastic wave equation. Figure 1(a) shows a schematic of the measurement apparatus for the repeated weak measurement of a 100 µm BEC with periodic boundary conditions in 1D.
The dispersive interaction of far detuned light and the atomic ensemble is described by the second quantized light-matter Hamiltonian where γ describes the coupling strength; D is the dimension;n(r) =Ψ † (r)Ψ(r) measures the D-dimensional atomic density in terms of atomic field operatorsΨ(r); andρ(r, z) =φ † (r, z)φ(r, z) gives the optical mode occupation in terms of optical field operatorsφ(r, z).
We consider a pulse of light incident on our Ddimensional system that therefore travels in a D + 1 dimensional space. Here, the vector r denotes the space containing the atoms and z is the perpendicular dimension in which the light travels. For an incident coherent state with amplitude α [i.e.,φ(r, z) |α = α |α ] and width z = cdt, we define a measurement strength parameter ϕdt 1/2 = γ|α|(2dt/c) 1/2 , a continuum analog to the expression in Ref. [19], where c is the speed of light in a vacuum.
Continuous weak measurements of strength ϕdt 1/2 output a measurement signal where n(r, t) is the atomic density and m(r, t) is the quantum projection noise associated with the measurement. The projection noise is characterized by Fourier domain Gaussian statisticsm k = 0 andm km * k = L D Θ(k 0 −|k|)dW k dW * k /2dt 2 , dW k is a Weiner increment with dW k = 0 and dW k dW * k = dtδ kk [25]. Fourier momenta greater than k 0 = 2π/λ are removed from the measurement signal via the Heaviside function Θ due to the inability of the physical measurement process to resolve information on length scales smaller than λ/2π.
Our protocols control the condensate by applying a feedback potential V [ε], generated using an estimator ε, defined in the following sub-section. The measurement and quantum control protocol is described by a stochastic equation of motion for Ψ(r, t) in D = 1 and 2 dimensions [20]: where are the contributions from closed system evolution, measurement backaction, and feedback, respectively. Figure 1(b-i) shows the middle 20 µm section of the postmeasurement density n of the BEC described by Ψ(r) in the above equations of motion.

C. Signal Filtering
Our simulations model a continuous measurement process, for which the measurement noise in Eq. (5) diverges as dt → 0. In reality, any physical measurement process has a non-zero detection time τ ; accordingly our simulations discretize time into steps of duration δt τ . This yields the discrete stochastic equation in Appendix A which quantifies the measurement strength in terms of κ = ϕ(δt) 1/2 . In general we model a physical detector's response with the convolution expressed in terms of a normalized filter kernel f (t). For example, a single measurement of fixed duration τ would be described by the box filter f box (t) = Θ(t + τ )Θ(−t)/τ . Here, we describe our measurement in terms of the low pass filter described by filter kernel f (t) = Θ(−t) exp(−t/τ )/τ . We term ε(x, t) as the estimator because it is a running average that estimates the true state of the system. Figure 1(b-ii) displays the estimator computed following the measurement process that resulted in the postmeasurement density shown in Fig. 1 This filtered measurement signal encourages us to introduce an optimal measurement strength κ * for which the filtered measurement result is equal to the postmeasurement density n |M (x), i.e. M(x) = n |M (x). As show in App. A, this leads to the optimal measurement strength which is exact only for the box filter. Reference [20] heuristically used this expression to obtain a nearoptimal measurement strength for the exponential filter, thereby coupling κ and τ . In Sec. III, our momentum dependent filtering protocols make τ a function of k and thereby are not described by a single optimal measurement strength. As a result, κ and τ are not coupled in the feedback cooling protocols presented here.

D. Feedback Cooling
Feedback cooling suppresses excitations, including those present prior to any measurements as well as those added by measurement backaction. Our feedback techniques apply a time-dependent feedback potential, derived from the measured density. Even in principle, this feedback cooling is limited because the measurement provides no information regarding the condensate phase, however this limitation can be mitigated with the enhanced feedback cooling techniques discussed later in Sec. III.
To generate a cooling potential V C (x) given a measurement record M(x), we take the stationary GPE for the post-measurement wavefunction ψ |M and, under the Thomas-Fermi (TF) approximation, substitute u 0 n |M → gM and compute the potential for which the ψ |M would be the ground state of the system: Reference [20] found that g = u 0 is the value of gain that results in the best cooling; this was expected as a result of the Thomas-Fermi approximation.
Using an estimator of the form shown in Eq. 10, the potential expression becomes Figure 1(b) shows this basic feedback cooling process using dark green arrows. The arrow between Figures 1(bii) and 2(b-iv) represents the generation of a feedback potential V [ε] according to Eq. (15). The arrow from Figure 1(b-iv) to (b-i) represents the application of V [ε] to the system via the Hamiltonian term in Eq. (9). In Sec. III, we build on the cooling capabilities of this basic feedback protocol by implementing momentum mode-specific feedback.

III. ENHANCED FEEDBACK COOLING
The goal of feedback cooling is to depopulate highmomentum modes (i.e., phonons in our BEC), which are excited during the measurement process due to measurement backaction. The previous cooling protocol [20] outlined in Sec. II D is applied uniformly to all momentum modes. In this reference, a Bogoliubov-de Gennes analysis indicated that momentum-dependent cooling protocols should give improved performance. In this section, we introduce momentum-dependent signal filtering and gain, describe several different protocols that improve feedback cooling, and discuss the merits of each.

A. Momentum-dependent Cooling
The component of the density distribution with momentum k > k 0 is inaccessible to the measurement process and therefore remains untouched by both measurement backaction and feedback cooling. However, the associated modes can become populated via scattering; this contributes to heating and prevents the system from reaching a steady state [20]. An optimal feedback cooling protocol should rapidly depopulate the momentum modes with k k 0 , so that modes with momentum k > k 0 cannot be populated by collisions.
Transitioning the discussion in Secs. II C and II D into Fourier space allows us to straightforwardly introduce momentum-dependent gain and filtering, e.g. g → g(k) and τ → τ (k). This results in a general cooling protocol with k-space estimator giving the feedback potential where F and F −1 denote the Fourier transform and its inverse, and the argument of the function delineate between spatial and spectral representations. The functions τ (k) and g(k) are selected based on physical considerations discussed below. The time constants τ (k) in the integral Eq. (16) used to compute ε(k, t) allow for a momentum-dependent compromise between reduced noise (longer averaging times) and rapid response to system dynamics (shorter averaging times). Here, we used where τ c is the filter time constant at k 0 , the maximum accessible momentum. As a result, the small k contributions to ε(k, t) have relatively longer filter times. We designed this function to be compatible with phonondynamics: a phonon of wavevector k will create a density modulation whose phase changes by 1 radian in a time 1/(ck). Therefore, averaging the information over longer times would not be useful in estimating the current density of the system. Similarly, the gain function g(k) in Eq. (17) includes the possibility that the system's density responds to applied potentials differently on different length scales. For example, low-k phonon excitations couple to density much more weakly than high-k particle-like modes; this motivated us to select a g(k) to decrease with increasing k. This is in contrast with the optimal protocol for adiabatic feedback cooling [20], which features increasing gain at large k. Our momentum-dependent gain function significantly attenuates the gain for momenta larger than a Gaussian width σ. Figure 1(b) displays this enhanced feedback cooling protocol using pink arrows. We show kernel of the estimator integral in Eq. (16) in Figure 1(b-iii).

B. Results
The momentum dependent cooling protocols outlined in Sec. III A cool the BEC more effectively than previous protocols according to several metrics. Figure 2 displays the differences in high-frequency density noise with and without momentum dependent cooling, while Figure 3 summarizes our results quantified in terms of total energy (top), condensed fraction (middle), and entropy (bottom) in both 1D and 2D. We also discuss the optimization of parameters τ C and σ.
The left column of Fig. 3 presents our 1D cooling results. In all cases the system was initialized at t = 0 in the ground state and approached a steady state as the cooling protocols were applied. The pink curves depict the time evolution of the previous best protocol [20], which was momentum independent. We minimized the final energy by tuning the gain g and filter time τ to obtain an optimal momentum independent protocol (blue), reducing the total energy by a factor of about 30. Interestingly, the momentum dependent gain protocol (light green) reduces the steady state energy, but performed worst as quantified by condensed fraction and entropy. The protocols with momentum dependent filtering reduce the total energy by an additional factor about 4 [ Fig. 3 (a)]. The momentum dependent filter (brown) had the most significant effect on cooling across all three metrics, leading to a higher condensed fraction and lower entropy [ Fig. 3 (b-c)]. In contrast, incorporating the momentum dependent gain parameter g(k) (dark cyan) given by Eq. (19) provided little advantage in cooling in terms of the total energy and lead to slightly worse out-comes in terms of the entropy and condensed fraction. This indicates that the physical intuition that shorterwavelength excitations would require reduced feedback strength is either invalid, or that our Gaussian ansatz was poorly motivated. Given the success of an overall gain derived from the Thomas-Fermi approximation, perhaps adjusting the k-dependent gain based upon a higher-order expansion would be better motivated.
With our previous momentum independent protocols, the system did not reach steady state-as quantified by condensed fraction, energy, and entropy-even at long times [20]. By contrast, all of our metrics become time independent for the momentum dependent cooling protocols, indicating that the system is in a steady state.
The 2D cooling results are presented in Fig. 3(d-f). Because there are no literature results for this type of feedback cooling in 2D, we do not plot the previous best protocol. The momentum dependent protocols reduce the total energy by a factor of slightly less than 4 compared to the 2D momentum independent protocol [ Fig. 3 (a)], very similar to the 1D results. However, the overall energy of these protocols is about a factor of 5 higher than their 1D counterparts. By contrast, the un-condensed fraction and entropy (e,f) are about 50 times lower than those achieved by the 1D protocols. We attribute this contrast to an overall decrease in the occupation of excited modes, but with an increase in energy deriving from the increased density of states for higher momentum states. This density of states effect is even evident in individual measurements, where each measurement in 2D adds about 2.5 times more energy than in 1D.
We quantify the overall behavior of the cooling protocols by fitting the time evolution of the energy, condensed fraction, and entropy to fitting functions of the form In this example for energy, E ∞ denotes the asymptotic value of the energy, τ E describes the nominal timescale for equilibration, and α E captures the rapidity of the approach to equilibrium. We optimize the feedback cooling protocols by selecting the gain and filtering parameters that resulted in the lowest steady-state energy E ∞ . The symbols g and τ f denote the momentum independent gain and filtering protocol parameters, respectively, while σ and τ C denote the momentum dependent gain and filtering protocol parameters, respectively. For each set of parameters, we simulate 16 trajectories in 1D and 5 in 2D over a period of 1 second. The results for each protocol are presented in Fig. 4, where a clear local minimum emerges in the final steady-state energy. In Fig. 4 we imposed an upper-limit of E ∞ /µ ∼ 10 1 ; trajectories with energies beyond this limit generally failed to reach equilibrium in 1 s and are therefore not well fit by Eqn. (20).
To reduce the parameter space for the full momentum dependent filtering and gain protocol, we selected the k = 0 gain g in Eq. (19) to be to the optimal gain for momentum dependent filtering. Numeric labels denote parameter scans for the following protocols: i. momentum dependent gain and filtering, ii. momentum dependent gain, iii. momentum dependent filtering, iv. momentum independent gain and filtering.

IV. QUENCH-COOLING INTO FAR FROM EQUILIBRIUM CONDENSED STATES
Section III presented an optimized protocol for continuously measuring a condensate while maintaining it at very low energy and entropy. In this section, we cool high energy, high entropy 1D systems with no initial condensate into final steady states with large condensed fraction and non-zero winding numbers. These integer valued winding numbers are described by quantized supercurrents. These condensed steady states are far from equilibrium, with energy in excess of the "critical" energy for condensation, as discussed below. We show that this results from significant occupation of high momentum states (above the imaging cutoff k 0 ) that, because this is a nearly integrable system, are only weakly coupled to the low momentum modes accessible to feedback cooling.

A. Quasi-thermal excited states
In our simulations, we introduce "thermal" energy by repeatedly adding low amplitude noise over the course of 37 ms and allowing 0.5 s for the system to equilibrate. The noise has a spectrally uniform distribution up to a cutoff k 0 = 2π/λ, modeling light-scattering from a nearresonant laser beam with wavelength λ. As indicated above, we create ensembles of nominally identical trajectories (1024 trajectories in 1D and 256 in 2D), with each trajectory having an independent noise realizations.
We find that the momentum distributions in 1D did not approach a thermal equilibrium state; this is the expected behavior because the 1D GPE is an integrable model. To facilitate equilbration, we break integrability with a disorder potential. After adding the desired energy, we adiabatically ramp on this potential over 2 s, allow the system to equilibrate for 5 s before adiabatically removing it over another 2 s. This process left the total energy unchanged but increased the entropy, giving the 1D momentum distributions in Fig. 5(a-i), for three values of ∆E chosen to give f c ≈ 1, 0.3, and 0. Our thermalization procedure yielded smooth distribution in the higher energy-middle and bottom-panels. By contrast, we can see that the lowest energy case (top) has a sharp peak at the k = 0 mode. Figure 5(a-ii) shows that in 1D the un-condensed fraction 1 − f c increases linearly for small ∆E and crosses over to the constant value of 1 (i.e., f c = 0) at a "critical energy" E c . The un-condensed fraction increases by no more than 0.1 at energies greater than E c . Finally, Fig. 5(a-iii) plots f c as a function of the ground state occupation probability P k=0 . We find that these quantities are nearly equal, consistent with the intuitive picture that condensate is associated with a macroscopic occupation of the ground (k = 0) state.
In Figure 5(b-i), we present the 2D momentum distributions projected onto the k x axis for the same approximate f c values as shown in 1D. As expected, we observe a sharply peaked distribution centered at k = 0 for low energy (top), and broad smooth distributions at higher energies (middle, bottom). No integrability-breaking disorder potential was required for thermalization in 2D.
Unlike the simple linear scaling of the un-condensed fraction with ∆E observed in 1D, Fig. 5(b-ii) shows that in 2D there exist regimes of different slope. In addition,    5(b-iii) shows an approximately quadratic dependence of f c on the ground state occupation P (k = 0). Although these behaviors are not a focus of our study, it is likely that as a function of increasing ∆E they result from the role of vortex anti-vortex pairs. For example, states with a small number of vortices will have a very large P (k = 0), but depending on the vortex configuration can have very little wavefunction overlap leading to a small value of f c .

B. Cooling from excited states
Here we add energy as described in Sec. II A, resulting in incoherent states above E c , and then apply the optimized cooling protocol to generate far from equilibrium condensed states. Figure 6(a) shows ∆E as a function of time with cooling starting at t = 0 s (green) starting with an added energy ∆E/E c = 2.05. For this amount of added thermal energy, our protocol quenches through the critical energy E c (red).
We observe that the cooling protocol is quite ineffective in energy removal for large added energy. To clarify this observation, we plot the energy divided into two parts: the kinetic energy contributions from momentum states above (dotted) and below (dashed) k 0 (the contribution from interaction energy is negligible). When the cooling process begins, the kinetic energy from states with k < k 0 immediately plummets to a fraction of the total energy as those modes are efficiently cooled. In stark contrast, the energy from states with k > k 0 slowly decreases over many seconds. Because these modes are above k > k 0 they cannot be directly cooled, and instead are slowly depopulated via scattering processes.
This effect is clearly visible in the inset to Fig. 6(c) that shows the final state momentum distribution with prominent shoulders at ±k 0 . This distribution further displays sharp k = 0 condensate peak surrounded by modes with occupation probabilities around 10 −6 . Relatively highly populated modes appear starting at k = k 0 . We attribute the decoupling of these groups of modes and the resulting far from equilibrium state to the integrability of the 1D GPE; essentially stating that in 1D binary collisions do not exchange momentum.
The non-equilibrium momentum distribution in Fig. 6(c) is not observable with our weak measurements due to the finite imaging resolution set by k 0 . However, time-of-flight imaging gives direct access to the momentum distribution, thus it could be observed experimentally.
We find that our quench-cooling protocol can also result in states with integer valued winding numbers associated with quantized supercurrents, with occurrence probability plotted in Fig. 7(a) fit to a Gaussian (with RMS width = 0.98). Figure 7(b) plots the the wavefunction phase as a function of position associated with winding numbers 0, 1 and 2. These states with non-zero winding number are still condensed, each with a single highly populated k-mode given by the the winding number; for example, a state with winding number 2 will be condensed into the k = 2 × 2π/L mode.
When [ Fig. 7(c), blue], indicative of a poor condensate. However, this is an artifact of the winding number distribution in Fig. 7(a), which leads to a formally fragmented condensate [26]: one with many macroscopically occupied modes. We can construct sub-ensembles for each winding number, and Fig. 7(c) shows that f c for these modes separately have f c > 0.9.

V. OUTLOOK
We presented a continuous weak measurement and feedback protocol that effectively removes high-energy excitations, up to the limit given by light-scattering, and preserves low-entropy BECs in both 1D and 2D. These results show a route toward an essential feature of any continuous feedback protocol, which is that the same quantum system can be subjected to repeated measurement without requiring an additional reservoir to remove the entropy added by the measurement process.
We furthermore demonstrated that these feedback cooling protocols can lead to new far from equilibrium condensed states that would be inaccessible in equilibrium systems. The observed separation of scales is reminiscent of condensation in pumped systems such as excition polariton condensates that can combine a nonthermal distribution of highly excited modes with a coherent condensate mode [27]. Interpreted in the context of traditional feedback theory, the observed momentum distributions are similar to feedback stabilized systems with Bode peaks located near the bandwidth limit.
The distribution of winding numbers in Fig. 7(a) for these states is reminiscent of the predictions of the Kibble-Zurek mechanism [28,29] describing quenches through phase transitions. In this work we did not study the width of these distributions as a function of quench speed or system size, as would be required to test if this feedback-cooling transition obeys Kibble-Zurek scaling. If it does, a potentially interesting direction for future study would be searching for non-mean field exponents using systems with novel or non-local feedback [18]. We note that the winding number is observable, and has for example been studied in superfluid ring experiments [30].
Additional work is needed to characterize analogous far from equilibrium condensed states in 2D systems. In 2D, quantized vortices-topological defects-are analogous to the winding number in 1D, but the large parameter space given by their number and positions within the condensate make the construction of low-entropy subensembles impossible and underlies the non-condensed Berezinskii-Thouless-Kosterlitz superfluid. rather than the Wiener increment. The Heaviside step function used in the manuscript is a special case of a dimensionless momentum-space apodizing functionf (k) withf (0) = 1 andf (k) =f (−k), implying the spatial correlation function m(r)m(r ) = f (r − r )/2.
Discrete temporal domain In numerical models we consider time points spaced by δt and generally define quantities on this grid as the average over the time interval δt, where δt is made small enough that dynamics during this interval can be neglected. Thus we have averages such as By contrast, random variables such as m t (r) require a slightly different average to give the temporal correlation m t m t = δ t,t /2. This results from the usual √ N scaling of N random events Together these lead to where κ = ϕ √ δt. Discrete spatial domain of length L For an accurate simulation, the process of discretizing the spatial domain into a grid with spacing δx requires that f (r) to be nearly constant over a range of δx; in our specific case this implies that 2π/δx k 0 . Applying our Fourier transform relations (see Appendix B for our conventions) again gives the discrete correlation function Continuous finite spatial domain of length L In this case the problem is a straightforward application of the suitable inverse Fourier transform under the assumption that the extent of the domain L is much larger than the characteristic width of f (r). This promptly gives m(k)m * (k ) = L Df (k)δ k,k /2, the spatial part correlation function following Eq. (5) in the main text.
Following such a finite-duration measurement the condensate density is n t (r) ≈ n t (r) [1 + 2κm t (r)] (A6) at lowest order in κ. Equating this expression with Eq. (A4) and assuming a constant density n 0 allows us to define the optimal measurement strength κ * = 1 2n 0 (A7) for which the final density and measurement outcomes are equal [20]. This concept can be extended to estimators of the form where K(t) is an arbitrary normalized filter kernel. For the special case of a box kernel with width τ the estimator is where m est,j is the filtered noise signal: a random variable again constructed to have variance of 1/2. The corresponding ideal measurement strength is The ideal measurement strength is reduced by a factor of δt/τ reflecting the fact that multiple measurements of length δt can be taken in the time τ . This yields the expected 1/ √ N scaling with the number of measurements taken in the time interval τ .
(DI) Discrete infinite domain with spacing δx. Similar to the CF cases, the pair f (r) = dk 2πf (k)e ikx ↔f (k) = δx x f (x)e −ikx now has a sum over x associated with the x integral of the CI case.
(DF) Discrete finite domain consisting of N sites. In this case we approximate the differential in both sums to arrive at Notice that this is not the standard textbook convention for which both expressions would have prefactors of 1/N 1/2 rather than 1/L and δx. Our motivation for making these choices is that the units of f (r) andf (k) are unchanged and the connection to continuum systems is direct and clear. Where there are no complex terms, so we have written Ψ * |m = Ψ |m . This formal expression can be further simplified on physical grounds. First, we ignore the nonolocal terms ∝ ϕ 2 2dt 2 m(r)m(r ) − dtδ(r − r ) . These terms vanish for perfectly uncorrelated white noise (in the spatial domain) and are small as k 0 increases, which sets the scale of spatial noise correlations. We can also replace m 2 (r) by it's correlation function m 2 (r); a reasonable approximation for fluctuating quantities.
The remaining local ϕ 2 term can be further simplified as (C8) where we use the continuum limit for the noise correlation function. In going from the second to the third line, we can make the approximation [Θ(k 0 − |k|) − 2] ≈ −Θ(k 0 −|k|), which is valid for large k 0 . The integral can then be carried out exactly. The final result is valid for dimensions D = 1, 2 but could be generalized to higher dimensions. Finally, we are left with the updated state Ψ |m (r) = 1 + ϕm(r) − ϕ 2 π k 0 4 D Ψ(r)dt.