Gravimetry through non-linear optomechanics.

Precision gravimetry is key to a number of scientific and industrial applications, including climate change research, space exploration, geological surveys and fundamental investigations into the nature of gravity.  A variety of quantum systems, such as atom interferometry and on-chip-Bose-Einstein condensates have thus far been investigated to this aim. Here, we propose a new method which involves using a quantum optomechanical system for measurements of gravitational acceleration. As a proof-of-concept, we investigate the fundamental sensitivity for gravitational accelerometry of a cavity optomechanical system with a trilinear radiation pressure light-matter interaction. The phase of the optical output encodes the gravitational acceleration g and is the only component which needs to be measured. We prove analytically that homodyne detection is the optimal readout method and we predict an ideal fundamental sensitivity of Δg = 10-15 ms-2 for state-of-the-art parameters of optomechanical systems, showing that they could, in principle, surpass the best atomic interferometers even for low optical intensities. Further, we show that the scheme is strikingly robust to the initial thermal state of the oscillator.

Introduction

The practise of measuring the gravitational acceleration g—also known as gravimetry—has led to important advances in both fundamental science and industry. For example, local gravity variations due to mass redistribution driven by climate change have been mapped with the GRACE satellite[1-3], and more recently, the Juno spacecraft mission reported the measurement of the gravity harmonics of Jupiter[4]. Furthermore, precise measurements of g can test for small deviations from Newtonian gravity on extremely small scales, which may provide indications of a deeper theory of quantum gravity[5]. In industry, precision accelerometry is extensively used in inertial navigation technologies and for conducting geological surveys.

While classical systems have long been utilised to perform accurate measurements of g, quantum systems offer several useful advantages, including reduced noise levels, a compact setup and most importantly an increased measurement sensitivity achieved through the power of coherence and interferometry. Over the past decade, a variety of quantum systems have been explored to this aim, in both theory and practice. The largest research effort to date has focused on atom interferometry[6-9], for which the highest achieved sensitivity currently stands at Δg = 4.3 × 10−9 ms−2[9]. A similar investigation has been carried out for both on-chip and fountain Bose–Einstein condensate (BEC) interferometry with best sensitivity Δg = 7.8 × 10−10 ms−2[10]. Finally, a proposal for using magnetically levitated spheres that predicts sensitivities of 2.2 × 10−9 ms−2 Hz−1/2 has been put forward in[11]. For comparison, the current commercial standard is set by the LaCoste FG5-X gravimeter which can achieve a measurement sensitivity of 1.5 × 10−9 ms−2 Hz−1/2[12]. More generally, the broader topic of using quantum systems to probe relativistic phenomena is currently being pursued with great interest (see for example[13-20]).

A key advantage to quantum systems are their interferometric properties. The following question arises: How can these interferometric properties be enhanced in order to improve the measurement sensitivity? One possibility is to place a quantum system in the form of a mechanical oscillator in an optical cavity, a research area known as quantum optomechanics[21]. See Fig. 1 for an illustration of a nanodiamond trapped in an optical cavity as an example of a class of optomechanical systems. The addition of the cavity allows for a strong coherent coupling between light and oscillator which, as we shall see, cancels out any initial thermal noise and fundamentally improves the measurement sensitivity of the device.

Fig. 1

An example of an optomechanical system. A nanodiamond is trapped and laser-cooled to milliKelvin temperatures in an optical cavity

Within classical optomechanics, the idea of gravimetry and accelerometry by optically detecting the mechanical oscillator has been experimentally realised by Cervantes et al.[22]. Other avenues, such as the detection of high-frequency gravitational waves through the driving of resonant mechanical elements was proposed also in[23]. In the related field of electromechanics, Schrödinger cat states and a Kerr nonlinearity have recently been found to be useful for the same applications[24]. However, the ensuing fundamental limits on the measurement sensitivity of gravimetry in the quantum regime of optomechanics using its trilinear radiation pressure interaction is yet to be investigated. Here we undertake this task and obtain some striking results: Firstly, it is possible, in principle, to surpass the sensitivity Δg that has been obtained in atom interferometers and other implementations to date. Secondly, due to the periodic decoupling of light and mechanics, the mechanical element does not require initial cooling to the ground state to improve the fundamental sensitivity of the gravimeter and, finally, the best possible sensitivity is achieved by a simple homodyne measurement of the cavity field, while only a low photon number in the cavity is required. That is, no measurement on the mechanical oscillator is required. Unlike the case of atomic interferometers, in optomechanics the interaction of light and matter is continuous, and we will see that our Hamiltonian cyclically entangles and disentangles the light and mechanics, leading to their decoupling. It follows that the experimental challenge will be to maintain the quantum coherence of the field and mechanics over the duration of each run of the experiment, which we set as one oscillation period of the mechanical element. This requirement, on which the plausibility of the scheme hinges, will be discussed in some detail.

Results

Optomechanics and Newtonian gravity

Let us begin by considering a general optomechanical system consisting of a mechanical oscillator coupled to a light-field in the cavity. The non-gravitational Hamiltonian that describes the dynamics of an optomechanical system is given by:[25,26]where are the annihilation and creation operators for the cavity field with frequency ωc, are the annihilation and creation operators for the mechanical oscillator with frequency ω, and k (usually denoted g in the literature, but this we shall here reserve for gravity) is a coupling constant that determines the interaction strength between the photon number and the position of the oscillator.

