{'fraction_non_alphanumeric': 0.0328273044274594, 'fraction_numerical': 0.009744872201284444, 'mean_word_length': 4.415093140059763, 'pattern_counts': {'":': 0, '<': 15, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Bloch Oscillations of Cold Atoms in a Cavity: Effects of Quantum Noise

We extend our theory of Bloch oscillations of cold atoms inside an optical cavity [ Venkatesh et al. Phys. Rev. A 80 063834 (2009)] to include the effects of quantum noise arising from coupling to external modes. The noise acts as a form of quantum measurement backaction by perturbing the coupled dynamics of the atoms and the light. We take it into account by solving the Heisenberg-Langevin equations for linearized fluctuations about the atomic and optical mean fields and examine how this influences the signal-to-noise ratio of a measurement of external forces using this system. In particular, we investigate the effects of changing the number of atoms, the intracavity lattice depth, and the atom-light coupling strength, and show how resonances between the Bloch oscillation dynamics and the quasiparticle spectrum have a strong influence on the signal-to-noise ratio, as well as heating effects. One of the hurdles we overcome in this paper is the proper treatment of fluctuations about time-dependent mean fields in the context of cold-atom cavity QED.


I. INTRODUCTION
When quantum particles in a periodic potential of period d are subject to a weak additional constant force F they do not uniformly accelerate like free particles, but instead undergo Bloch oscillations [1] at an angular frequency given by: Bloch oscillations (BOs) of cold atoms in optical lattices were first observed in 1996 by uniformly accelerating the lattice [2]: in a frame co-moving with the lattice the atoms experience a constant force. At about the same time, the accelerating lattice method was used to observe Wannier-Stark ladders [3], which are a different aspect of the same "tilted lattice" physics. The method has subsequently been employed to realize beam splitters for atom optics capable of large momentum transfers, see, e.g. [4].
In gravity-driven BOs the lattice is held fixed in space but oriented vertically so that gravity provides the force F g = mg on the atoms (of mass m). From Eq. (1), a measurement of ω B corresponds to a measurement of the applied force F if we know d/ . This Bloch oscillator may be viewed as an interferometer in momentum space [5] and has been experimentally demonstrated by a number of groups [6][7][8][9]. For example, the experiment [9] used gravity-driven BOs of strontium atoms to measure the local acceleration due to gravity at the level of ∆g/g = 5 × 10 −6 . Like any interferometer, long coherence times are crucial for precision measurements and in [9] the BOs were coherent over 7 s, corresponding to ≈ 4000 oscillations. This remarkable degree of coherence was greatly facilitated by the choice of strontium atoms, which have very weak s-wave scattering, and thus dynamical instabilities normally associated with superflow in lattices [10,11] were highly suppressed. Variations on this scheme that improve the visibility of the BOs, including frequency [12] and amplitude [13,14] modulation of the lattice, have allowed for the measurement of gravity at the level of ∆g/g = 10 −9 . In these latest experiments the BOs were coherent for over 20 seconds.
The experiments referred to above all involve destructive measurements of the BOs due to the nature of the imaging process of the atoms, whether it be in situ or by a time-of-flight technique after the lattice has been switched off [14]. Therefore, a precision measurement of ω B by the above methods requires that the experiment be re-run many times, each run being for a slightly different hold time, so that the oscillations can be accurately mapped out. This not only takes a long time, but also requires that the initial conditions be recreated as faithfully as possible for each run.
In [15] we proposed a scheme for continuous (i.e. nondestructive) measurements of BOs based upon placing the atoms inside a Fabry-Perot optical resonator which would allow for an estimate of ω B from the data acquired over a single run. A related scheme has also been independently proposed for ring cavities [16]. The periodic potential is now provided by the standing wave of light which forms inside the cavity when it is pumped by a laser. Orienting the cavity vertically, the atoms execute BOs along the cavity axis as depicted in Fig. 1. The enhanced atom-light coupling inside a high-Q cavity means that the oscillating atoms imprint a detectable periodic modulation on both the phase and amplitude of the light which can be seen either in transmission or reflection. Thus, the measurement is performed upon the light leaking out of the cavity rather than directly upon the atoms.
The strong atom-light coupling that can be realized in cavity-QED stands in contrast to the case of optical lattices in free space where the atoms exert only a tiny backaction upon the light. The optical dipole interaction between a single cavity photon and a single atom is characterized by the Rabi frequency g 0 = (µ/ ) ω c /(ε 0 V ), where ω c and V are the frequency and volume of the relevant cavity mode and µ is the atomic transition dipole moment. Defining the cooperativity C ≡ g 2 0 /(2γκ), where 2γ is the spontaneous emission rate of the atom arXiv:1212.3594v3 [quant-ph] 30 Jul 2013 in free space and 2κ is the energy damping rate of the cavity, 1/C is the number of atoms required to strongly perturb the light field. The normal mode splitting that results from strong coupling has been directly observed in a number of cold atom optical cavity experiments [17][18][19][20]. In the experiment [20], which was performed with a Bose-Einstein condensate, the cooperativity was C = 145. Even more pertinently, these systems have been used to detect the presence of single atoms [21][22][23], as well as to follow their dynamics in real time [24,25]. The collective dynamics of ultracold atomic gases have also been tracked using cavities [26][27][28]. The key experimental steps necessary for the continuous monitoring of BOs in a cavity have, therefore, already been demonstrated.
The drawback with any continuous measurement scheme is measurement backaction. In cavities this backaction typically takes the form of cavity photon number fluctuations which lead to random force fluctuations on the atoms, as is evident in the erratic nature of the single atom trajectories seen in the experiments [24,25] referred to above. In the many atom context, quantum measurement backaction generally manifests itself in a heating of the atom cloud (although under some circumstances it can lead to cooling [29]). In the cavity-optomechanical regime (where the collective motion can be modelled as a harmonic oscillator of angular frequency ω) the heating rate is expected to be R = (x zpf / ) 2 S F F (−ω) [30], where x zpf is the zero-point fluctuation and S F F is the spectral density of the force fluctuations (which is directly proportional to the cavity photon number fluctuations). This heating rate is in agreement with observations when convolved with technical fluctuations [31].
In the system considered in this paper (see Fig. 1), we can divide the backaction into two types. One type comes from the fact that the atoms sit in an optical lattice whose depth is periodically modulated in time at the frequency ω B due to the effect of BOs. This backaction is a classical effect in the sense that it occurs even when the light field is treated classically (no photons). The nonlinearity that arises from this backaction can lead to swallowtail loops in the atomic band structure [32,33] that mimics the effects of direct atom-atom interactions [34][35][36][37]. These loops are the counterpart in the atomic wave function of optical bistability in the light [26][27][28]38]. The second type of backaction arises only when the fluctuations due to the discrete photon nature of the light field are taken into account and is related to the heating effect mentioned above. The characteristic frequency of these latter fluctuations is κ which is much larger than ω B .
The first type of backaction was analyzed in our previous paper [15] where our main aim was to show that, despite the self-generated time modulation of the intracavity optical lattice, the Bloch acceleration theorem still applies and the BO frequency is not modified (although harmonics can be generated). This latter result is clearly very important if the cavity BO method is to be used for precision measurements and may be viewed as a conse-  [15]. A dilute cloud of cold atoms undergoes BOs in the combined intracavity lattice potential and the acceleration due to gravity. The transmitted light field's intensity and phase are modulated at the Bloch frequency. A in-situ precise measurement of the Bloch frequency (and hence the force) can be performed by detecting the transmitted light.
quence of the fact that the formula (1) does not depend on the depth of the lattice, only its spatial period. An estimate of the effects of the second type of backaction was also given in our previous paper, but this estimate was obtained under the assumption that the photon number fluctuations were purely due to the photon shot noise found in a coherent state of light. This ignores the correlations that build up between the atoms and the light inside the cavity and our main aim in this paper is to solve the dynamics of the coupled photon and atom fluctuations systematically from first principles and thereby capture these correlations. This will allow us to properly determine the sensitivity of the measurement of the Bloch frequency to quantum fluctuations.
The plan of this paper is as follows: in Section II we introduce the physical system, the associated hamiltonian, and the equations of motion. We then review in Sections III and IV the meanfield approximation and the associated numerical results which were the focus of our previous paper [15], before introducing in Section V the main model to be treated in this paper which adds quantum fluctuations. This is an elaboration of the linearization approach presented in, e.g. [39][40][41], to include a time-dependent meanfield component (due to the BOs). The fluctuations correspond to quasiparticles (excitations out of the meanfield), and their spectrum is analyzed in Section VI and then used to help interpret the numerical results for the quantum dynamics presented in Section VII. We also develop a simple rate equation picture, valid in the weak coupling regime, to help us understand the rate of quasiparticle excitation. Following this we change gears slightly and apply the above results to investigate how quantum fluctuations affect a precision measurement of ω B by calculating the signal-to-noise ratio (SNR). We present the theory lying behind these calculations in Section VIII and in Section IX we examine the results, paying particular attention to whether or not there is an optimal value for the atom-light coupling parameter β = N U 0 /κ. We also present results illustrating the dependence of the SNR on other system parameters such as the number of atoms and the intracavity lattice depth. We summarize our results and give some further perspective in Section X. We have also provided three appendices that give details omitted from the main text: the first derives an approximation wherein the cavity field is assumed to be in a coherent state and the atomic fluctuations about the meanfield are treated as independent oscillators, the second discusses the effects that BOs have on cavity cooling, and the third discusses our approach to calculating two-time correlation functions.

II. HAMILTONIAN AND EQUATIONS OF MOTION
Our system consists of a gas of N bosonic atoms inside a vertically oriented Fabry-Perot optical cavity. A single cavity mode of frequency ω c is coherently pumped by a laser with frequency ω p that is detuned from both the atomic and the cavity resonance frequencies. This sets up a standing wave mode along the cavity axis of the form cos(k c z), where k c = ω c /c. The relevant frequency relations are characterized by the two detunings where ω a is the atomic transition frequency. In the dispersive regime, the occupation of the excited atomic state is vanishingly small and it can be adiabatically eliminated. A one-dimensional hamiltonian for the atomcavity system in the dispersive regime can then be written as [42,43] whereΨ(z, t) andâ(t) are the field operators for the atoms and the cavity photons which obey the equal time bosonic commutation relations [Ψ(x, t),Ψ † (x , t)] = δ(x− x ), and [â(t),â † (t)] = 1, respectively. The single atom dispersive light shift has been denoted by U 0 ≡ g 2 0 /∆ a .
The hamiltonian has been written in a frame rotating with the pump laser frequency ω p , and this leads to the appearance of the two detunings. The first term is just the free evolution of the cavity mode. The second term represents the laser coherently pumping the cavity at rate η, and the third term describes the atomic part of the hamiltonian. The first two terms of the atomic part represent the kinetic energy and a light induced potential energy. This latter term can either be understood as the atom moving in a periodic potential with average amplitude U 0 â †â or, if combined with the first term in the hamiltonian, as a shift in the resonance frequency of the cavity due to the coupling between the atom and the field. The third term in the atomic part provides the external force that drives the BOs. We assume this force arises from the vertical orientation (z increases in the downward direction) of the cavity and is given by We have not included direct atom-atom interactions in the hamiltonian (3) because under realistic experimental conditions they are three orders of magnitude smaller than the recoil energy E R ≡ 2 k 2 c /(2M ) which characterizes the single-particle energy (kinetic and potential) of an atom in an optical lattice. Consider, for example, the meanfield interaction energy per particle E int /N = (2π 2 a s N/M ) |Φ(r)| 4 d 3 r for a cloud of N = 5 × 10 4 87 Rb atoms trapped in a 178 µm long cavity [28]. Here a s = 5.3 nm is the s-wave scattering length. We take the normalized 3D wave function Φ(r) to be the product of a ground band Bloch wave that extends 178 µm along z and a gaussian 25 µm wide in the transverse plane. Then, evaluating the Bloch wave for a lattice which is 3 E R deep and made from 780 nm light (456 wells are occupied), we find the ratio E int /N : The interactions can be tuned to smaller values still using a Feshbach resonance: the experiment [44] increased the dephasing time of BOs from a few oscillations to 20,000 using this technique. The fact that the atoms all interact with a common light field whose magnitude is modified by the sum of their individual couplings gives rise to a nonlinearity (the classical backaction referred to above) that is in some ways analogous to that due to direct interactions [32,45], but in other ways differs and can lead to novel behavior [42,43,46]. Natural units for the length and energy in cavity-QED are given by 1/k c and the recoil energy E R , respectively. From here on we scale all lengths by 1/k c and consequently define x ≡ k c z. We scale frequencies by the recoil frequency ω R ≡ E R / and time by 1/ω R and retain the same symbols for the scaled variables. The Heisenberg-Langevin equations of motion for the light and atomic field operators in the scaled variables are [42] is the dimensionless form for the external force. The operatorξ(t) is the Langevin term and is assumed to be Gaussian white noise with the only non-zero correlation being Mathematically, the Langevin noise terms are necessary in order to preserve the commutation relation [â(t),â † (t)] = 1 in an open system. Physically, their origin is vacuum fluctuations of the electromagnetic field that are transmitted into the cavity via the mirrors and they thus only appear in the equations for the light field. Nevertheless, the noise is conveyed to the atomic dynamics by the atom-light coupling.

III. MEANFIELD DYNAMICS: THEORY
The approach we follow in this paper is based upon a separation of the field operators into meanfield and quantum parts:â In the meanfield approximation the light is assumed to be in a classical state with amplitude α(t) = â(t) , where |α(t)| 2 corresponds to the average number of photons in the cavity, and the atoms are assumed to all share the same single-particle wave function ϕ(x, t) = Ψ (x, t) / √ N . The equations of motion for the meanfield amplitudes α(t) and ϕ(x, t) are where the second equation has the form of a Schrödinger equation. The expectation value that appears in the first of these equations provides the time-dependent coupling between the atomic probability density and the cavity mode function. Multiplying this integral is the collective atom-cavity coupling parameter N U 0 . When measured in units of the cavity linewidth we denote this parameter by β We illustrate the effect that β has on the meanfield dynamics in Figs. 2 and 3 below.
In reference [32] we studied the influence the classical backaction nonlinearity has upon the band structure of atom-cavity systems. The band structure is given by the steady state solutions [α = 0, ϕ(x, t) = ϕ(x) exp(−iµt/ )] of the coupled equations of motion (6a) and (6b) in the absence of the external force f . It is straightforward to see that, despite the nonlinearity, exact solutions of the steady state problem are given by Mathieu functions (like in the linear problem of a quantum particle in a fixed cosine potential). Mathieu functions are Bloch waves and so can be labelled by a band index b and quasimomentum q [47] where U q,b (x + π) = U q,b (x) has the same period as the lattice. In the reduced zone picture q is restricted to lie in the first Brillouin zone −1 < q ≤ 1. Substituting the Bloch wave solution into the equations of motion yields the steady state equations where the subscript ss denotes "steady state". Solving these equations one obtains a band structure analogous to that in the linear case but with the striking difference that the nonlinearity can lead to swallowtail loops in the bands. It is important to appreciate that this band structure is not for the atoms alone, but for the combined atom-cavity system. For example, the eigenvalue µ is actually a chemical potential rather than the band energy (for the underlying energy functional with the light adiabatically eliminated see [32]), and another difference from the linear case is that the lattice depth s = U 0 |α ss | 2 is not fixed, but instead depends on the values of {b, q}. So, for example, the lattice depth changes during a BO as q is swept along the band.
The external force f breaks the spatial periodicity and means that Bloch waves are replaced by Wannier-Stark states as the stationary solutions of the equations of motion (in fact, in finite systems the Wannier-Stark states are resonances rather than true eigenstates [48]). The spatial periodicity can be restored by applying the unitary transformationφ(x, t) = exp(−if tx)ϕ(x, t) which removes the f x term appearing in the hamiltonian in the Schrödinger equation (6b) and introduces a shift f t into the momentum operator We denote the frame resulting from this transformation as the transformed frame (TF), and the original frame as the lab frame (LF). Let us now consider the dynamics under the influence of the force term. We take the initial atomic state ϕ(x, t = 0) = ϕ(x, t = 0) to be a Bloch state in the ground band with quasimomentum q = q 0 . In the adiabatic approximation the atoms remain in the ground band but the force causes the quasimomentum to sweep periodically through the first Brillouin zone in accordance with the Bloch acceleration theorem as can be seen by comparing Eqns. (9b) and (10). In fact, a careful analysis [49] shows that Eq. (11) holds even when adiabaticity is broken and interband transitions are allowed providing these transitions are "vertical", i.e. they conserve q. This standard approach to BOs remains valid even when the lattice depth is modulated in time, as takes place in cavities, because amplitude modulation does not break the spatial periodicity of the potential and so cannot change q [15]. We therefore find that at any later time t, the exact atomic meanfield can be expressed as In general U(t) is in a superposition of bands and so is no longer the steady state solution of Eqns. (9a) and (9b), although it does retain its Bloch form. The advantage of the TF is that the quasimomentum is frozen at its initial value and we havē so that it is only the spatially periodic function U(t) that evolves in time. From the point of view of numerical computation this allows us to work with a basis of periodic functions (we normalize our wave functions over one period of the lattice). At any given time a relatively small number of basis functions can accurately describe the atomic meanfield state and this greatly reduces the numerical effort in the calculation of BOs. By working in terms of Bloch waves, our approach is predisposed towards treating wave functions which are localized in momentum space rather than coordinate space. This choice is sensible because momentum space is a natural setting for BOs as is evident from Eq. (11). This is also in line with existing experiments demonstrating cold atom BOs in free space optical lattices [2, 6-9, 12, 13], where the initial state is generally a fairly narrow wavepacket in momentum space. In this paper we shall therefore restrict ourselves to states that are completely localised in quasimomentum (δ-function wave packet).

IV. MEANFIELD DYNAMICS: RESULTS
We now present our numerical results for the meanfield dynamics. The initial state at time t = 0 is taken to have quasimomentum q = 0, and be given by the solutions α ss and U 0,0 (x) of the meanfield steady state equations   Under the action of the external force the atoms begin performing BOs, which for atoms in extended Bloch states gives rise to a breathing motion of the atomic density distribution on each lattice site [15]. The classical backaction imprints an oscillation on the amplitude and phase of the light field at the Bloch frequency ω B . In Fig. 2 we plot the time-dependence of the intracavity lattice depth s(t) = U 0 |α(t)| 2 seen by atoms, which is proportional to the number of cavity photons |α| 2 . The experimental signature of the BOs is the photon current transmitted by the cavity, and this is given in the meanfield approximation by κ|α(t)| 2 , and hence is directly proportional to s(t).
The size of the backaction is controlled by the collective coupling β = N U 0 /κ, as is apparent from the different curves in Fig. 2. As β is increased the change in the lattice depth over a Bloch period increases and hence the visibility or contrast of the BOs as measured by a photon detector outside the cavity increases also. We define the contrast as Each curve in Fig. 2 has a different pumping strength η in order to maintain the same minimum lattice depth of 3E R . If the lattice becomes too shallow interband transition rates (e.g. due to Landau-Zener tunnelling around the band edges) become so high that the atoms effectively fall out of the lattice. On the other hand, if the lattice becomes too deep the contrast decreases (see Fig.  9a below and also Fig. 5 in [32]). A depth of 3E R gives a reasonable compromise. Therefore, although in the rest of this paper we will examine the effects of changing the various system parameters, we will always maintain the minimum lattice depth at 3E R (except in Fig. 9a and Fig.  14). This also allows us to make comparisons between the effects of different parameter values upon, e.g. the quantum fluctuations, whilst keeping the atomic meanfield dynamics as similar as possible.
As the coupling β is increased other effects appear apart from an increase in the contrast. These effects are visible in Fig. 3 (see also Fig. 2 in [15]). In Fig. 3a we see that small-amplitude fast oscillations of the lattice depth appear on top of the basic BO. Referring to the Fourier transform of s(t) plotted in Fig. 3b, we see that the basic BO dynamics is governed by the fundamental ω B and its low lying harmonics, whereas the fast oscillations are clustered around the tenth harmonic (see inset) and include a continuum of frequencies with some peaks at half harmonics. In this context it is important to bear in mind that the band gaps change continuously in time as q is swept through the Brillouin zone and so a range of frequencies is to be expected.
The precision to which ω B can be measured in the scheme proposed in this paper depends upon the contrast. From the results shown in Fig. 2 it may therefore seem that in order to make the most sensitive measurement possible one should choose β to be as large as possible. However, this is false for two reasons. One is the effect of quantum fluctuations due to measurement backaction which is also controlled by β and will be the focus of Section VIII. Another reason, which enters even at the meanfield level, is the possibility of bistability in cavity photon number for large values of β (when the pumping is sufficiently large). In [28] this bistability was studied experimentally in a uniform unaccelerated condensate, which in our language has a quasimomentum q = 0. In [32] we studied this problem theoretically and generalized it to include finite q: we showed that bistability arises from the appearance of swallowtail loops in the bands. In the semiclassical picture of a BO the quasimomentum scans adiabatically through the entire band and so when it encounters a swallowtail loop the system can follow a branch that suddenly terminates at some later time, leading to fundamentally nonadiabatic behavior [15,50]. Hence, in a scheme to measure BOs, it would be better to be in a parameter regime where the cavity is not bistable for any value of q. In Fig. 4 we plot the pump strength required to maintain the lattice depth at a mini- mum value of 3E R as a function of β. The red (solid) and blue (dot-dashed) lines enclose the values of η for which the steady state photon number in the cavity displays bistability for at least some values of the quasimomentum. We see that for β values as large as 25 (at the fixed detuning ∆ c = −0.75κ) one can avoid bistability and get large contrast in the lattice depth evolution. Having emphasized that our choice for the pumping strength η is guided by the tradeoff between contrast and bistability according to Fig. 4, let us now explain how we chose the rest of the system parameters used in the calculations. There are three parameters we hold constant throughout this paper; the first is the cavity damping rate κ = 345 ω R which is the value realized in the experiment [28]. As a guide to the magnitude of the atomic recoil frequency ω R used as the frequency unit, we note that for 87 Rb atoms in 780nm light ω R = 2π × 3.8 kHz. The second constant parameter is the Bloch frequency ω B = ω R /4. The gravitational force on 87 Rb atoms in a 780nm lattice provides a Bloch frequency very close to this value. Finally, unless specified otherwise, we keep the meanfield atom number fixed at N = 5 × 10 4 and vary U 0 in order to vary β. This last choice is motivated by a scaling symmetry of the meanfield equations [Eqns. (6a)-(6b)], which also holds for the quantum operator equations [Eqns. (22a)-(22b)] below; solving the coupled equations for the set of parameters {U 0 , N, η, ∆ c , κ} is exactly the same as solving them for {U 0 r, N/r, η/ √ r, ∆ c , κ}, where r is some positive scaling factor. In both the scaled and unscaled versions the lattice depth s(t) = U 0 |α(t)| 2 is maintained at the same value. Thus, each specific calculation one performs represents a family of parameters. Our choice for N keeps the atomic density dilute enough in a typically sized cavity that the approximation of ignoring collisional atom-atom interactions remains valid. There is some latitude in the choice of ∆ c , but the contrast one obtains at a given value of β is larger for ∆ c closer to the cavity resonance.
On the other hand, one also has to make sure that the effective cavity detuning is less than zero so that we are in the cavity cooling regime for the fluctuations [40] (see Section VI). We set ∆ c = −0.75κ since we find that it maximizes the contrast for the coupling value of N U 0 /κ = 1. We will examine the effect of changing the number of atoms N and the minimum lattice depth when we examine the signal-tonoise ratio in Section VIII.

V. QUANTUM DYNAMICS: THEORY
The approach we take to quantum dynamics is based upon a linearization about the meanfield solution, retaining the quantum operators δâ and δΨ only to first order in the equations of motion. This corresponds to the Bogoliubov level of approximation [51,52], suitably generalized to describe coupled atomic and light fields. A new feature of our problem in comparison to previous linearization-based treatments of cavity-QED systems, e.g. [39,41,53], is that our meanfield is time-dependent because of the BOs. This means that the fluctuation modes, which must be orthogonal to the meanfield mode, also evolve in time (not just their occupations).
Linearizing about the meanfield solution may appear to be an innocent strategy, but, as is well known from the theory of Bose-Einstein condensation, care must be taken with such U(1) symmetry breaking approaches because they introduce a macroscopic (meanfield) wave function with a particular global phase at the cost of particle number conservation [54]. In particular, when performing a linearization about the condensate there is always a trivial fluctuation mode parallel to it with zero frequency (the "zero mode") which corresponds to unphysical fluctuations of the global phase. These issues are even more acute when the condensate is time-dependent and the boundary between condensate and fluctuation is further blurred [55].
The zero mode problem can be handled by only including fluctuations that are at all times orthogonal to the meanfield. We achieve this by applying the projector P (t) [55,56] so that One consequence of this is that the commutator between atomic fluctuations is given by [57] δΨ Unlike the usual bosonic commutator for the fluctuation fieldδΨ, this is time dependent. Next, we transform the atomic fluctuation operator from the LF to the TF which simplifies the calculation for the same reasons as mentioned in Section III for the meanfield. Since only "vertical" fluctuations between bands can occur, both the fluctuations and the meanfield have the same quasimomentum (which in the TF is frozen at its initial value), and so we can expand the fluctuations and meanfield in the same basisφ This makes the numerics a little easier. Note that we have set the initial quasimomentum in these equations to q 0 = 0 without loss of generality. Meanwhile, back in the LF, the quasimomentum evolves according to the Bloch acceleration theorem given by Eq. (11).
We can now write down the coupled equations of motion for the cavity and atomic fluctuation operators in the TF as The structure of these equations is such that without the Langevin termξ(t) the operators δâ andδΨ ⊥ would be fixed at their initial values and so the quantum parts of the fields would remain zero for all time. The Langevin fluctuations appear as an inhomogeneous term in the cavity field equation and act as a source that drives the evolution of δâ which in turn drives the evolution ofδΨ ⊥ via the atom-cavity coupling.
As pointed out in [41], the dynamics of the complex valued operators in the above equations can be solved either by separating out their real and imaginary parts (optomechanics approach) or by simultaneously solving the equations for the hermitian conjugates of the operators (the Bogoliubov-de Gennes approach). We choose the latter. Collecting the fluctuations into the column vectorR(t) = δâ δâ †δΨ ⊥δΨ † ⊥ T , and the noise operators that act as source terms into the column vector where T denotes transposition, we obtain the operator matrix equation where we have introduced the operators i.e. V is an integral operator that acts on a function g(x).
Since they fall on the off-diagonals, the terms involving V and W couple the cavity and atom fluctuations. Observe, however, that in the linear approximation used here the atomic fluctuation operatorsδΨ ⊥ (x, t) are not directly coupled to the cavity fluctuation operators δâ(t) because this would lead to terms which are of second order. Rather, the coupling between the two sets of quantum fields is mediated by the meanfields α(t) andφ(x, t).
The matrix M(t) is non-normal, i.e. it does not commute with its Hermitian adjoint and its left and right eigenvectors are not the same. However, it does have the following symmetry property: a linear transformation T that swaps the first and second, and simultaneously, the third and fourth rows, produces a matrix which is proportional to the complex conjugate of the original [40] T .M.
This symmetry, which is a general feature of Bogoliubovde Gennes type equations [11], implies that the eigenvalues (and the associated eigenvectors) occur in pairs of the form ±ω n + iγ n i.e. with the same imaginary parts but with real parts of opposite sign. We shall explore the spectrum of the fluctuation matrix M further in the Section VI. We also note that when written in matrix form the role of the projection operator becomes clear since one can immediately see that the vectors (0 0φ(x, t) 0) T and (0 0 0φ * (x, t)) T span the zero eigenvalue subspace of the matrix M and the trivial fluctuations live in this subspace.
The time evolution of the fluctuation operators is given by solving Eq. (23a). However, measurable observables are given by expectation values and correlation functions of these operators rather than by the operators themselves. To this end we consider the covariance matrix C(t) associated with the vectorR Particular cases of C jk (t) (or more precisely, its sum) include the total number of photonic and atomic fluctuations The latter correspond to the number of atoms excited out of the meanfield component (i.e. the atomic depletion). To obtain the time evolution of the covariance matrix, consider the formal solution to Eq. (23) [58] (28) where G(t) is a matrix satisfying: We drop the dependence of G on the initial time for notational convenience in what follows. Inserting this formal solution in Eq. (25) we find Using the property of the Langevin noise terms given in Eq. (4c), we can simplify Σ(t) as Our main numerical task is thus to solve the matrix differential equation given by Eq. (29). In addition, the matrix elements of M(t) have to be computed from the meanfields {α(t),φ(x, t)} obtained by solving the coupled equations Eq. (6a) and Eq. (6b). These latter equations are simply a set of ordinary differential equations that we solve using an adaptive time-step Runge-Kutta scheme.
We then solve the matrix differential equation for G(t) using the same time grid as the meanfield solution. For the matrix differential equation, and the associated solution for the covariance matrix C(t), we can again use a Runge-Kutta algorithm or exponentiate the fluctuation matrix M(t) over the (small) time step intervals [59]. As a check on the results we can use the fact that the elements of the covariance matrix C(t) have to obey the commutator relations Eq. (18) for the operators making upR(t). For example, when the atomic operator is expanded as in Eq. (20), the expectation value of the commutator relation Eq. (18) gives The left hand side of this equation gives the difference between certain entries of the covariance matrix and we can calculate its expected value (the right hand side) from the meanfield solution. The degree of agreement between the two sides provides a measure of the accuracy of the fluctuation calculation. In general, we find that the accuracy can be increased by taking smaller time steps.
In closing this section, we would like to point out that the meanfield solution already includes Landau-Zener type tunnelling that causes the coherent excitation of higher bands. By contrast, the effect of Langevin fluctuations ξ is to incoherently populate different bands. Within the linear approximation used here, the depletion of the atomic meanfield by quantum excitations is not self-consistent, i.e. the meanfield is always normalized to N atoms, whatever the number of depleted atoms δN . The linearized equations are only valid when δN << N , as expected from a Bogoliubov-type approach. The red (solid) and blue (dotted) lines correspond to hybridised atom-cavity modes and generally have non-zero imaginary parts except at certain special points such as at the band center and edges where they can become marginally stable and decouple from the cavity. The green (dash dotted) line corresponds to a cavity-like mode, i.e. the real part of its frequency is close to the effective detuning frequency ∆ eff c (q), and the imaginary part is close to −κ.

VI. SPECTRUM OF ELEMENTARY EXCITATIONS
Before presenting the results of the combined meanfield and quantum dynamics (see the next section), we shall first examine the excitation spectrum of the atomcavity system. The excitation spectrum gives insight into the dynamics, heating effects, and will also be of use in explaining resonances that affect the signal-to-noise ratio of the Bloch frequency measurement, a topic we will discuss in Section IX.
We first note that there are two distinct types of excitation, and hence spectra. The coupled atom-cavity band structure discussed in Section III refers to meanfield excitations which are labelled by a band index and a quasimomentum. They involve every atom and photon responding identically since, by the nature of the meanfield approximation, they are assumed to be described by a single wave function ϕ(x, t) and the coherent amplitude α(t), respectively. On top of these, there are also elementary excitations or quasiparticles whose energies are the complex eigenvalues ±ω n + iγ n of the matrix M(t) given in Eq. (23b). A clear description of the difference between the meanfield and the quasiparticle spectra for a BEC in a (non-cavity) optical lattice can be found in [60] and references therein. In the atom-cavity system the quasiparticles correspond to single quanta of the combined fields and are thus polaritons. The fact that they come in pairs can be interpreted as an analogue of particles and antiparticles [11].
Whereas the meanfield band structure is always real, the quasiparticle energies have an imaginary part which comes from the leaking of the cavity field out of the cavity. If γ n < 0, we have dynamical stability and |γ n | can be interpreted as the lifetime of the quasiparticle. This damping effect has potentially very important applications in cavity-assisted cooling [39,61]. If, on the other hand, γ n > 0 we have dynamical instability and heating.
In general, the elementary excitations have a band structure all of their own, i.e. the solutions of the Bogoliubov-de Gennes equations take the form of Bloch waves with a band index and quasimomentum that can differ from that of the meanfield solution about which we are linearizing. However, as discussed in Section III, here we only allow excitations that preserve the quasimomentum (vertical transitions), and thus our quasiparticles have the same quasimomentum as their parent meanfield solution. Some examples of the quasiparticle band structure are plotted in Fig. 5 (see also Fig. 12a in Appendix B).
The eigenvectors of M can be classified into three kinds: cavity-like modes, hybridised atom-cavity modes, and marginally stable modes [40]. The cavity-like modes (depicted by the green dash-dotted lines in Fig. 5) are close to being pure cavity field modes with only a small atomic component. Hence, their eigenvalues have a real part with magnitude close to the effective detuning ∆ eff c = ∆ c − N U 0 cos 2 (x) , and an imaginary part approximately equal to −κ. The hybridized modes (depicted by the red solid and blue dashed lines in Fig. 5) have some atomic and some cavity field properties, whereas the marginally stable modes are purely atomic in nature with zero cavity component. As we shall demonstrate below, the marginally stable modes occur at the points q = 0 and q = ±1, i.e. at the band center and edges, and their name derives from the fact that their imaginary part is zero.
The properties of the hybridised and marginally stable modes are determined by the sign of ∆ eff c . When ∆ eff c < 0 we find γ n < 0 and we are on the cooling side of the effective resonance. On the contrary, when ∆ eff c > 0 we find γ n > 0 and we are on the heating side. A calculation of the dynamics on the heating side is not stable since the linearization will fail after a short time due to the exponentially growing number of quasiparticles. Thus, the calculation of the spectra serves a very useful purpose: it guides our choice of ∆ c so as to ensure that we are always on the cooling side of the resonance.
In Fig. 5a the red (solid) and blue (dashed) lines give the magnitudes of the real parts of the frequencies of the two lowest quasiparticle eigenmodes as a function of quasimomentum. The magnitudes of the imaginary parts are plotted in Fig. 5b. Notice that the imaginary part of one of the modes goes to zero at the band edges (red solid line) and the other goes to zero at the band center (blue dashed line). This implies that these excitations are only marginally stable at those specific values of the quasimomentum. The vanishing of the imaginary part of the frequencies at these points can be understood as follows: a Bloch wave with q = 0 (q = 1) is even (odd) about the center of a single cell [0, π] of the cos 2 (x) potential. Since the atomic meanfield solutionφ(x) is a Bloch wave it has well defined parity at these points. The same is also true for the atomic part of quasiparticle eigenmodes δΨ ⊥ (x) [and of courseδΨ † ⊥ (x)] of the fluctuation matrix M, for these are also Bloch waves. In this case the integral dxφ * (x) cos 2 (x)δΨ ⊥ (x) will sometimes vanish identically because the integrand can contain functions with opposite parity. Examination of the fluctuation matrix M given in Eq. (23b) shows that it is exactly this integral that controls the weight of the cavity part of the quasiparticle eigenmodes, and so at q = 0, ±1 we can have undamped quasiparticles with γ = 0. For other values of quasimomentum the meanfield wave function has no particular parity and there are no marginal modes.

VII. QUANTUM DYNAMICS: RESULTS
In this section we present results from the numerical solution of the quantum equations of motion. We assume that at t = 0 the fluctuation fields corresponding to δâ andδΨ ⊥ are in their vacuum states and expand the atomic part in the basis given in Eq. (21). This allows us to construct the covariance matrix Eq. (25) C(t = 0) which we then evolve to later times using Eq. (29). In order to perform this task we need the fluctuation matrix M(t) as a function of time which in turn requires the meanfield solution {φ(x, t), α(t)} as input. We therefore solve the meanfield dynamics on the same discretized time grid in parallel with the computation of Eq. (29).
Once we have computed C(t), we can use it to calculate the physical quantities of interest, such as the number of atomic excitations δN (t), as defined in Eq. (27). This should not be confused with the number of quasiparticles, which are generally made up of both atomic and cavity field components. If it were not for the Langevin noise, the evolution would be perfectly coherent and δN would be zero. However, the presence of Langevin noise in the electromagnetic field generates atomic excitations via the atom-cavity coupling. In Fig. 6 we plot the fraction δN (t)/N as a function of time for five Bloch periods for two different coupling values. The red (solid) curves are given by a full solution of the quantum equations, whereas the black (dashed) curves are made with a coherent state approximation for the cavity field, which will be outlined below and is discussed in more detail in Appendix A. The gradient of the curves in Figs. 6 gives the heating rate and we note from Fig. 6b that the coherent state approximation slightly underestimates the true heating rate for atoms. The behaviour of δN (t)/N over longer times (40 Bloch periods) is shown in Fig. 7a and the equivalent quantity for the photons is shown in Fig.  7b. We see that while the number of photons excited out of the meanfield saturates due to the damping by photon loss from the cavity, the atoms maintain a finite heating rate over all times we have investigated. This is perhaps surprising because we are on the cooling side of the resonance [see Eq. (15)] all the time (despite the modulations in the effective cavity detuning due to the BOs). In the inset in Fig. 7a we show the case without BOs, and as can be seen, we recover the cooling. The presence of BOs clearly counteracts the cooling to some degree and prevents δN (t)/N from reaching a steady state. This residual heating effect is analyzed in detail in Appendix B, but it turns out to be due to the transport of quasiparticles to higher energy states by Landau-Zener transitions that are driven by the BOs. A finite heating rate implies that at long propagation times the validity of our linearized approach will break down because δN/N will no longer be small. However, for all the times and coupling strengths considered in this paper we have δN/N < 1/100.
In order to gain further insight into the dynamics, let us develop a semi-analytic model that we can compare against the exact results: it will allow us to see when atom-light correlations are important. The model makes two approximations: firstly we treat each eigenstate of the instantaneous meanfield Hamiltonian Eq. (10) as an independent oscillator mode uncoupled from the other modes, and secondly we approximate the state of the  Fig. 2. The atom heating rate in these figures oscillates because it is lower at the Brillouin zone edges than at the center. Referring to Fig. 5 we see that at the zone edges the quasiparticle mode with the smallest real part (red solid curve) becomes marginally stable, i.e. the cavity light field part and the atomic part decouple.
light inside the cavity as a coherent state. Coherent states have a noise spectrum that corresponds to the vacuum and so neglect correlations with the atoms. In fact, the second approximation follows naturally from the first as we show in Appendix A. The results of the approximate model are the black dashed curves in Fig. 6. The agreement with the exact results at weak coupling (β = 0.1) is excellent, but begins to break down over time at stronger coupling (β = 1), thereby revealing the dynamic generation of correlations. In fact, as mentioned in the Introduction, the heating rate of a cloud of cold atoms inside a cavity has been measured by Murch et al [31], and they found it to be consistent with the predictions of vacuum noise. The new feature in our problem is that the effective cavity drive detuning ∆ eff c (t), which appears in the phase terms in Eq. (35), is Bloch periodic due to the meanfield dynamics.
To motivate the coherent state approximation consider the exact solution to the first order inhomogeneous differential equation for cavity field fluctuations Eq. (22a), which can be formally written as whereX(t) ≡ dxφ * (x, t) cos 2 (x)δΨ ⊥ (t) + h.c., and A(t) = −∆ eff c (t) − iκ, as above. We see that the cavity field fluctuation has two distinct contributions: the first term depends on the Langevin noise which accounts for vacuum fluctuations, whilst the second term depends on the state of the atoms. The coherent state approximation consists of dropping the latter term in favour of the former to give In writing δâ(t) in this way we have taken advantage of the fact that the cavity decay rate κ is much faster than the frequency ω B at which ∆ eff c (t) evolves, and so the integrand is appreciable only for times t − t κ −1 during which ∆ eff c is a constant and can be evaluated at time t. The regime of validity of the coherent state approximation can be estimated from its derivation which requires r ≡ √ N U 0 |α(t)|/ √ 2κ << 1. Note that r 2 = βs(t)/2. In our earlier discussion (Section IV) of desirable parameters, we stipulated a minimum lattice depth of s(t) ∼ 3 ω R , which implies that the validity of the coherent state approximation here is contingent upon β 1, i.e. this is a weak coupling approximation.
The assumption of uncorrelated vacuum noise is a standard one in the field of cavity optomechanics [30,[62][63][64][65]. The paradigmatic example is a cavity with one end mirror attached to a spring or cantilever i.e. a harmonic oscillator driven by radiation pressure. Although ultracold atoms in a very shallow lattice in a cavity can be mapped onto this system [28,31,41,66,67], that is not the case here because the atomic Bloch states do not map faithfully onto a single harmonic oscillator. Nonetheless, we have obtained our approximate model by applying a similar philosophy by mapping onto a collection of independent oscillators (the eigenstates of H). The coherent state approximation for the atomic excitation occupation δN (t) that is plotted as the black (dashed) curves in Fig. 6 is the sum over the occupation numbers of these independent oscillator modes δN (t) = j δN j (t). The details of the mapping are presented in the Appendix A and here we shall only sketch out the main idea which is to consider the noise as a perturbation to the oscillator dynamics, and then use Fermi's golden rule to calculate the noise induced transition rates amongst each oscillator's states. This leads to a rate equation describing the occupation number dynamics for each oscillator [63] which is Eq. (A9) in Appendix A. In this expression Γ uj and Γ dj are the transition rates "up" and "down" for the jth oscillator and they are proportional to S F F (−ω j ) and S F F (ω j ), respectively, where S F F (ω) is the spectral density of force fluctuations (shot noise power spectrum). Thus, each oscillator is driven and damped by vacuum noise, with the rates of driving and damping being time dependent (due to the meanfield BO dynamics).
In the next two sections we examine the effects of the fluctuations upon a precision measurement, i.e. how the fluctuations put a limit on how large a value of β can be chosen for a precision measurement.

VIII. SIGNAL-TO-NOISE RATIO: THEORY
We now explore how the inclusion of quantum noise affects the precision measurement proposal in [15]. Recall the basic idea shown schematically in Fig. 1: a cloud of cold atoms undergoes Bloch oscillations (e.g. due to gravity) inside a Fabry-Perot cavity, and the light field transmitted through the cavity is measured in order to determine the Bloch frequency. In order to quantify the measurement performance we will compute the signal-tonoise ratio using standard input-output theory [68].
Let us consider a double sided cavity with mirrors with matched reflectivities providing equal amplitude damping rates of κ/2. The quantum part of the input fields for both the top (driving side) and bottom (detection side) mirrors is given by the electromagnetic vacuum. Since we are not going to consider classical fluctuations of the driving laser we do not include a classical laser field contribution in the input field , but introduce it via the hamiltonian in Eq. (3). In our consideration of system dynamics in earlier sections we implicitly assumed a single sided cavity giving an amplitude damping rate of κ, and associated with this decay is a Langevin noise term √ 2κξ(t). In a double sided cavity we have two independent noise terms of the form √ κ{ξ t (t),ξ b (t)}. It can be shown that the dynamics of the intracavity system (both meanfield and fluctuations) are independent of whether we assume a double sided or single sided cavity as long as we divide the net damping equally among the two mirrors (provided they have matched reflectivities). The transmitted light field is the output field at the bottom mirror which is related to the input field at the bottom mirror asâ whereâ in andâ out in this equation refer to the fields at the bottom mirror. The transmitted photon current is given by the operatorÎ out (t) =â † out (t)â out (t), where againâ out refers to the field leaving the bottom mirror.
An experimentally straightforward method for measuring the Bloch frequency consists of recording the transmitted photon current using a photodetector. It is useful to consider the Fourier transform of the data [69] N (ω, T ) = T 0 dt cos(ωt)Î out (t) (38) and define the signal-to-noise ratio for the measurement as where ∆N 2 (ω, T ) ≡ (N − N ) 2 . Thus, the SNR is the ratio of the spectral density of the photon current to its variance and provides one measure of the sensitivity of the scheme. Let us first evaluate the SNR for a classical cavity field a(t) = α(t). In this case one finds that the signal ampli-tude and variance are given by In order to obtain an approximate magnitude for the SNR we further assume that the detection rate goes as ≈ R(1+ cos[ω B t]) [15], where is the contrast parameter defined in Eq. (14). Setting the classical photon current κ|α(t)| 2 in the above formulae equal to this detection rate gives Despite appearances, this result does include quantum noise to a certain degree because without the Langevin operators the variance given in Eq. (41) would have been zero, i.e. even when the cavity field is classical the output field contains a quantum partâ out = √ κα(t) −ξ b (t).
Thus, the above calculation includes detector shot noise, also known as measurement imprecision [64], but neglects the effect of quantum fluctuations on the coupled dynamics inside the cavity, i.e. quantum measurement backaction. Note that this is a different approximation from the coherent state approximation used in Section VII, where quantum fluctuations were included in the cavity dynamics by using a Glauber coherent state, i.e. a state with vacuum noise, for the cavity field, albeit one whose fluctuations are unaffected by the presence of the atoms. The SNR given by Eq. (42) predicts that the sensitivity of the scheme can be increased indefinitely by increasing the mean total number of photons collected RT and also the contrast . The former effect is the standard one expected from the general theory of measurements with uncorrelated fluctuations. The latter is intuitively plausible too, but, however, can not be the whole truth because, as stated above, it neglects the effect of measurement backaction upon the dynamics which is expected to become important at larger values of β. When fluctuations are includedâ(t) = α(t) + δâ(t), and the mean signal amplitude is given by In fact, this is not so very different from the meanfield photon number given by Eq. challenging to evaluate numerically not only because the fluctuations occur on time scales κ −1 much shorter than the BOs, but also because they require the storage and manipulation of data at two times. Furthermore, the continuous driving by the BOs means that the correlations are not stationary in time, i.e. they do not just depend on t 1 − t 2 , and this forces us to calculate the SNR in parallel to the system dynamics starting at t = 0. Unfortunately, due to limited computing power, we have only been able to track the SNR over ten Bloch periods which is certainly shorter than the coherence time of the BOs for the parameters we use. An actual experiment would, of course, not suffer from this limitation and would benefit from running until the BO coherence time is reached. The main steps of our algorithm for calculating the two time correlations are provided in Appendix C.

IX. SIGNAL-TO-NOISE RATIO: RESULTS
We now show how the SNR depends on the various system parameters. Due to the size of parameter space, this will not be an exhaustive study, but rather an ad hoc choice that nevertheless we hope is experimentally relevant. We begin by looking at the SNR as a function of the coupling parameter β = N U 0 /κ. In Fig. 8a we plot the SNR evaluated at ω B for an integration time of 10 Bloch periods. We change β by increasing U 0 but also change η to maintain the same minimum lattice depth of 3E R throughout. The results without measurement backaction (i.e. the dynamics in the cavity is purely meanfield) are plotted by the red (dots) curve which monotonically increases until about β = 12. The initial increase of the SNR with β is in line with expectations based on Eq. (42). The turnover of the red curve near β = 12 is in a sense an artefact that arises from having evaluated our SNR at ω B : it so turns out from the meanfield solution that for β > 12 the fraction of the power in the fundamental of s(ω) begins to decline and is diverted to higher harmonics. However, there is no real reason other than simplicity to only consider SNR(ω B ) (any harmonic of ω B gives information about the applied force and inclusion of all of them in the data analysis would extract the maximum possible information from the measurement). The full calculation including measurement backaction is plotted by the black (crosses) curve. The first thing to notice is that measurement backaction always lowers the SNR. Secondly, the full SNR monotonically increases only until β ≈ 7, and thereafter suffers from dramatic dips which we explain below as being due to resonances with quasiparticle excitation energies. These two observations are the main results of this paper. In Fig. 8b we plot the SNR as function of the total integration time T for three values of β, two before the first dip in the SNR and one in it. This plot further illustrates that for β > 7 there is a dramatic lowering of the SNR. The BO dynamics are also clearly visible due to the fact that the contrast is periodically growing and shrinking as the lattice depth grows and shrinks. In Figs. 9a and 9b we show how the SNR depends on other parameters, namely the lattice depth and total number of atoms. In particular, in Fig. 9a we plot the SNR as a function of the minimum lattice depth in the cavity for the coupling value N U 0 /κ = 1. The lattice depth is changed by increasing η. The red (dots) curve gives the SNR calculated using only meanfield dynamics plus the effect of shot noise at the detector and justifies the comment made in Section IV that for larger lattice depth the contrast decreases. The SNR calculation including measurement backaction fluctuations is given by the black (crosses) curve, and has the same qualitative behaviour but is somewhat lower. Fig. 9b plots the SNR as a function of N , where for different values of N we keep N U 0 /κ = 1 constant by scaling U 0 . We also scale the pump strength η to maintain the same intracavity lattice depth s(t) in all the cases. As we pointed out in Section IV, this method of scaling the system variables leaves the form of the meanfield and fluctuation equations unchanged. The only quantitative change is that the meanfield cavity field solution α(t) is scaled by the same √ r factor as the pumping. This leads to a linear scaling of the SNR as a function of N (with and without fluctuations) as shown in the plot. It is interesting to note that the rate of increase is different for the calculation including fluctuations compared to that without. Clearly there is a gain in the SNR with N . Finally, we shall explain the physical origin of the complicated series of dips in the SNR when β > 7 that are seen in Fig. 8a. Consider the spectrum of quasiparticle excitations about the adiabatic meanfield solution introduced in Section VI. For the example shown in Fig. 5a, the smallest excitation frequency occurs at the band edge q = ±1 and the largest at q = 0. As β is increased in the usual manner (holding the minimum lattice depth constant), the q-dependence of the quasiparticle spectrum evolves, as shown for the quasiparticle mode ω 1 in Fig. 10a. Thus, the range of frequencies (i.e. across the entire Brillouin zone) contained in ω 1 also evolves with β and is shown in Fig. 10b. If the meanfield dynamics happens to contain any frequencies that fall in this range there is clearly the possibility of a resonance, exciting quasiparticles and lowering the SNR. This is exactly what happens as can be seen from Fig. 11 which plots the total power in the harmonics of ω B that fall in the frequency range covered by ω 1 . The two peaks in Fig.  11 at β ≈ 8 and β ≈ 12 coincide exactly with the dips in Fig. 8a. Referring back to the inset in Fig. 3b, which was deliberately evaluated at β = 7.75 for this very purpose, we can see the part of the meanfield spectrum that falls in the range spanned by ω 1 . In the absence of BOs the quasiparticle excitation ω 1 is very narrow, with a width given by the imaginary part γ 1 evaluated at q = 0. However, the BO dynamics effectively broadens the resonance by orders of magnitude to that shown in Fig. 10b and this has a dramatic effect on the SNR.

X. DISCUSSION AND CONCLUSIONS
In this paper we have extended our previous analysis of Bloch oscillations of ultracold atoms inside a cavity to include the effects of quantum noise in the electromagnetic field. The quantum noise originates from the open nature of the cavity and can be interpreted as a form of quantum measurement backaction because it perturbs the dynamics. The magnitude of the backaction is controlled by the dimensionless atom-light coupling parameter β = N U 0 /κ and we find that it can strongly affect the sensitivity of a measurement of the Bloch oscillation frequency ω B and therefore the determination of the magnitude of the external force F driving them.
Our treatment is based upon the coupled Heisenberg equations of motion for the atoms and light which we linearize about their meanfield solutions, i.e. a Bogoliubov level approximation. We solve the time-dependent meanfield level dynamics exactly and hence coherent effects such as Landau-Zener tunneling between bands are fully taken into account. A spectral decomposition of the meanfield solution shows that it is dominated by ω B and its first few harmonics, but as β is increased spectral power begins to spread to higher frequencies.
Quantum noise is introduced via Langevin operators which act as inhomogeneous source terms in the Heisenberg equations. These terms excite quasiparticles (quantized excitations with a mixed atom-light character) out of the meanfield. In the standard situation [39,41,61] where there is no external force then if the system is started off with no quasiparticles their number initially grows in time but eventually saturates due to competition between cooling and heating processes (provided we are in the cooling regime ∆ eff c = ∆ c − N U 0 cos 2 (x) < 0 which means that the quasiparticle energy has a negative imaginary part). By contrast, in this work we have found that the presence of an external force, and hence BOs, profoundly changes this behaviour so that following some initial transients the heating rate settles down to a constant value even when we are nominally in the cooling regime. Nevertheless, for the parameter regimes we tested the heating rate was modest and the fraction of the atoms excited out of the coherent meanfield over the lifetime of the simulation was always less than 1% even for quite strong coupling.
In order to gain some insight into the numerical calculations we used Fermi's golden rule to develop a semianalytic model for the heating rate in terms of a simple rate equation for the number of atomic excitations. In so doing we approximated the cavity light field by a coherent state whose quantum fluctuations are the same as those of the vacuum. This is a common approximation in cavity optomechanics but ignores the quantum correlations that build up between the atoms and the light. Comparing this with the exact numerical results for the number of atomic excitations, we infer that the field is close to a coherent state for small β, but differs from it as β is increased, as expected. Furthermore, this comparison allowed us to see the dynamic generation of atom-light correlations.
The above calculations can be applied to the estimation of the signal-to-noise ratio for a continuous measurement of ω B . For example, we find that the SNR decreases with intracavity lattice depth, and increases with the number of atoms. Our principal result, however, concerns the dependence upon β. We find that the SNR can be severely reduced due to resonances between the quasiparticle spectrum and the Bloch oscillating meanfield for certain ranges of β. Indeed, the SNR behaviour depicted in Fig. 8 is much more complicated than that found in the standard example of a quantum limited position measurement of a harmonic oscillator, e.g. the end mirror of a resonant cavity [64]. In that system, the SNR is determined by the competition between the "measurement imprecision" (detector shot noise), which decreases with increasing measurement strength, and the measurement backaction, which increases with increasing measurement strength, and correlations between the two can be ignored to a good approximation. This leads to a smooth curve (see Fig. 5 on p.1171 of the review [64]) with a single maximum at the measurement strength where the two effects are equal. This is where the measurement should be performed for maximum sensitivity. By contrast, in our case we have a cloud of atoms occupying Bloch states in an optical lattice and thus our system does not correspond very well to a single harmonic oscillator (except in the limit where the lattice is extremely weak so that the atoms are predominantly in a state which is uniform in space [28], but then Landau-Zener tunneling will be so severe that the atoms will quickly fall out of the lattice when the external force F is applied). Add to this the fact that our system is driven by an external force and so scans through the entire Bloch band in a time-dependent fashion, leading to the possibility of resonances, and it is not surprising that our resulting SNR in Fig. 8a does not have a simple maximum as a function of β. However, we can make the parameter dependent statement that it seems safest to choose β < 7 which lies below the point where the resonances set in (and for β > 25 we find optical bistability which will destroy the Bloch oscillations [32]). The resonances only occur in the calculation when quantum measurement backaction is included and so provide a salutary example of when the latter is important. Nevertheless, away from the resonances the SNR for this continuous measurement is large and is in pretty good agreement with an approximate calculation based upon purely meanfield dynamics in the cavity plus detector shot noise.
In comparison to previously studied cold atom cavity-QED systems, or even cavity optomechanical systems, a new feature of our Bloch oscillating system is the timedependence of the meanfield. Apart from the resonances discussed above, this also has implications for the computational scheme we use to calculate the results. For example, all the fluctuation modes should be orthogonal to the meanfield mode as well as to each other, and hence they also evolve with time. Furthermore, the twotime correlation functions that are needed to calculate the signal variance that enters the SNR are not stationary in time, meaning that a large amount of data must be stored. This is especially true because the Bloch period is roughly three orders of magnitude larger than the quantum fluctuation timescale 1/κ and hence the calculation of the SNR over even a few Bloch periods is quite intensive in the regime where the coherent state approximation breaks down. In non-cavity BO experiments it has been shown that coherent dynamics can run for thousands of Bloch periods [9]. In a continuous measurement scheme, such as that proposed here, the quantum measurement backaction reduces the coherence time but unfortunately we have been unable to go much beyond ten Bloch periods with our numerical computations of the SNR and thereby find this coherence time for our scheme (we have, however, given an estimate in [15] based upon the idea that the spontaneous emission rate sets the upper limit on coherent dynamics). Nonetheless, our short time calculations illustrate quantitatively that it may be advantageous to remain at small β and integrate for longer times.

XI. ACKNOWLEDGEMENTS
This research was funded by the Natural Sciences and Engineering Research Council of Canada. We thank A. Blais, E. A. Hinds, J. Goldwin, J. Larson, and M. Trupke for discussions related to this work and D. Nagy for helpful correspondence.

Appendix A: Coherent State Approximation
In this appendix we provide details of how the coherent state approximation introduced in Section VII can be used to derive a simple rate equation for the occupa-tion numbers of atomic fluctuation modes. Working in the TF, consider the atomic fluctuation operatorδΨ ⊥ (t). Rather than expanding it in plane waves like in Eq. (21), let us instead expand in the instantaneous eigenbasis ν k (x, t) of the time-dependent meanfield hamiltonian H(t) [given in Eq. (10)] Substituting the decomposition Eq. (A1) into the equation of motion Eq. (22b) we obtain We see that the dynamics of the δb j (t) are coupled amongst themselves: this is obvious from the second term on the right hand side, but also occurs due to the third term as can be seen from Eq. (34). In order to obtain a description in terms of independent oscillators the contribution from these two terms must vanish, and we will now examine when this happens. We begin with the second term (with the time derivative) on the right hand side of Eq. (A3). It can be shown that [70]: In general, contributions to the above overlap element are suppressed for states well separated in energy due to the denominator. Also, we will show below that for k = j the element vanishes. Hence the dominant contribution comes from adjacent levels i.e. k = j ± 1 and is given by where the second line is obtained by taking a derivative of the instantaneous hamiltonian and realizing that, due to the opposing relative parity of adjacent states, the term in the overlap integral due to the potential is zero. The above term can be neglected if the Bloch frequency is small compared to the energy gap ∆ ± . To proceed further we will assume that this is the case but in the next appendix we will see that this cannot be guaranteed in general. Specifically, this approximation is most likely to break down at the times when the quasimomentum comes close to the center or the edge of the Brillouin zone where there are avoided crossings. Thus, this calculation will be valid only for short times (since at longer times the system will have repeatedly gone through such crossings) and/or at parameter regimes where the gaps are large compared to Bloch frequency. Coming back to the case when k = j we have ν k (t)|ν k (t) = 1 and so i.e. the derivative is purely imaginary. For the time dependent hamiltonianH(t) the potential term cos 2 (x) has an inversion symmetry about x = 0 and we can always choose the instantaneous eigenbasis ν k (x, t) to have real coefficients when expanded over plane waves. As a result the above term goes to zero and the second term in Eq. (A3) can be excluded. Turning now to the third term on the right hand side of Eq. (A3), we can see from Eq. (34) that it does not couple the different modes δb j (t) if the light field fluctuations are independent of the atomic fluctuations i.e.
which is exactly the coherent state approximation.
Having now seen the conditions under which the fluctuations in the instantaneous eigenmodes ofH(t) become independent, let us assume that these conditions are fulfilled so that the fluctuations obey the uncoupled equa-tions of motion These equations describe the atomic fluctuation dynamics in terms of a collection of independent oscillator modes that are acted upon by the shot noise forceF(t). As described in [63], we can now use Fermi's golden rule to derive a rate equation for each of the oscillator occupation numbers δN where the damping and diffusion rates are These depend on the spectral density (power spectrum) of the shot noise force In the above expressions the shot noise spectrum is evaluated at the shifted oscillator frequencies defined by ω j = E j (t)−µ(t) with the instantaneous chemical potential µ(t) = φ(t)|H(t)|φ(t) . This shifting helps in removing the slow time dependence of the couplings u j (derived from the meanfield Bloch oscillations). Since the damping and diffusion rates for the different oscillators are not the same, it is in general not possible to write down an equation similar in form to Eq. (A9) for the total δN (t), and we have to settle instead for δN (t) = j δN j (t).
Appendix B: Absence of cavity cooling in the presence of Bloch oscillations In this appendix we analyze the long-time behaviour of the number of atomic fluctuations δN . We do this in order to understand the apparent absence of a cavity cooling effect in the results shown in Fig. 7a. In the standard case where there is no external force [40,41], cavity cooling occurs when the effective detuning ∆ eff c ≡ ∆ c − N U 0 cos 2 (x) is negative. This ensures that the quasiparticle energies have a negative imaginary part γ n < 0 which implies dynamical stability as explained in Section VI. Under these circumstances δN reaches a steady state and the heating rate vanishes as shown in the inset in Fig. 7a. This is, however, not what we see in the presence of an external force as shown in the main body of Fig. 7a where the heating rate settles down to a constant nonzero value. The external force must therefore disrupt the cooling mechanism, and in this appendix we shall see that indeed the periodic driving due to the BOs drives the quasiparticles to higher energy states thereby heating the system.
The heating rate is given by the change in the occupation numbers of the various quasiparticle states as a function of time. These states are nothing but the instantaneous eigenvectors of the fluctuation matrix M(t) introduced in Section VI [40] M(t) r (n) (t) = (ω n (t) + iγ n (t)) r (n) (t) , and have a mixed atom-photon character. However, the fluctuation matrix M(t) is non-normal and so its left and right eigenvectors are not the same. The left eigenvectors l (n) are defined as The left eigenvectors can be used to define the quasiparticle mode operatorρ n (t) corresponding to the nth mode asρ where the bracket on the right in the above equation denotes a scalar product andR(t) = We therefore see that the required quasiparticle occupation numbers ρ † n ρ n (t) as a function of time can easily be computed from the numerical solution of the covariance matrix C(t) [Eq. (30)] once the eigenvectors l (n) (t) are obtained. Before we look at the results, we should first comment on the relation between the quasiparticle occupation number and the atomic fluctuation number δN (t). As mentioned in Sec. VI, quasiparticle modes come in three types and the most relevant ones are the hybridized atom-light modes which have the strongest atom-light coupling and tend to lie lowest in the spectrum. Since the hybridized modes contain both atomic and light components, their occupation number is not exactly equal to the atomic fluctuation occupation number. Nonetheless, in this system the atom-light entanglement is not very large [41] and the total quasiparticle occupation number closely tracks the atomic fluctuation number (as we have verified). Moreover, to establish a connection with the calculation in Appendix A, we note that for small N U 0 /κ the atomic part of the hybridised quasiparticle modes are very close to the higher band eigenstates of the instantaneous meanfield hamiltonianH(t). Thus, the mode occupation of the oscillators in Appendix A can be roughly mapped to the quasiparticle occupation numbers here.
In order to understand the occupation number dynamics, consider first the real part of the quasiparticle spectrum plotted in Fig. 12a  Since the gaps in the spectrum in (a) are smaller than the Bloch frequency the level populations are partially exchanged at the avoided crossings: the gaps become smaller higher up in the spectrum and indeed we see that the exchanges between higher lying states are almost complete. The system parameters are the same as the case with N U0/κ = 1 in Fig. 2 avoided crossings every half Bloch period which alternate between being with the band above and below. On the scale of the plot, the gaps at the crossings are not discernible but for the present parameters it turns out that even the gap between the lowest two bands is smaller than the Bloch frequency (recall that in this paper we have set ω B = 0.25ω R ) and the magnitude of the gaps gets smaller as we go higher up in the spectrum. During the course of Bloch oscillations these avoided crossings are repeatedly traversed at the Bloch frequency and consequently the occupation number dynamics at the avoided crossings are increasingly non-adiabatic as we go up in the spectrum due to Landau-Zener transitions. For example, in Fig. 12a at t = T B /2 the green (dot-dashed) curve of the third band approaches the blue (dashed) curve of the second band and as a result the populations of the two levels are almost completely exchanged as can be seen at the corresponding time in Fig. 12b. We therefore have the following picture: the occupation number of a given quasiparticle band increases either by direct scattering out of the meanfield due to quantum noise or by upcoming quasiparticles from the immediately lower band by a Landau-Zener transition. The occupation decreases due to quasiparticles scattering back into the meanfield [the hermitian conjugate term to the excitation processes in Eqns. (22a) and (22b)], or due to the finite lifetime of quasiparticles associated with cavity decay at rate κ as described by the A(t) term in Eq.
(22a), or due to Landau-Zener transitions to the next higher band. This has to be contrasted with the dynamics without BOs where the quasimomentum is fixed at q = 0 and the fluctuations occupy a stationary quasiparticle ladder. Without Landau-Zener transitions there is no directed transport of quasiparticles up the ladder and cooling effects, due to the finite quasiparticle lifetime 1/γ n , have time to act.
In Fig. 7a notice that the linear behaviour is established at later times for larger N U 0 /κ. In order to understand why this happens we explore the quasiparticle number dynamics for N U 0 /κ = 5 in Fig. 13, i.e. a factor of 5 greater than in Figs. 12a and 12b. From the inset we can immediately see that the two lowest quasiparticle bands are well isolated (by more than ω B ) from the rest of the ladder. As a result, the occupation numbers in these modes evolves in an adiabatic manner, in contrast to the situation for N U 0 /κ = 1. In fact, over the times plotted in Fig. 13, the blue (dashed) band reaches a steady average occupation number. But the higher quasiparticle energy levels represented, for instance, by the green (dotted) and black (dot dashed) lines have smaller gaps and behave akin to Fig. 12b because they are rapidly emptied by Landau-Zener transitions. Another relevant observation comes from Fig. 7b, where we see that the fluctuation photon number reaches its quasi-steady state around the same time as the atomic fluctuation number begins to exhibit linear growth. This can be understood shows the real part of the qp energy measured in units of ωR as a function of time over a single Bloch period TB for the lowest three bands. Since the gap between the lowest two bands (red (solid) and blue (dashed) lines) and the rest of the spectrum is larger than the Bloch frequency, their dynamics is decoupled from the rest. System parameters are as in Fig. 2.
now in the light of the above discussion since the lowest quasiparticle modes are coupled most strongly to the light field. The red (solid) band in the inset of Fig. 13 has two minima and demonstrates how for larger N U 0 /κ the quasiparticle bands can be strongly modified from the single particle (linear) band structure. We conclude this appendix by examining another way to control the band gaps in the quasiparticle spectrum and as a result the time taken for the linear increase behaviour (denoted by t l henceforth) to set in. In Fig. 14 we plot the atomic fluctuation number as a function of time for β = N U 0 /κ = 1 and three different initial meanfield lattice depth values that are set by the pump strength. Since the initial atomic state has q = 0, the initial lattice depth is the minimum lattice depth over the Bloch pe-riod. Furthermore, we are at relatively small β, and so the linear band picture holds good and one can anticipate that t l increases with lattice depths due to the widening of band gaps. In the inset in Fig. 14 we plot t l as a function of the initial lattice depth for a range of values at N U 0 /κ = 1. As expected, we see a general trend of increasing t l for larger lattice depths. We have identified t l from the numerical simulation for atomic fluctuation number by requiring that the average change in the rate of increase of δN (t) over a Bloch period converge to three When the intracavity light field is written asâ(t) = α(t) + δâ(t), the signal variance is given by: ∆N 2 (ω, T ) = κ T 0 cos 2 (ωt) |α(t)| 2 + δâ † δâ(t) dt + 2κ 2 Re T 0 dt 1 dt 2 cos(ωt 1 ) cos(ωt 2 )α(t 1 )α(t 2 ) δâ † (t 1 )δâ † (t 2 ) + κ 2 T 0 dt 1 dt 2 cos(ωt 1 ) cos(ωt 2 )α * (t 1 )α(t 2 ) δâ(t 1 )δâ † (t 2 ) + T 0 dt 1 dt 2 α(t 1 )α * (t 2 ) cos(ωt 1 ) cos(ωt 2 ) δâ † (t 1 )δâ(t 2 ) − 2(κ) 3/2 Re T 0 dt 1 dt 2 α * (t 1 )α(t 2 ) cos(ωt 1 ) cos(ωt 2 ) ξ b (t 1 )δâ † (t 2 ) θ(t 2 − t 1 ) . (C1) In this appendix we provide details of how we numerically compute the signal variance (and hence the SNR). The important extra computational step compared to the covariance matrix calculation in Eq. (29) is the evaluation of the two time correlations such as δâ † (t 1 )δâ † (t 2 ) . In the vector notation for the fluctuations, the two time correlations are elements of the correlation matrix Λ(t 1 , t 2 ) = R (t 1 )R T (t 2 ) . The time evolution for the correlation matrix is given by: Let us consider the case when t > t 0 . Then the last term in above equation gives a correlation between the Langevin operators at some future time t and the system fluctuation operators at t 0 . Due to the delta-correlated nature of the Langevin noise this term will be zero. This means that Eq. (C2) becomes homogeneous and we can solve it with the initial condition at t = t 0 , Λ(t 0 , t 0 ) = C(t 0 ). Also note that the time evolution operator for the numerical evolution in Eq. (C2) is same as the one for the covariance matrix [denoted by G(t) in Eq. (29)], which is an expression of the quantum regression theorem [68]. A separate computation for t < t 0 is not needed since they are related to the elements of Λ(t, t 0 ) with t > t 0 by complex conjugation. For example: δâ † (t 0 )δâ(t) = δâ(t)δâ † (t 0 ) * .
We can evaluate the correlation ξ b (t 1 )δâ † (t 2 ) using a similar approach as above for the time evolution of the vector ξ b (t 1 )R(t 2 ) . In this case the initial condition for the evolution is ξ b (t 0 )δâ † (t 0 ) = √ κ/2. Since the evolution operators for the correlation matrix and covariance matrix evolution are the same the calculation can be performed without additional computational cost. The main difficulty in computing the signal variance arises from the fact that the two time correlation functions are not stationary. As a result, in order to evaluate the integrals in Eq. (C1) the correlation matrix needs to be computed for all values of 0 < t 1 , t 2 < T . This is the memory intensive step in the computation and we simplify the situation by performing the correlation matrix computation over a coarser grid than the one used in the numerical solution of Eq. (29). This is justified since we find typically the correlation matrix elements do not change significantly over the very short time steps chosen in the solution of Eq. (29). Moreover, for the results presented in Section IX, we have taken care to check that the numerical solutions converge to a value independent of the size of the coarse grid. The necessity of evaluating two-time correlators over a two dimensional time grid is the main limiting factor to the maximum integration time for the SNR calculations. Another point to bear in mind is that for β values larger than the ones that we have presented here we have found that the size of the coarse grid needs to be essentially matched with the size of the finer computational grid over which Eq. (29) is solved. As a result the calculation for strong coupling becomes very memory intensive indeed.