In order to introduce a coupling to a gravitational potential in the Hamiltonian, we add a term of the form . Here, m is the mass of the mechanical oscillator, g is the gravitational acceleration, is the position operator acting on the mechanical oscillator and θ is an angle from the vertical axis that we include in order to describe inclined systems, similar to[27]. Note that while the mass m appears as a coupling in the Hamiltonian, we will later see that measurements of g are mass-independent, which is what we expect from the equivalence of inertial and gravitational mass. With the addition of Newtonian gravity, the Hamiltonian of the system thus becomes

Optomechanical couplings

While we will keep the subsequent discussion general, let us here provide three examples of common optomechanical systems and their respective coupling constants. For a a Fabry–Perot cavity with a mechanical oscillator mirror, k takes the form[25,26]where L is the length of the cavity and m is the mass of the mirror. A levitated nano- or micro-crystal (e.g. a diamond or silicon bead), on the other hand, has a k given by[28,29]where ε0 is the permittivity of free space, Vc is the cavity mode volume, and kc is the wavevector of the laser, given by 2π/λ, where λ is the laser wavelength. P = 3Vϵ0(ϵ – 1)/(ϵ + 2) is the polarizability of the levitated object of volume V and ϵ is the relative electric permittivity. Alternatively, we can also consider a BEC trapped in a cavity. Here, the collective motion of the ensemble acts as the massive oscillator. For this system, the coupling constant is given by[30,31]where N is the number of atoms in the ensemble, g0 is the single-atom cavity QED coupling rate, M = Nm is the collective mass of all the trapped atoms with individual mass m, k is the wavevector of the laser and Δca = ωp – ωc with pumping frequency ωp. We will return to these expressions when computing the fundamental sensitivity limits for each system in the latter part of the paper.

System dynamics

In order to simplify the time evolution operator corresponding to the above Hamiltonian, we rescale by dividing all terms by the oscillator frequency ω. As a result, the time parameter t now represents the labframe time multiplied by ω, such that the oscillator has undergone a full oscillation cycle at t = 2π. The operator can then be written in the following decoupled form (see ref. [26] for details of the derivation in the absence of gravity):where r = ωc/ω, η = 1 – e−, , and . As a rule, we will denote any dimensionless quantity with a bar. For time-dependent variables, such as dissipation rates, this means they have been rescaled with respect to ω.

We now assume that the cavity field mode and the mechanical oscillator are initially in coherent states ǀ〉C and ǀ〉O respectively. For laser light injected into the cavity, this is the natural assumption. The oscillator, on the other hand, will in reality be initialised as a thermal state, which corresponds to a random coherent state ǀ〉O according to a thermal distribution. However, by starting out with a coherent state we will later argue that the gravimetric phase accumulated by the light does not depend on O so that our procedure works equally well for an arbitrary thermal state. A formal proof of this statement can be found in Supplementary Note 1. The initial state at t = 0 is then given by , and under it gives us the following statewhere τ = t – sin t, and ǀϕ(t)〉O are coherent states of the oscillator given by . In the derivation of this state, we have adopted a rotating frame for the cavity field, thus ignoring the free evolution induced by the term exp{–ira†a}.

The state in Eq. 7 show us that light and mechanics will entangle and disentangle periodically, with maximum entanglement occurring at t = π. At t = 2π, the oscillator state returns to regardless of the values of and β, and therefore by extension a thermal state also returns to its initial state because it will undergo the same compact evolution. This means that the initial oscillator state does not impact the fundamental sensitivity of this scheme. As already mentioned, a formal proof of this can be found in Supplementary Note 1. Most importantly however, at t = 2π the cavity state is completely decoupled from the oscillator, meaning that all information about g is transferred to the phase of the cavity state. As a result, any measurement scheme only needs to consider the cavity state after one oscillation period, meaning that direct or indirect access to the oscillator state is not required. This will greatly simplifies an experimental implementation, as measuring the oscillator state is generally difficult. This convenient property arises from the interferometric properties of the oscillator; its quantum nature allows it to acts as an interferometer to ensure that any initial thermal noise is removed from the cavity field, and thereby our scheme does not require cooling of the oscillator to a pure ground state. In other words, our results are valid for both coherent and thermal states. Note however that decoherence ensuing from damping to the oscillator motion during the state evolution will adversely affect the final measurement sensitivity and cause the oscillator state to grow increasingly mixed. We will not consider this kind of decoherence in this work, and instead assume that the mechanical element remains coherent over one oscillation period.

We can visualise some of the dynamics of the state in Eq. 7 by computing the expectation values of the field quadratures and [32]. We focus on the cavity state, which we obtain by tracing out the oscillator. The traced-out cavity state is given by

For decoherence-free evolution, the trajectories traced out by the system in phase space can be seen for different values of in Fig. 2a with and Fig. 2b with . Both figures are plotted with α = β = 1. We observe that the system performs increasingly complex trajectories for larger values of .

Fig. 2

Trajectories in phase space. Figure showing the position and momentum quadratures of the cavity state ρ that starts as a coherent state for different values of . As grows, we see the system perform increasingly complex trajectories in phase space. Parameters used here are α = β = 1, and a and b

We noted above that the light and mechanics periodically entangle and disentangle during its evolution. In order to see this more clearly, we can compute the linear entropy S(t) for the traced-out cavity state ρc(t) in Eq. 8. The linear entropy is defined as

The linear entropy tells us about the entanglement between the cavity and oscillator states. The results can be found in Fig. 3a for pure state and in Fig. 3b for states undergoing decoherence with photon dissipation rate . We see that S(t) increases until the state is maximally entangled at t = π. While a pure state completely decouples the light and mechanics at t = 2π for any values of and , a decohering state becomes increasingly mixed and does not return to its original state.

Fig. 3

Linear entropy of the traced-out cavity state. Plots showing the linear entropy S(t) for a free evolution and b noisy evolution with photon dissipation rate . Both cases have been plotted for values

Quantum metrology

We now come to our main results which concern the use of optomechanical systems as gravimeters. The question we wish to answer is: what is the best fundamental sensitivity Δg with which an optomechanical system can measure the gravitational acceleration g? Here, Δg denotes the standard deviation of a gravimetric measurement. We can directly predict Δg from the system’s dynamics by calculating the Fisher information IF(t) which provides a natural lower bound on the variance Var(g) of an unknown parameter, in our case g. This relationship is captured by the Cramér-Rao inequality[33-35]where N is the number of measurements. Thus if we maximise IF(t), we minimise the measurement spread of g.

Quantum Fisher information

The Fisher information comes in two forms: the measurement-specific classical Fisher information (CFI) and the quantum Fisher information (QFI). The QFI, which we denote HQ(t), is computed by optimising over all possible positive-operator valued measures (POVMs) and their resulting CFI[36]. Thus HQ represents the ultimate bound on obtainable information from a system, but it does not reveal which specific measurement is required to achieve it. For a general mixed quantum state ρ(t) the QFI is given by , where is the symmetric logarithmic derivative. In general, it is difficult to obtain analytically, especially for noisy systems. There are however methods for finding a noisy bound on the Cramér-Rao inequality[37]. A similar method for many-body systems was proposed in[38], and numerical methods were shown to be effective for a class of specific systems[39]. We shall not be using these methods here, as we shall instead investigate specific measurements for the noisy scenario to better approximate an experimental setting. This will later allow us to prove the optimality of the homodyne measurement.

Let us start by deriving a fundamental bound to the sensitivity. We specialise to the simpler case where the state ρ(t) is pure. Setting ρ(t) = ǀΨ(t)〉〈Ψ(t)ǀ the QFI becomeswhere we have used the notation ∂ = ∂/∂.

At first glance, the QFI of the global system might not seem very relevant as the mechanical part of the optomechanical system cannot easily be measured directly. However, we recall that the coherent state returns to at t = 2π, so that all information about g is transferred to the phase of the pure, decoupled cavity state. Since the decoupling time does not depend on β, this is also the case for a thermal state that may be written as a statistical mixture of coherent states (see Supplementary Note 1 for a proof of this statement). Calculating the QFI for this state will therefore provide an experimentally accessible notion of the fundamental sensitivity of the device. We find the following expression for HQ(t) at t = 2π:

Note that the mass term m is cancelled by the appearance of √m in the coupling constant , so that the final accelerometry measurement will be mass-independent. We also note the strong dependence on and ω, and that the expression scales linearly with the number of photons .

To find HQ(t) for the global state at any time t we resort to numerical calculations. We consider the case to allow for future comparisons with subsequent numerical evaluation of the CFI, which will be restricted to the same narrow parameter range. The resulting HQ(t) as a function of t can be found in Fig. 4a with and   β = 1 for various values of . We note that HQ(t) reaches its maximum value at t = 2π, which means that Eq. 12 returns the largest possible value during one oscillation period for any choice of system.

Fig. 4

Fisher information for measurements of gravitational acceleration g. Plots of the quantum Fisher information (QFI) and classical Fisher information (CFI) for measurements of g. a shows the QFI HQ(t) for 1, 4 and 9 photons with rescaled couplings . b shows the dimensionless CFI with and without photon dissipation rate  for a momentum measurement (λ = 1/2) and . c shows how the peak of the CFI for a momentum measurement narrows with increasing and with constant . d shows the IF(t) for a momentum measurement of photons that leak from the cavity with environmental coupling and . All plots use parameter   β = 1

Classical Fisher information

Let us now consider a specific measurement of g. The CFI IF(t) determines the minimum standard deviation of a parameter estimator once we have chosen a single specific measurement with POVM elements . The CFI is given by the expressionwhere is a conditional probability distribution.

We now consider a general homodyne measurement on the traced-out cavity state ρC. For notational convenience, we use a general Hermitian operator , where λ denotes a label that rotates between the field quadratures[40]. Any two operators that differ by λ = π/2 form a conjugate pair which satisfies the position-momentum commutator relation. In the following, we shall refer to the choices λ = 0 and λ = π/2 as a position and momentum measurement respectively. In order to calculate IF(t) we must find the probability distribution , where are the eigenstate of . While the position eigenstates themselves are not proper vectors, we can make use of a standard result from the quantum harmonic oscillator: [40], to writewhere H(x) are the Hermite polynomials of order n. These probabilities in turn gives rise to a CFI of the formwhere

Timescales of measurements

Let us analyse the expression for IF(t). We immediately note that any terms in the sum with n = n′ do not contribute to the Fisher information. The remaining behaviour of IF can be inferred from the second exponential in f, namely as this will dominate the entire expression for large . If we simplify the expression in the exponential, we find that it is equal to

For and large , the first term will dominate, and the exponential will be small for any t that is not a multiple of 2π. In other words, the Fisher information for a homodyne measurement becomes significant only when light and mechanics are completely decoupled. Figure 4c shows how the CFI for a momentum measurement (with λ = π/2) for and becomes increasingly narrow as grows larger. For clarity, we have rescaled IF with in the plot. Note that for small we still find large IF at times t ≠ 2π.

We saw earlier that the QFI scales with , which mean that the scheme favours systems with a large single-photon coupling. We shall soon show that the CFI coincides with the QFI at t = 2π, but in the meantime we must explore what the narrowing of the CFI at t = 2π entails for our measurement scheme. The narrow peak of the CFI will require the homodyne measurement to be performed within an increasingly narrow time-window. We can estimate the timescale in question by finding the full-width-half-maximum (FWHM) of the peak. To do so, we consider only the dominant first term for small perturbations in t around t = 2π, thus . That brings the first term into the form , which is now a Gaussian distribution. For a Gaussian function with , the FWHM is given by . In our case, we find . We already noted that terms with (n − n′) will not contribute to the CFI, and any term with will just cause the peak to narrow further. Thus we only consider the terms with , leaving us with , and so we conclude that any measurement must be performed roughly on a timescale of .

Optimality of homodyne detection

Let us see if we can simplify the expression for IF(t) even further and whether it bears any semblance to the QFI. At t = 2π, ϕ(2π) = β and η = 0. Then setting and to integer values causes IF(2π) to lose all dependence of . The coefficients reduce to and . We now consider the generating function for the Hermite polynomials . Taking the derivative of both sides results in , which we can use to show that Eq. 15 reduces to the compact expression

This expression coincides precisely with the QFI in Eq. 12 for complementary choices of λ and α. To better see why, we rewrite the term in the brackets as . We now note that when λ = 0, only contributes to the CFI, whereas at λ = π/2, only contributes. For both of these specific choices of λ, and when matched by α being either entirely real or imaginary, the CFI coincides precisely with the QFI in Eq. 12 because the term in the brackets reduces to or , respectively. We conclude that the homodyne measurement saturates the QFI limit up to a phase dependence of α, which can always be accounted for by changing the quadrature of the homodyne measurement. Note, however, that at other times than t = 2π. the homodyne measurement will be zero for all choices of λ and α, and so it only saturates the QFI when the light and mechanics have decoupled.

Finally, the absence of from IF(2π) is not a problem for sensing g—it just means that the sensitivity at times t = 2π is independent of the actual value of g. Numerical analysis suggests that larger values for causes the CFI to oscillate increasingly quickly before reaching its maximum value (see Supplementary Note 2). The optimality of the homodyne detection for sensing within our scheme is greatly advantageous as it is a routine measurement which is easy to accomplish. It has in fact also been shown to be an optimal measurement[41] in other contexts.

Decoherence

The calculation above is valid for pure states, but in practice every measurement will suffer various forms of decoherence. We will here investigate the effects of decoherence on the CFI for a narrow parameter range, as realistic parameters are very difficult to simulate numerically. We shall later use these results as indications of the behaviour of realistic systems.

There exists a large variety of decoherence effects for optomechanical systems, such as decoherence due to photons leaking from the cavity, or phonons gradually being lost from the mechanical element. The latter manifests as a gradual damping of the oscillator motion, which moves the state towards a mixture in the coherent state basis[42-45]. This problem has previously been treated analytically, which is possible because the decoherence operators commute with the Hamiltonian. Thus we refer to these works and will not treat the mechanical decoherence here. Instead, we make the assumption that the phonon decoherence is negligible over one oscillation period of the oscillator.

The effect of photons leaking from a cavity on a state ρ(t) can be modelled using a Lindblad master equation of the formwhere denotes the anti-commutator, are Lindblad operators and is the decoherence rate κ with respect to the rescaled time t. This equation cannot easily be solved analytically since the operator does not commute with the Hamiltonian in Eq. 2. Some solutions have been found for specific cases, for example when assuming that the photon leakage occurs only during the injection of the state into the cavity. The decoherence can then be modelled as a series of beamsplitters[46]. We will not consider these modifications here, but instead solve the Lindblad master equation numerically and compute the Fisher information IF(t) for the resulting mixed state.

In all subsequent numerical evaluations, we will set and α = 1 (note the choice of , which will optimise the CFI for λ = π/2). Larger values will cause the system to quickly grow numerically unstable due to the inclusion of non-linear terms such as (a†a)2 in the evolution in Eq. 6. While is experimentally achievable with the right choice of parameters, we can justify setting by noting that it physically corresponds to a heavily inclined cavity with θ ≈ π/2. Since we are interested in the general behaviour of the CFI under decoherence, we will here be working with the dimensionless Fisher information  (see the Methods section). Thus these numerical investigations should only be seen as a indication as to how decoherence will affect IF(t), and not as predictions for the sensitivity of a realised device. We shall later extrapolate from these results to make a prediction for realistic systems.

The behaviour of can be found in Fig. 4b for a momentum measurements with λ = π/2. A measurement with λ = 0 and α = i would show the same results. We note that larger values of do affect the CFI adversely, but setting implies that about 10% of the pure-state CFI is still accessible.

Measurements of leaking photons

In practise, a homodyne measurement is performed by monitoring and measuring the photons that continuously leak from the cavity. Aside from the experimental considerations, such a scheme also negates part of the photon dissipation considered above. We briefly estimate the CFI obtained through such a setup by using a simplified model where a pure vacuum state of the environment ǀ0〉E is added to our original state ǀΨ(t)〉CO, giving us the combined initial state . We then add a rotating wave interaction term to the Hamiltonian in Eq. 2, of the formwhere γ is the interaction strength and and are the creation and annihilation operators of the environment. The effect of this interaction Hamiltonian is to couple the cavity state to the environment which causes information about g to slowly leak out from the cavity into .

As before, we evolve the full state for a single-photon and with parameters . To maximise the CFI, we choose and λ = π/2. The results can be found in Fig. 4d for a rescaled coupling strength , where . As evident from Fig. 4d, we suffer a 10-2 reduction in the information that can be extracted from the system. Note also that the behaviour of IF(t) for this scenario will most likely also resemble a delta function centred around t = 2π for realistic parameters. Additional plots for this simplified model can be found in Supplementary Note 3.

Ideal sensitivities

In this section we shall first calculate the ideal Fisher information for the three optomechanical systems considered above, and then discuss the experimental challenges and advantages to an optomechanical gravimeter. As we here calculate the fundamental sensitivity, which is unlikely to be realised, we will only concern ourselves with order-of-magnitude estimates. These results are meant to showcase the potential of optomechanical systems, and to do so we have chosen state-of-the-art parameters that have been implemented in a variety of systems. For discussions of an experimental implementation including noise, see the Discussion.

Starting with the Fabry–Perot cavity system, we choose a fully vertical cavity with θ = 0 and use the following state-of-the-art experimental parameters: We choose a mass m = 10−6 kg, oscillator frequency ω = 103 Hz, cavity frequency ωc = 1014 Hz, cavity length L = 10−5 m and a photon number of . For these values, the rescaled coupling constant in Eq. 3 becomes , which gives us a Fisher information of IF = 1.58 × 1028 m−2 s4. This implies a sensitivity of Δg ≈ 7.96 × 10−15 ms−2.

Next, we look at a levitated micro-object confined in an ion trap interacting with an optical cavity, as demonstrated very recently in refs[47,48]. Again setting θ = 0 for maximal effect, we use mass m = 10−14 kg, oscillator frequency ω = 102 Hz, cavity frequency ωc = 1014 Hz, volume V = 10−18 m3, cavity mode volume Vc = 10−14 m3, electric permittivity ε = 5.7 for nanodiamonds, laser wavelength λ = 1064 × 10−9 m and a photon number of . From these values we obtain , which leads to IF  =1.15 × 1029 m−2s4. This gives us a final sensitivity of Δg ≈ 2.94 × 10−15 ms−2 for levitated nanospheres.

Finally, let us also consider cold atoms trapped in a cavity. Based on[30], we choose the following parameters: a wavelength λ = 780 nm, implying ωc = 1015 Hz, a single-atom coupling of g0 = 107 Hz, an atomic oscillation frequency ω = 102 Hz, a single-atom mass m = 10−25 kg, a detuning of Δ = 1011 Hz, and a laser wavevector of k = 108 m−1. With N = 105 atoms trapped in the cavity, we find that and IF  = 1.58 × 1019 m−2 s4, giving a sensitivity of Δg ≈ 2.5 × 10−10 ms−2. The reason for this disparity seems to be that the polarisability of the collection of cold atoms is not high enough to match the polarisability exhibited by the nanosphere. The number of trapped atoms can hardly match the number of atoms in a single nanosphere. One would either have to increase the number of atoms trapped in the cavity or increase the single-atom coupling strength to increase the Fisher information.

Comparison of theoretical results

Let us briefly compare the results obtained here with the performance of other quantum systems in the literature. In Table 1 we have listed a variety of experimentally implemented gravimeter systems with their best achieved sensitivity to date. Table 2, on the other hand, lists the ideal fundamental limits to sensitivities calculated in this work and others. The values for Δg and are presented in units of ms−2 and ms−2 Hz−1/2, respectively. The last column in Table 1 lists the integration time for each experiment, whereas in Table 2 the last column lists the experimental cycle time set by the oscillation frequency of the system in question. For atom interferometry, it is suggested in[49] that sensitivities of Δg ~ 10−12 ms−2 might be achieved, and a study of the fundamental limits has very recently been presented in[50].

Table 1

Comparison between gravimetry sensitivities obtained by various experimental systems

Experiments
System	Δg	Δg/Hz	Int. time
LaCoste FG5-X[12]	1 × 10−9	1.5 × 10−7	6.25 h
Atom intf[9].	5 × 10−9	4.2 × 10−8	100 s
On-chip BEC[10]	7.8 × 10−10	5.3 × 10−9	100 s
Optomech. accel[22].	3.10 × 10−5	9.81 × 10−7	10−3 sa

These include the commercial LaCoste FG5-X, atom interferometry, gravimetry through on-chip Bose–Einstein condensate (BEC) and classical optomechanical accelerometry. The second column lists the sensitivity Δg in ms−2 and the third column lists the -noise in ms−2Hz−1/2. The last column indicates the integration time needed to achieve each sensitivity

aThis value was provided to us by the authors of ref.[22]

Table 2

Comparison between sensitivities obtained by theoretical predictions for a variety of systems

Theoretical predictions
System	Δg	Δg/Hz	Cycle time
Magnetomech[11]	2.2 × 10−7	2.2 × 10−9	10−4 s
Fabry–Perot optomecha	10−15	10−16	10−3 s
Levitated optomecha	10−15	10−16	10−2 s
Cold atomsa	10−10	10−11	10−2 s

These include magnetomechanics, a Fabry–Perot optomechanical system, a levitated nanosphere optomechanical system and trapped cold atoms. The second column lists the sensitivity Δg in ms−2 and the third column lists the -noise in ms−2Hz−1/2. The last column indicates the cycle time or oscillation frequency ω for each system

aValues calculated in this work are denoted

Discussion

In this work, we investigated a new scheme for measurements of the gravitational acceleration g using a compact cavity optomechanical system with the usual trilinear optomechanical coupling to the cavity field. We derived a fundamental limit to the sensitivity Δg by computing the QFI and showed that the optimal sensitivity is achieved by a homodyne detection scheme performed on the cavity state at time  t = 2π. That is, no direct measurement of the mechanical oscillator is required. Using the expression in Eq. 13 and state-of-the-art experimental parameters, we predict a upper bound on the sensitivity of order Δg ~ 10−15 ms−2 for both a Fabry–Perot cavity and a levitated microsphere cavity, and Δg ~ 10−10 ms−2 for trapped cold atoms. These values compare favourably to all other currently available experimental and theoretical gravimetry proposals (see Tables 1 and 2). Furthermore, the quantum nature of the oscillator ensures that any thermal distribution in its initial state does not affect the fundamental sensitivity. However, as our scheme relies on superpositions involving distinct coherent states, we require thermal decoherence during one period of the oscillator motion to be negligible, which we estimate requires a Q-factor of at least 106 for the case of a Fabry–Perot cavity (see below). To explore the effects of photons leaking from the cavity, we numerically explored a narrow parameter range with , which physically corresponds to a nearly horizontally aligned cavity. We found that this form of decoherence does affect the system’s performance, but not severely. Finally, we briefly investigated what proportion of Δg we retain by performing measurements on the photons that leak from the cavity. Using a simplified noise model, we found a reduction of 10−2 in the resulting Fisher information. Given these results, we believe that there is significant potential in the use of quantum optomechanical systems for measurements of gravity and acceleration.

Let us now address some of the experimental challenges related to this scheme. Due to measurement inefficiencies and additional sources of decoherence not considered here, the final performance of optomechanical systems will naturally be expected to be lower than the values presented in Table 2. While we have shown that the initial optomechanical state does not need cooling to the ground state, thermal noise due to external influences during the evolution will gradually decohere the oscillator motion. We estimate that in the case of a Fabry–Perot cavity cooled to a temperature of milliKelvin, a number of ħω/(kBTth) = N phonons are present in the system at any time. Here, kB is Boltzmann’s constant and Tth is the system’s temperature. To retain coherence throughout the evolution, we require that , where κ is the phonon dissipation rate. In other words, the timescale of phonon decoherence κ must be much less than the characteristic timescale of the system. With ω = 1 kHz, as we assumed for Fabry–Perot cavities, we find N = 105 and κ = 10-2 Hz. A cavity which achieves such a decoherence rate must have a mechanical Q-factor of at least Q = ω/κ~106 to retain coherence, a regime which is not unprecedented.

Next, let us discuss which parameters yield the best sensitivities. Firstly, we note that the QFI in Eq. 12 ultimately scales with . In addition to the factor in the denominator, we acquire an extra from the rescaled coupling constant . The final factor of ω comes from the dependence of ω in k2. Given this scaling, we require ω to be as small as possible. At the same time, we also require the photon dissipation rate κ to be low. From our simulations, we saw that we require at least . This combination is difficult to achieve as low ω means the cavity must remain coherent over longer timescales. Therefore, the main experimental challenge of this scheme is to reduce ω and κ at the same time. Taking our numerical results as guidance, we essentially require that , which is nothing but the resolved side-band regime[28].

In the above, we used state-of-the-art parameters to calculate the ideal QFI for a variety of systems. However, as we just saw, the photon dissipation rate κ must be very low for these sensitivities to be achieved, and this has not yet been experimentally demonstrated for the parameters we used in Table 2. As technology improves we expect that this to be possible in future experiments, but for now, let us estimate the sensitivities that could be achieved today already. One of the best coherence times to date was demonstrated in[51], which achieved a cavity linewidth of κ = 660 Hz. To achieve a rescaled photon rate of for this system, we let ω = 6600 Hz and use L = 9.4 cm as reported in the paper. We keep m = 10−6 kg (since the QFI is ultimately independent of mass) and let ωc = 1014 Hz as before. Because the oscillation frequency ω is rather high, we choose to calculate IF for the Fabry–Perot cavity with a mechanical mirror, as this system performed slightly better for higher ω. The resulting coupling constant is , and the Fisher information is IF ≈ 2.16 × 1015 m−2 s4. This leads to Δg ≈ 2.15 × 10−8 ms−2. If we now assume that decoherence causes a similar proportion of the Fisher information to dissipate at these parameters compared to the ones chosen in our numerical simulations, we see that we retain about 10% of the pure-state Fisher information. Using this assumption, we find Δg ≈ 6.80 × 10−8 ms−2 and a -noise of 8.37 × 10−10 ms−2/. This is directly comparable with the values in Table 1, and so we believe that this scheme could be experimentally realised today, although the experimental challenges are of course substantial.

Let us briefly discuss ways in which we can decrease κ further and how this might affect the Fisher information. A heuristic estimate for κ can be given by considering the number of times per second that a single photon traverses the cavity. Each time the photon is reflected at the mirror, it has a T = 1—R chance of being transmitted instead of reflected. Here, T is the proportion of transmissions and R is the proportion of the number of reflections. The photon bounces off a mirror c/L times per second, where c is the speed of light. Thus we can take the dissipation rate to be κ = Tc/L, which means that increasing L decrease the photon dissipation rate κ, as the photon is effectively spending longer inside the cavity. However, increasing L also decreases the single-photon coupling constant, as we saw in the calculation above. This is true for all couplings we quote here, but it is perhaps most clearly seen for the case of the mechanical mirror and a Fabry–Perot cavity, with kFP given by Eq. 3. kFP scales with L−1, and so do the other couplings, through their dependence on the cavity volume Vc or the single-photon coupling g0. We recall that the Fisher information depends on , which means that it ultimately scales with L−2. Thus, changing L by an order of 10 will decrease the Fisher information by an order of 102. This contributes to the challenges of realising this scheme. However, it is important to note that there are realistic ways of increasing L without changing the single-photon coupling: One such method was explored in[52], where L was increased by adding an optical fibre to the cavity.

Furthermore, in the above we proved the optimality of a homodyne detection scheme, but we also found that such a measurement must be performed within a rather narrow temporal window, of timescale 1/k. Let us here estimate how quickly these measurements have to be performed based on the values we calculated for the coupling constant k. The nanospheres displayed the highest single-photon coupling with  Hz for the choice of ω = 102 Hz. Thus any homodyne measurement must be performed within 10−5 s, so we require at most microsecond precision, which is perfectly achievable. In comparison, we calculated for the levitated microsphere, which allows for a very comfortable ≈0.19 s window.

In spite of these challenges, optomechanical systems come with a number of advantages. They can remain stationary while performing the measurement, in contrast to on-chip BECs or BEC fountains which need to be launched, and the short cycle time of optomechanical systems allows for a large number of measurements to be performed very quickly. An additional point which we did not elaborate on above is that the spatial resolution of optomechanical systems will be extremely high since the oscillator is displaced only by a minuscule distance. As a result, it will be possible to determine very fine local variations in g, something which is not possible using larger systems. The scheme presented in this work also allows for the creation of macroscopic spatial superpositions, which, as pointed out in[11], is of great interest to testing gravitational collapse models (see for example[53-55]).

Before we conclude, let us now briefly discuss the underlying physical differences between atom interferometry and optomechanical systems for the purpose of gravimetry. We estimate that the QFI for atom interferometry is given by up to an unknown geometric factor, where n is the number of photons that deliver a momentum kick to the atoms, T is the total time over which the gravimetric phase is accumulated, and kc is the laser wavevector. See Supplementary Note 4 for the derivation. If we compare this to the Fisher information for optomechanical systems, we find that the Fabry–Perot cavity has a QFI that is larger by an enhancement factor , with c being the speed of light. This is due to the cavity confinement, whereby each photon interacts with the oscillator c/(2Lω) times per oscillation cycle, which is also the time period over which the gravimetric phase is accumulated. For the levitated nanosphere, we find a , where, again, for a micro-object containing ~1013 atoms, the polarizability P is much higher than that of a single atom. In practice, however, both of the enhancement factors will be damped by a factor ~1/(ωT)4 with respect to atom interferometry as the time of atomic interferometry T is typically larger than the time 1/ω of our scheme. Thus the sensitivity Δg is seen to improve by a factor of in our optomechanical scheme with respect to atomic interferometers. As n increases, the differences level out. However, different initial states in an optomechanical system, such as the superposition of two Fock states, will also scale with n2 (see Supplementary Note 5). Strictly speaking, the enhancement is valid for when the cavity field remains coherent for the time 1/ω over which our phase accumulation, i.e., κ«ω (the resolved side-band regime). However, our numerical results indicate that even in the presence of finite decoherence, say, κ~0.1×ω, the Fisher information is lowered only by a factor of about 10 compared to the case of loss-less cavities. Finally, we can also compare the treatment presented in this work to a position measurement of a classical oscillator that has been displaced due to gravity. While a classical treatment of the problem returns a preliminary measurement sensitivity similar to what we have derived in this work, it fails to take into account effects such as radiation pressure and the full quantum nature of the cavity field. Most importantly, a classical treatment does not utilise the coherent nature of the oscillator, which as we saw above negates any initial thermal noise in the state, and does not allow for the inclusion of other quantum states, such as squeezed states. After completing this work, the authors became aware of similar work carried out by Armata, Latmiral, Plato and Kim [58].

Methods

Dimensionless Fisher Information

For any numerical calculations, we must ensure that the numerical quantities are dimensionless. In order to calculate IF(t) under decoherence, we separate Eq. 13 above into a dimensionful and dimensionless part by writing

Here, is a dimensionful prefactor. The remaining integral

is now dimensionless and is what we will evaluate numerically. A final estimate for Δg can then be obtained by multiplying the value for by , but as this is only a rescaling we chose to only present the results for for clarity. Note that for the choice of , this dimensional prefactor is equal to .

Numerical methods

To evolve the system, we use the Python library Qutip[56] and a 4th order Runge–Kutta–Fehlberg method[57] for verification.

If we wish to compute the CFI for states undergoing decoherence, Eq. 15 is no longer valid and we must evolve the state numerically. We do so by computing the dimensionless part in Eq. 20 for a mixed state ρ(t). The probability distribution p(x|g) is easy to obtain numerically, since any operator has a finite matrix representation from which we can obtain the eigenstates and use these as our POVM elements. For example, we can define a position operator as a finite-dimensional matrix and solve for its eigenstates.

To obtain we must also compute the derivative . This can be done in a number of ways. The simplest one is to use a higher-order method of the central difference theorem. We obtained good and accurate results with the 4th order five-point method. For a function f(x) with parameter x and step-size h, the first derivative with this method is given by

As this method requires five data point, it is an expensive computation. This was our preferred numerical method as computing the CFI can still be done within reasonable timescales using the optimised master equation solver provided by the Qutip library. It does however contain two different sources of numerical errors: errors in the solver and errors that originate from the cut-off in the numerical derivative.

To verify that the error introduced by the numerical differentiation is not severely affecting the results, we used an additional method which provides an exact result. We reproduce it here for completion and in the hope that it might benefit others. We start by noting that as long as the POVM element П does not depend on , we can write the derivative as

Note that we are differentiating with respect to instead of g and that we have suppressed the dependence of t for clarity. This statement also holds for subsystems of ρ(g), which we can see by noting that the derivative distributes over a joint separable system as

Performing a measurement with П that only acts on subsystem A then gives

The second term reduces to zero because . While we have shown this for separable states, the same argument can be extended to entangled states by linearity.

In order to obtain the evolution for this state, we must now solve a modified version of the master equation. That is, given the Lindblad equation in Eq. 19,where denotes the anti-commutator, we now differentiate both sides with respect to to obtainwhere we have again used the notation . A more complicated form is obtained if the Lindblad operators depend on , which here is not the case. In coupled form, we can write

The system can be solved using any standard higher-order method, such as the family of Runge–Kutta ODE solvers. Note that the Qutip Master Equation solver cannot be used as Eq. 27 is not in standard Hamiltonian form.

Once the time-evolved state has been obtained, we proceed as usual to compute the probability distribution with the set of POVM elements {П}. With this method, we avoid round-off errors that appear in the five-point numerical derivative above.

Numerical stability

Let us make a few remarks regarding the numerical stability of the simulation. We start by considering the nature of coherent states and how they are represented numerically. Coherent states have support on infinite Hilbert spaces, whereas numerically we must work with finite matrices. It is therefore necessary to introduce a cut-off in the dimension used to represent the state. This leads to a gradual loss of coherence as information is pushed beyond the cut-off. In other words, numerically we use a finite Hilbert space , meaning that we truncate the space by letting , where . Furthermore, the appearance of in causes the system to become anharmonic and numerical instabilities grow fast for Hilbert spaces with small dimension N < 50.

The amount of information lost when using smaller Hilbert spaces is difficult to assess, since any good ODE solver will preserve the purity of the state throughout the simulation. Rather, the loss of information can be noted as a gradual deterioration of the trajectory in phase space, with the effect that states fail to return to their original position in phase space at t = 2π. That is, we require that and for the simulation to be deemed stable. This is however only true for noiseless evolution.

The system dynamics depend strongly on the dimensionless constants and . Larger and will cause the system to evolve more rapidly, as evident from their appearance in the phase of the state in Eq. 7. This in turn causes the numerical inaccuracies to accumulate more rapidly. When computing the CFI for mixed states, we restrict our investigations to the parameter range for precisely this reason.

Finally, it should be noted that we have not provided a full error estimate for any of the results computed here. However, since we are only interested in the general behaviour of the CFI, any small inaccuracies in the numerical estimates will not matter for the results obtained in this work.