Analogue quantum simulation of the Hawking effect in a polariton superfluid

Abstract

Quantum effects of fields on curved spacetimes may be studied in the laboratory thanks to quantum fluids. Here we use a polariton fluid to study the Hawking effect, the correlated emission from the quantum vacuum at the acoustic horizon. We show how out-of-equilibrium physics affects the dispersion relation, and hence the emission and propagation of correlated waves: the fluid properties on either side of the horizon are critical to observing the Hawking effect. We find that emission may be optimised by supporting the phase and density of the fluid upstream of the horizon in a regime of optical bistability. This opens new avenues for the observation of the Hawking effect in out-of-equilibrium systems as well as for the study of new phenomenology of fields on curved spacetimes.

Graphic abstract

1 Introduction

Analogue gravity is a type of analogue quantum simulation that enables the laboratory study of quantum field theories on curved spacetime via the equivalence between the kinematics of excitation in material systems and of massless fields in astrophysics [1, 2]. For example, a one-dimensional trans-sonic fluid flow that is a fluid whose flow velocity goes from being sub- to super-sonic forms an acoustic horizon where the flow velocity of the fluid equals the speed of sound [3]. In this configuration, excitations of the acoustic field behave as though they propagated on an effectively curved-spacetime whose properties are directly controlled by the geometry of the fluid flow [4], thus simulating a Black Hole. Importantly, quantum fluctuations of the vacuum scatter at the acoustic horizon, causing the emission of correlated waves on both sides of the horizon by the Hawking effect (HE) [3, 5].

Experimental evidence for correlated emission by the HE was recently reported in analogue gravity setups based on classical [6, 7] and quantum fluids [8]. While the thermal fluctuations of classical fluids overpower quantum fluctuations at the horizon such that vacuum emission cannot be observed there, this can be done with quantum fluids. Vacuum emission would yield a non-separable state at the output [9,10,11,12,13,14,15,16,17,18], whose degree of entanglement could be quantified from the density and correlation spectra [19, 20].

Although most theoretical work on the HE has been dedicated to analogues based on atomic Bose–Einstein condensate (BEC) [21,22,23,24,25,26,27,28,29,30], correlated emission with comparable properties may also be observed in quantum fluids of microcavity polaritons [31,32,33] where a sonic horizon has already been experimentally realised in one- and two-dimensional flows [34, 35]. In both quantum fluids, the HE manifests itself as the emission of collective Bogoliubov excitations that propagate in opposite directions on either side of the horizon. However, to date, theoretical works on polariton analogues have predicted a signal that appears hardly measurable because of its low strength (correlations on the order of $5\times 10^{-5}$ [31,32,33]) and short propagation length (about ${12}\,\upmu {\hbox {m}}$ from the horizon for the experimental configuration of [34]).

The main difference between quantum fluids of atoms and of polaritons is that the latter are intrinsically out of thermal equilibrium. Radiative and nonradiative dissipative processes in microcavities must be compensated for by optical pumping, and so the non-equilibrium state is not determined by thermodynamic equation conditions [36]. The radiative decay of polaritons does not solely render real-time, in-situ diagnosis of the fluid properties possible (a notable experimental simplification compared with Bose–Einstein condensates of atoms), it also is at the origin of a unique phenomenology in the collective dynamics. Specifically, in the regime of high fluid density of interest to horizon physics, a gap may open between the dispersion relation and the frequency of the pump [37, 38]. Here we show how this affects the strength of the HE and how, in turn, two-point correlations can be used as a diagnostic for out-of-equilibrium effects.

Fig. 1

Waterfall geometry with a polariton fluid. a Sketch of the system: the pump spot is structured in a two-steps profile: on the first step (far left), the fluid is set above the bistable regime ($|F_p|$=9). On the second step, the fluid density is supported near the turning point of the bistability loop ($|F_p |$=1.2, point C in c). The angle $k_p$ of the pump field creates a transsonic fluid flow across an attractive obstacle at $x=0$. b Bistability loop. c Top: pump profile near the defect ($|F_p |$=1.2). The pump strength is set to zero ${10}\,\upmu {\hbox {m}}$ before the defect. Bottom:Spatial properties of the fluid when setting the fluid density and phase in the upstream region as close as possible to point C of the bistability loop. Black, pump intensity; red, fluid velocity; blue, speed of excitations

In this paper, we explore the parameter space of quantum fluids of polaritons and identify a regime that is specifically favourable to the formation of correlations by the HE. The hydrodynamics of the fluid are controlled by its density and phase profiles, which are in turn connected with the optical bistability of the system (the hysteresis cycle of its polariton-density-to-optical-power relationship) [39], and so we investigate vacuum emission by the HE from this perspective. We study the influence of the fluid density on either side of the horizon on vacuum emission. In doing so, we explain how out-of-equilibrium effects in the configurations considered in [32,33,34] limit the emission of Bogoliubov excitations, and we show how to engineer the fluid such that the strength and spatial extension of the correlation signal become amenable to experimental detection. We obtain an order of magnitude increase in the length and a factor 4 in the strength of quantum correlations from vacuum fluctuations compared to [32,33,34]. Fine control upon the working point provides us with a better understanding of the influence of the properties of the quantum fluid of polaritons on emission by the HE in systems out of equilibrium. Our results open the way to the experimental observation of the HE in polaritonic systems.

2 Polariton fluid and analogue gavity

In order to reproduce the kinematics of a scalar quantum field near an event horizon, we use the model of the so-called waterfall geometry illustrated in Fig. 1 in laboratory frame coordinates x and t. This flow profile is realised in a polariton wire, a laterally confined microcavity in which the polariton dynamics are effectively one-dimensional. The microcavity is pumped with a continuous wave, coherent pump laser incident at a given angle with respect to the normal to form a stationary flow along the wire. The light field is structured in a step-like intensity profile (black line in Fig. 1c). The first step of high intensity excites the fluid in the nonlinear regime while the second step supports the fluid density at a tunable working point. As in [33, 34], the cavity features an attractive defect (a localised ${1}\,\upmu {\hbox {m}}$ long broadening of the wire) placed after the region where the pump lies.

Our system is a one-dimensional quantum fluid of exciton-polaritons whose flow velocity goes from being sub- to super-sonic, thus forming a sonic horizon where the local flow velocity of the fluid equals the local speed of sound. For now, we begin with the theoretical description of a homogeneous quantum fluid and of the propagation of quantum fluctuations of its phase and density.

Exciton-polaritons are quasi-particles resulting from the interaction of light with matter in a semiconductor microcavity. Photons emitted by a laser are sent into a cavity formed by two Bragg mirrors, wherein their dispersion is the usual Fabry–Perot dispersion, giving them an effective (very low) mass. These trapped photons create excitons, which are bound electron-hole pairs, in the semiconductor microcavity. Strong coupling between the photons and excitons trapped in quantum wells gives rise to two eigenstates for the total Hamiltonian, known as the lower polariton (LP) and upper polariton (UP) branches, separated by the Rabi splitting. Furthermore, the Coulomb interaction between excitons results in an effective non-linearity for exciton-polaritons (polaritons). The dynamics of the mean-field are governed by a generalised Gross–Pitaevskii equation, which leads to Euler and continuity equations describing the system as a quantum fluid. Historically, polaritons have first been described as two-dimensional quasi-particles [37], although the theory may be reduced to one-dimensional cavities called wires [32,33,34], as in the present case.

In our case, all energies involved are small compared to the Rabi splitting so the exciton-polariton system can be described by the mean field approximation [36]. At this level the system is described by a single scalar field $\Psi$, the field of lower polaritons, whose dynamics are governed by the driven-dissipative Gross–Pitaevskii equation (GPE)

$$\begin{aligned} i\partial _t \Psi (x,t)= & {} F_p(x,t)\nonumber \\{} & {} +\left( \omega _0 \!-\! \frac{\hbar }{2m^*} \partial _x^2 \!+\! V(x) \!+\! g|\Psi (x,t)|^2 \!-\! i\frac{\gamma }{2}\right) \Psi (x,t).\nonumber \\ \end{aligned}$$

(1)

$F_p$ is the field of the pump laser, $\omega _0$ is the frequency of the lower polaritons at the bottom of the branch, $m^*$ is their effective mass, V is the ‘external potential’ (that is controlled via the cavity geometry), g is the effective non-linearity, $\gamma$ is the loss rate. The field $\Psi (x,t)$ is written in the laboratory frame.

In order to make the link between this Gross–Pitaevskii equation and the description of the polariton ensemble as a fluid, we perform the Madelung transformation: we write the field of lower polaritons as $\Psi (x,t) = \sqrt{n(x,t)} e^{i\theta (x,t)}$ and insert this expression into Eq. (1). We take the real and imaginary parts to arrive at the Euler and continuity equations for the polariton fluid [36]:

$$\begin{aligned}&\partial _t \theta + \frac{m^* v^2}{2\hbar } + \frac{\hbar }{2m^*}\frac{\partial _x^2\sqrt{n}}{\sqrt{n}} + V + gn + \frac{\Re {F_pe^{-i\theta }}}{\sqrt{n}}= 0,\nonumber \\&\partial _t n + \partial _x (n v) = \gamma n - 2\Im {F_pe^{-i\theta }}\sqrt{n} \end{aligned}$$

(2)

with

$$\begin{aligned} v = \frac{\hbar }{m^*}\partial _x \theta . \end{aligned}$$

(3)

The first and second equations of (2) correspond to the Euler and continuity equations of atomic Bose–Einstein condensates (BECs), albeit with terms coming from coherent pumping, losses and quantum pressure. We see that the properties of the fluid depend on two parameters, namely its density $n=|\Psi |^2$ and phase $\theta$. The spatial variations of the phase are encapsulated in v, which we identify from Eq. (3) as the flow velocity of the fluid. Via n and $\theta$, the geometry of the flow may be all-optically controlled by engineering the profile and phase of the pumping laser.

Let us first consider the case of a non-dissipative, non-driven fluid (leave out $\gamma$ and $F_p$). We also neglect the quantum pressure so Eq. (2) become

$$\begin{aligned} \begin{aligned}&\partial _t \theta + \frac{m^* v^2}{2\hbar } + V + gn = 0\\&\partial _t n + \partial _x (n v) = 0. \end{aligned} \end{aligned}$$

(4)

We investigate small excitations in this system, i.e. low-k phononic modes on top of the fluid (see Sect. 2.1 for details). The linearisation of Eq. (4) around a background state yields a wave equation that is strictly isomorphic to the wave equation of a massless scalar field $\Psi _1$ (the Klein–Gordon equation) on a 1+1D curved spacetime,¹$\Delta \Psi _1\equiv \frac{1}{\sqrt{-\eta }}\partial _\mu \left( \sqrt{-\eta }\eta ^{\mu \nu }\partial _\nu \Psi _1\right) =0$, with the metric tensor

$$\begin{aligned} \eta _{\mu \nu }=\frac{n}{c^2}\begin{pmatrix} -\left( c^2-v^2\right) &{} -v\\ -v &{} 1 \end{pmatrix}, \end{aligned}$$

(5)

$c\propto \sqrt{n}$ the ‘speed of sound’ and $\eta =\textrm{Det}(\eta _{\mu \nu })$. The various components of this ‘acoustic metric’ [3] are given by the fluid velocity. Notably, there is an event horizon where $v=c$ (the time component $\eta _{00}$ of the metric goes to zero) [4]. So optically controlling the flow velocity and the speed of sound permits the engineering of an effectively curved spacetime for acoustic waves. Crucially, this is also valid for the quantized acoustic field, such that correlated pairs are emitted at the horizon by the Hawking effect [5, 21].

In paragraph 2.5 we will show that the driven-dissipative nature of polaritons does not modify the kinematics of excitations. On the contrary, these additional interesting properties can be harnessed to reveal new effects.

Fig. 2

Bogoliubov dispersion. a Dispersion in the co-moving frame of the fluid at point C in Fig. 1. Real part of the dispersion (11) for a pump vector $k_p = {0.25}{\,\upmu {\hbox {m}}^{-1}}$. Blue, $\omega _+$, positive-norm modes; orange, $\omega _-$, negative norm modes. Red dashed lines, speed of sound. The dispersion in the laboratory frame is obtained by Doppler-shifting the fluid frame dispersion: b Laboratory-frame subsonic dispersion exactly at point C (Eq. (11)). c Laboratory-frame subsonic dispersion away from point C on the upper branch of the bistability loop (Eq. (13)). d Laboratory-frame supersonic dispersion of a ballistic fluid (Eq. (14)). Blue, positive-norm branch; orange, negative-norm branch. Filled dots, local modes of positive group velocity; circles, local modes of negative group velocity. $\omega _\textrm{min}$, lower frequency of the gapped positive-norm branch; $\omega _\textrm{max}$, upper frequency of the negative-norm branch

2.1 Bogoliubov excitations

We now study the propagation of small fluctuations such as quantum fluctuations in this fluid using the linearisation of the Gross–Pitaevskii equation with the Bogoliubov method.

Assuming the pump beam to be monochromatic with frequency $\omega _p$, $F_p(x,t)=F_p(x)e^{-i\omega _pt}$, we write the polariton field as $\Psi (x,t)=\Psi (x) e^{-i\omega _p t}$. In the steady state the Gross–Pitaevskii Equation (1) becomes

$$\begin{aligned}{} & {} \left( \omega _0 - \omega _p - \frac{\hbar }{2m^*} \partial _x^2 + V(x)\nonumber \right. \\{} & {} \quad \left. + g|\Psi (x)|^2 - i\frac{\gamma }{2}\right) \Psi (x) + F_p(x) = 0, \end{aligned}$$

(6)

We first consider a configuration where the wire is pumped with a spatially homogeneous and monochromatic pump of incident wavevector $k_p$. The phase gradient of the fluid is then set by, and equal to, $k_p$ while its density is homogeneous. The steady-state GPE (6) simplifies to

$$\begin{aligned} \left( g|\Psi |^2 -\Delta _p - i\frac{\gamma }{2}\right) \Psi + F_p = 0, \end{aligned}$$

(7)

where $\Delta _p$ is the effective detuning defined as the difference between the pump energy and that of lower polaritons,

$$\begin{aligned} \Delta _p = \omega _p - \omega _0 - \frac{\hbar k_p^2}{2m^*}. \end{aligned}$$

(8)

Bogoliubov excitations are then mathematically obtained by linearising the GPE (6) around a background: $\Psi \rightarrow \Psi + \delta \Psi$, and $\Psi ^* \rightarrow \Psi ^* + \delta \Psi ^*$. The dynamics of the excitations $(\delta \Psi , \delta \Psi ^*)$ is given by the Bogoliubov matrix:

$$\begin{aligned} i\partial _t\begin{pmatrix}\delta \Psi \\ \delta \Psi ^* \end{pmatrix}= \mathcal {L}\begin{pmatrix}\delta \Psi \\ \delta \Psi ^*\end{pmatrix}. \end{aligned}$$

(9)

We go to the reference frame co-moving with the fluid via a Galilean transform ($x \rightarrow x -vt$). In the special case of a homogeneous system where the interaction energy matches the detuning, $gn=\omega _p-\omega _0-\frac{\hbar k_p^2}{2m^*}$, the Bogoliubov matrix $\mathcal {L}$ can be written in this frame as

$$\begin{aligned} \mathcal {L} = \begin{pmatrix} gn + \frac{\hbar k^2}{2m^*} - i\gamma /2 &{} gn\;e^{2ikx} \\ -gn\;e^{-2ikx} &{} -gn - \frac{\hbar k^2}{2m^*} - i\gamma /2 \end{pmatrix}. \end{aligned}$$

(10)

where k is the wavenumber of Bogoliubov excitations in the co-moving frame. Upon diagonalization, we retrieve the Bogoliubov dispersion relation in this co-moving frame, which relates k to the frequency of the Bogoliubov excitations $\omega$:

$$\begin{aligned} \omega ^{\Delta _p=gn}_\pm = \pm \sqrt{\frac{\hbar k^2}{2m^*}\left( \frac{\hbar k^2}{2m^*} + 2gn\right) } - i\gamma /2. \end{aligned}$$

(11)

Figure 2a shows the real part of Eq. (11), the dispersion curve, which, in this case where $\Delta _p=gn$, is identical to that of atomic BECs. There are two branches $\omega ^{\Delta _p=gn}_\pm$ of the dispersion, which are symmetrical around the point $\omega =0,\,k=0$. At low k, the dispersion curve has a linear slope: $\omega ^{\Delta _p=gn}_\pm \xrightarrow [k\rightarrow 0]{} \pm c_B k$, with $c_B$ is the speed of excitations in the fluid, which is also the speed of sound in the fluid.

$$\begin{aligned} c_B = c_s = \sqrt{\hbar gn/m^*} \end{aligned}$$

(12)

At large k, the dispersion is that of free massive particles, $\omega ^{\Delta _p=gn}_\pm \xrightarrow [k\rightarrow \infty ]{} \hbar k^2/2m^*$. There, $|{\partial \omega ^{\Delta _p=gn}_\pm }/\,{\partial k}|> c_s$ the gradient of the dispersion curve is larger than the speed of sound.

2.2 Acoustic horizon

The polariton ensemble behaves as a fluid whose dispersive properties depend on its density and velocity, and we will use these properties to study the waterfall geometry, which is the basis of our black hole analogue model. The waterfall separates two regions of fluid density n and phase $\theta$. The polariton fluid flows across a defect close to $x=0$ in the positive x direction at velocity $v = \frac{\hbar }{m^*}\partial _x \theta$, as shown in red in Fig. 1c. The speed of excitations in the fluid depends on their wavenumber. As mentioned above, Bogoliubov excitations have a group velocity $\nicefrac {\partial \omega }{\partial k}$. However, the relevant case for our model is when the excitations have a small wavenumber k and their speed, given by Eq. (12), does not depend on k and is proportional to the square root of the fluid density (in blue in Fig. 1c).

As can be seen in Fig. 1c, the speeds of the fluid and of the excitations, calculated from the GPE equation Eq. (1), vary a lot when the flow goes across the defect, which creates an attractive potential $V(x)=-{0.85}\,\upmu \hbox {e}\hbox {V}$ with a Gaussian distribution centred at $x=0$. Upstream of the defect, for $x<{-2}{\,\upmu {\hbox {m}}}$ we have $v<c_B$, the flow is subsonic while downstream of the defect, for $x>{-2}\,\upmu {\hbox {m}}$, $v>c_B$, the flow is supersonic, so there is an acoustic horizon at $x_H={-2}\,\upmu {\hbox {m}}$. The region where the flow is subsonic is outside the horizon, the region where the flow is supersonic is inside the horizon.²

In the next section, we will focus on the properties of the fluid in the waterfall configuration and the kinematics of collective (or Bogoliubov) excitations therein. As can be seen in Fig. 1c, the flow velocity v is rather flat, except in the near-horizon region. In particular, v spikes around $x=0$ because of approximate conservation of the flow current there (under dissipation). Beside this narrow feature, the flow velocity may be treated as homogeneous on either side of the region around $x=0$. Meanwhile, $c_B$ is flat in the pumped region, slightly bumps between the edge of the pump and the defect as well as just after the defect (under the approximate conservation of the flow current mentioned above) and then decreases as the polariton population decays. For the sake of simplicity, we may consider that the speed of excitations is homogeneous on either side of the region around $x=0$.

We will first discuss in more detail the dispersive properties of the fluid as a function of the parameters of the polariton field and identify the plane waves of this field before explaining how to construct “global modes” of the inhomogeneous system (including the horizon). These are the modes in which correlated emission from vacuum occurs. We then show the influence of out-of-equilibrium physics on the kinematics of Bogoliubov excitations on either side of the horizon.

2.3 Dispersion relations and optical bistability

In this section, we treat the system as composed of two regions homogeneous in their properties on either side of the horizon at $x_H={-2}\,\upmu {\hbox {m}}$. The properties of the polariton field in the upstream region are set by those of the pump, and they are given by Eq. (6), with no external potential ($V(x)=0$). By linearizing Eq. (6) in these conditions, we obtain the dispersion relation in the laboratory-frame:

$$\begin{aligned} \omega _{\pm } (k) = \pm \sqrt{\left( \frac{\hbar \delta k_p^2}{2m^*} - \Delta _p + 2 gn\right) ^2 - (gn)^2} + v \delta k_p - i\gamma /2\nonumber \\ \end{aligned}$$

(13)

with $\Delta _p :=\omega _p - \omega _0 - \nicefrac {\hbar k_p^2}{2m^*}$, the effective detuning between the pump energy $\hbar \omega _p$ and that of polaritons flow. $k_p$ is the wave-number of the pump field, $\delta k_p=k-k_p$, and g the interaction strength. The term $v\delta k_p$ is due to the Doppler effect in the laboratory frame [37]. Note that unlike the configuration considered in the previous section, we do not assume that the interaction energy matches the effective detuning, and the dispersion curve is thus modified.

On the other hand, in the downstream region, the pump field is zero, so the polariton fluid propagates ballistically there [42] and the dispersion is

$$\begin{aligned} \omega _\pm (k) = \pm \sqrt{\frac{\hbar \delta k_0^2}{2m^*}\left( \frac{\hbar \delta k_0^2}{2m^*} + 2gn\right) } + v\delta k_0 - i \gamma /2, \end{aligned}$$

(14)

with $\delta k_0=k-k_0$ where $k_0 = \nicefrac {m^* v}{\hbar }$ denotes the wave-number of the ballistic fluid.

In the upstream region, due to presence of the pump laser, the relationship between the density of polaritons n and the intensity of the pump laser $|F_p|^2$ is obtained from squaring Eq. (7)

$$\begin{aligned} \left( \left( gn - \Delta _p\right) ^2 + \frac{\gamma ^2}{4}\right) n = |F_p|^2 \end{aligned}$$

(15)

In the case where the energy of the laser is above that of the lower polaritons, $\Delta _p > \gamma \sqrt{3}/2$, this equation has several solutions in n for a given value of $F_p$ [43]. This degeneracy of fluid densities is due to optical bistability [39]. As shown in Fig. 1b, there is a hysteresis relationship between the fluid density n and the pump strength $|F_p|$, which comes from the bistability loop (cf “Appendix A”).

The previous equation can be written in terms of the velocity of excitations $c_B$.

$$\begin{aligned} \left( \left( \frac{m^* c_B^2}{\hbar } - \Delta _p\right) ^2 + \frac{\gamma ^2}{4}\right) \frac{m^* c_B^2}{g\hbar } = |F_p|^2. \end{aligned}$$

(16)

This shows that bistability has a critical influence on the propagation of excitations of the fluid, including Bogoliubov excitations.

Under optical bistability, the shape of the dispersion curve (the real part of the dispersion relation (13)) depends on the working point along the loop, i.e. on $\Delta _p$. In the special case where $\Delta _p = gn$ (point C in Fig. 1b), one recovers the behaviour given in Eq. (11) and the dispersion curve has a linear slope at low $\delta k_p$ where interactions are phononic, while the dispersion at high k is that of free massive particles, $\omega ^{\Delta _p=gn}_\pm \xrightarrow [k\rightarrow \infty ]{} \hbar k^2/2m$.³ Because of the linear, sound-like, dispersion of excitations at short k, point C is referred to as the “sonic point” of the bistability loop.

Operation at $\Delta _p < gn$ is also possible, in which case the linear behaviour at short k disappears and the dispersion becomes quadratic. As we will show in Sect. 2.5, this bears consequences on the emission and propagation of Bogoliubov excitations in the fluid. For now, we consider that $\Delta _p = gn$. In the downstream region where the fluid is ballistic, interactions at low $\delta k_0$ are phononic, as shown in Eq. (14).

2.4 Modes of the system

In the configuration of Fig. 1a, the fluid flow is transsonic: it goes from being sub- to supersonic with a sonic horizon ($v=c_B$) at $x={-2}\,\upmu {\hbox {m}}$. We plot the dispersion curve in the laboratory frame (the real part of Eq. (13)) of the subsonic fluid flow upstream in Fig. 2b/c, and of the supersonic fluid flow downstream in Fig. 2d. Blue (orange) curves correspond to $\omega _+$ ($\omega _-$) solutions of Eq. (13). In the subsonic case of Fig. 2b, the sonic behaviour close to $\omega =0$, with a slope shifted by the Doppler effect, can easily be seen. In Fig. 2d the fluid is supersonic and the shape of the dispersion curve of the excitations in the laboratory frame is changed a lot. The sonic behaviour at $\omega =0$ manifests itself by a discontinuity in the slope of the dispersion curve (14) in Fig. 2d.

Note that in the rest frame of the fluid, the $\omega _+$ ($\omega _-$) modes, which are the positive (negative) norm modes [46] have positive (negative) energies. However, in the laboratory frame, the Doppler effect modifies the shape of the dispersion relation. For subsonic fluid flows, the $\omega _-$ branch is at negative laboratory frame energies. For supersonic flows, part of the negative norm $\omega _-$ branch is pulled up to positive laboratory frame energies (up to a maximum energy which we denote by $\omega _\textrm{max}$) while part of the positive norm $\omega _+$ branch is pulled down to negative laboratory frame energies. This has a critical effect on the transmission of the excitations at the horizon.

Now that we have described the dispersive properties of the transsonic fluid, we consider the kinematics of Bogoliubov excitations therein. Because of the steady-state condition of the system, these are plane wave modes. Eq. (13) is a fourth-order polynomial, so there are four (positive laboratory-frame frequency) solutions to the equations of motion in each spatial region on either side of the interface. These solutions are found at the intersection point of line of constant $\omega$ with the dispersion branches at positive energies in Fig. 2c and d. Although these solutions share the same $\omega$ (which manifests energy conservation in the laboratory frame), they have distinct k, i.e. they are local modes of the homogeneous system. For $\omega >0$ in the upstream region there are two propagating modes of positive norm and two modes of complex $\omega$ and k, which are exponentially growing and decaying modes. For $\omega <\omega _\textrm{max}$ in the downstream region, there are four propagating modes, two of which have positive norm while the other two have negative norm. For $\omega >\omega _\textrm{max}$, there are two propagating modes of positive norm and two exponentially growing and decaying modes.

Local modes in a homogeneous region may be sorted by their respective group velocity $v_g = {\partial \omega _\pm }/\,{\partial k}$: those which have positive group velocity propagate rightwards (towards positive x) while those which have negative group velocity propagate leftwards. In the supersonic region downstream, because of the superluminal behaviour of the dispersion relation for large k, there are leftward propagating modes. This is a specificity of analogue systems based on quantum fluids (be they atomic or polaritonic [23, 40, 45]).

From the local modes, one may construct modes of the whole transsonic fluid—the global modes (GMs) [46]. These are solutions to the equation of motion that are valid in both regions on either sides of the interface. GMs correspond to waves scattering at the interface, and they describe the conversion of an incoming field to scattered fields in both regions. Although the system is driven-dissipative, GMs may be constructed similarly to those of conservative systems, i.e. as superpositions of the plane wave solutions in the two homogeneous regions on either side of the interface [40]. GMs are identified via their defining local mode depending on its group velocity $v_g$: an in GM describes the scattering of an incoming plane wave to various outgoing plane waves while an out GM describes a single outgoing plane wave resulting from the scattering of various incoming waves.

The horizon separates two regions of differing properties, so the vector bases of the in and out GMs are different. As a result, scattering at the horizon mixes the annihilation and creation operators of the field (as described by a Bogoliubov matrix [20, 47]). The matrix giving the scattering of the operators of the in-going modes a into the operators of the outgoing modes b can be expressed as [23]

$$\begin{aligned} \begin{pmatrix}b_u \\ b_{d1} \\ b^{\dag }_{d2} \end{pmatrix} = {\textbf {S}} \begin{pmatrix}a_u \\ a_{d1} \\ a^{\dag }_{d2} \end{pmatrix} \end{aligned}$$

(17)

where the indices u and d indicates the upstream and downstream positions. The ${\textbf {S}}$ matrix elements are determined from the transmission/reflection coefficients of an in-going mode into an outgoing mode. Because of their negative norm, the d2 modes must be quantized using creation operators $a^{\dag } _{d2}$ and $b^{\dag } _{d2}$ rather than annihilation operators. This is at the origin of the Hawking radiation. As a result quantum fluctuations in the in GMs $\{u_{in},d1_{in},d2_{in}\}$ will be converted into pairs of real excitations in the out GMs $\{u_\textrm{out},d1_\textrm{out},d2_\textrm{out}\}$: the Hawking radiation will occur in correlated pairs $u_\textrm{out}-d2_\textrm{out}$ (Hawking-partner) on top of the classical background formed by the mean-field of the polariton fluid.⁴

2.5 Effects of out-of-equilibrium physics

So far we have discussed the dispersive properties of the fluid when $\Delta _p=gn$, that is, when operating at the sonic point C of the bistability loop. In the regime $\Delta _p<gn$, the microcavity acts as an optical limiter [37]: as can be seen in Fig. 1b, the growth of the speed of excitations $c_B$ with the pump strength on the upper branch of the loop is sub-linear. While the fluid is stable in this regime, a gap opens between the $\omega _\pm$ branches of the dispersion curve, see Fig. 2c. We mark the bottom of the $\omega _+$ curve as $\omega _\textrm{min}$.

This behaviour is markedly different from that observed in systems like quantum fluids of atoms. There, the oscillation frequency of the condensate wavefunction corresponds to the chemical potential. Instead here it corresponds to $\omega _p$. The opening of the gap illustrates how tuning $\Delta _p$ gives access to a unique phenomenology of collective dynamics. Specifically, the linear behaviour at short k disappears as soon as the gap opens and the dispersion is always quadratic, meaning that the polariton ensemble cannot be superfluid. As we will see in Sect. 3.1, this departure from superfluid propagation modifies the density of the fluid in the region $x<0$ as a function of the pump strength and profile.

In Sect. 2.4, we have seen that the Hawking Effect consists in the mixing of in GMs of opposite sign of norm at the horizon, $u_{in}$ from $x<0$ and $d1_{in}$ and $d2_{in}$ from $x>0$. The downstream modes only exist over the limited interval $0<\omega <\omega _\textrm{max}$. When the gap between $\omega _-$ and $\omega _+$ opens, $u_{in}$ only exists for $\omega>\omega _\textrm{min}>0$ so the frequency interval for scattering is reduced to $\omega _\textrm{min}<\omega <\omega _\textrm{max}$. We will show in the simulations (cf Sect. 3.2) that the strength of vacuum emission by the HE is thus decreased.

In brief, when $\Delta _p<gn$, out-of-equilibrium physics manifests itself in the opening of a gap between the branches of the dispersion relation and a modification of the shape of the dispersion to a purely quadratic form. This affects the generation of Bogoliubov excitations as well as their propagation in the fluid.

3 Emission by the Hawking effect

We now perform calculations with the cavity parameters of [34]: $\hbar \gamma =$ 0.047 meV, $\hbar g =$ ${0.0003}\,{{\hbox {meV}}\, \upmu {\hbox {m}}}$ (same as in [32]), $m^* = 3 \cdot 10^{-5} m_e$. Importantly, $\omega _p-\omega _0 =$ 0.49 meV was kept constant throughout.

We study vacuum emission via non-local correlations in the fluid density [21], which we quantify with the normalised spatial correlation function

$$\begin{aligned} g^{(2)}(x,x') = \frac{G^{(2)}(x,x')}{<n(x)> <n(x')>}. \end{aligned}$$

(18)

$G^{(2)}(x,x')$ is the two-point correlation function of the field (cf “Appendix B”).

In Fig. 3 we plot the operation point on the bistability loop of the fluid on either side of the horizon (solid line, upstream, dashed line, downstream), the pump profile (black line) and ensuing properties of the inhomogeneous fluid—characterised by its velocity (red line) and the speed of excitations $c_B$ (blue line)—as well as the resulting density-density correlations (18). We are interested in the fluid properties and their influence on correlated emission. Note that we plot $c_B$ for better comparison between the working point along the bistability loop and the local properties of the inhomogeneous fluid.

Fig. 3

Correlated emission as a function of the trans-sonic flow profile. Left column, bistability loop: solid line, upstream region; dashed line, downstream region. Dots, corresponding pump strength. Middle column, pump strength $|F_p|\, [\textrm{ps}^{-1}.\upmu \textrm{m}^{-1/2}]$: solid black, upstream; dashed black, downstream. Velocities [$\upmu {\hbox {m}} \, {\hbox {ps}}^{-1}$]. blue, speed of excitations $c_B$; orange, fluid flow velocity v set by $k_{p,u}$. If the pump strength is zero in the downstream region (c, d and e), the fluid propagates ballistically (free evolution with no support) there. Right column, spatial correlation function $g^{(2)}(x,x')-1$ (Eq. (18)), colour scale from $-1.25\times 10^{-4}$ to $1.25\times 10^{-4}$

3.1 Fluid configurations

We consider various flow profiles on either side of the horizon. The bistability of the fluid on either side of the horizon may be tuned by controlling the wave-number of the fluid in either region by means of the pump ($k_{p,u}$ or $k_{p,d}$ in the up- or downstream region, respectively), see “Appendix A”. The fluid density may be supported on the higher branch of the bistability loop by means of the pump intensity. In order to explore the full parameter space we have computed 36 correlation spectra. Not all combinations are interesting, though, so in Fig. 3 we present five configurations that give typical behaviours. We consider two types of settings: in rows (a) and (b), there is an additional beam that supports the polariton density in the downstream region, while in rows (c), (d) and (e) the pump strength is set abruptly to zero at $x=-{7}\,\upmu {\hbox {m}}$ and so the fluid is left free to evolve from that point on (across the defect into the downstream region). Row (a), the fluid density is set near (but not at) the sonic point in both regions with a jump in $k_p$ at the interface ($k_{p,u} = {0.25}\,\upmu {\hbox {m}}^{-1}$, $k_{p,d} = {0.55}\,\upmu {\hbox {m}}^{-1}$), the pump strength is set to 0 at $x=-{7}\,\upmu {\hbox {m}}$ the horizon; row (b) the fluid density is set near to and away from the sonic point on the upper branch of the bistability loop in the up- and downstream region, respectively, with a jump in $k_p$ at the interface ($k_{p,u} = {0.25}\,\upmu {\hbox {m}}^{-1}$, $k_{p,d} = {0.58}\,\upmu {\hbox {m}}^{-1}$); rows (c), (d) and (e), the fluid density is set gradually closer to the sonic point on the upper branch of the bistability loop in the upstream region ($k_{p,u} = {0.25}\,\upmu {\hbox {m}}^{-1}$),. In configuration (e) the fluid density is supported as close as possible to the sonic point.

In all configurations, the fluid builds up in the region $-{10}\,\upmu {\hbox {m}}<x<x_d$: a relatively small amplitude bump in the density forms before the defect (the fluid velocity slightly dips in that region). On the other hand, while the density of the fluid is mostly flat in the configuration of Fig. 3a, in that of Fig. 3b its amplitude undulates widely over ${100}\,\upmu {\hbox {m}}$ downstream of the horizon before flattening down. This illustrates how attempting to force the fluid properties to a working point away from the sonic point after it has propagated across an obstacle destabilises it. Meanwhile, the fluid density decreases exponentially (with a high, ballistic, wavenumber) in configurations (c), (d) and (e) where there is no pump in the downstream region.

Given the variety of fluid properties and the possible fast variations within, the description of the system as two homogeneous media adopted in Sect. 2 and amenable to analytical solutions is not valid everywhere. Instead we must calculate the flow profile and quantum fluctuations at all points. To this end, we use the Truncated Wigner Approximation (see “Appendix B”) to evolve the wave function and obtain the properties of the fluid at all points in the cavity as well as the dynamics of the Bogoliubov excitations therein. This numerical method is adapted to analogue systems based on atomic as well as polaritonic quantum fluids [21, 43]. Here, it enables the study of vacuum emission on highly varying backgrounds. All maps result from a statistical average over $10^6$ Monte-Carlo realisations.

3.2 Correlation diagrams

In all configurations, correlations may be sorted by the spatial region in which the involved modes propagate, which correspond to four quadrants in the plots. The South West quadrant ($x<0$, $x'<0$) corresponds to correlations in the upstream region; the South East and North West quadrants correspond to correlations across the horizon in the up- and downstream regions; the North East quadrant corresponds to correlations in the downstream region. All configurations have some common traces, which are most visible in Fig. 3e: (i) anti-correlations along the $x=x'$ diagonal that indicate anti-bunching under repulsive polariton interactions; (ii) a negative moustache-shaped trace in the upstream-downstream region that indicates correlations across the horizon between Hawking radiation and modes in the downstream region; (iii) oblique interference fringes localised along the $x=0$, $x'>0$ half line (and, symmetrically, the $x'=0$, $x>0$ half line), and another series (iv) of fringes localised along the $x=0$, $x'<0$ half line. While traces (i) and (ii) are generic features of the HE in quantum fluids, see, e.g. [21, 23, 30, 32, 33], the fringes (iii) and (iv) are new. They indicate correlations between the propagating modes $u_\textrm{out}$ and $d1_\textrm{out}$ and a mode bound to the horizon. This coupling manifests the excitation of quasi-normal modes of the acoustic field under vacuum-driven perturbations [48]. (Note that this was not seen in [32], where a repulsive defect was used, see “Appendix D”)

Configuration 3e leads to a longer Hawking moustache (ii) than obtained in other models for a quantum fluid: about ${35}\,\upmu {\hbox {m}}$- and ${105}\,\upmu {\hbox {m}}$-long in the up- and downstream regions, respectively. It also features stronger correlations than obtained in previous polariton analogue studies, with a maximal relative amplitude of $1.5\cdot 10^{-4}$. In comparison, the Hawking moustache is shorter and weaker in all other configurations: we observe the influence of the regime of density of the fluid properties (working point on the bistability loop on either side of the horizon) on vacuum emission at the horizon and propagation in either region thereafter. We remark that supporting the density of the fluid in the downstream region as in configuration (a) does not aid emission while coming at a higher technical cost (due to the precise matching of the pump wavenumber across the defect). Finally, although emission occurs in all configurations, as a comparison between configurations (c), (d) and (e) shows, its strength is limited by the opening of the gap in the dispersion upstream, which reduces the interval $\{\omega _\textrm{min},\omega _\textrm{max}\}$ over which the HE occurs. This is why the emission strength in [32] was lower than in Fig. 3e although a similar pump profile was used.

In brief, we have established that operating such that the fluid is as close as possible to the turning point C of the bistability loop upstream of the horizon and letting the fluid propagate ballistically downstream enhances the emission and propagation of Bogoliubov excitations. In that regard, operating with a flat pump profile whose spatial extension is well controlled is better than with a Gaussian profile.

4 Discussion

We showed how engineering the density of a quantum fluid of polaritons can enhance the emission and propagation of paired Bogoliubov excitations in a transsonic flow. Our work sheds light on the interplay between optical bistability and parametric amplification in fluids of light. The bistable behaviour of a system can thus be exploited to study field theoretic effects like the HE in the laboratory. Specifically, we have found that fine control over the fluid properties may be achieved with a step-like pump profile.

Here we observed the generation and propagation of paired Bogoliubov excitations of the quantum fluid on either side of a sonic horizon when supporting the density of the fluid at various points in the bistable regime. Support of an inhomogeneous fluid density and velocity may be achieved by changing the wave-number of the pump. In an experiment, this is easily implemented with high spatial resolution (limited by diffraction) thanks to spatial light modulators [49,50,51,52,53]. We found that letting the fluid flow ballistically across an attractive defect so as to form a horizon yields Hawking correlations of the order of $2\times 10^{-4}$ over more than ${100}\,\upmu {\hbox {m}}$. Thus we obtained a strong increase of the quantum correlations by the HE: these are a fourfold enhancement in the strength and tenfold enhancement in the propagation length of correlations compared with previous results in quantum fluids of light. Furthermore, we showed how moving away from this optimal configuration reduces the signal strength and length—two effects directly linked to the kinematics of Bogoliubov excitations in the out-of-equilibrium fluid. In this way, our work demonstrates that the correlation traces are a diagnostic for the influence of out-of-equilibrium physics on mode conversion in inhomogeneous flows. For example, in the paper [48], we have explained how a dissipative quench of a mode bound to the horizon yields novel, local correlation traces between a local and propagating modes, i.e. a quasi-normal mode of the acoustic field.

Given the strength of the Hawking correlations reported in this paper, these should be amenable to experimental detection with state-of-the-art apparatus. As such, our methods open the way to the theoretical and experimental study of the quantum statistics of the HE in driven-dissipative systems: for example, one could calculate (and observe) the Hawking correlations in reciprocal space [29], thus gaining frequency-resolved information on them [30], which could in turn be used to measure entanglement in the Hawking emission [19, 20].

Appendices

Appendix A The physical system

In this appendix we present the field theory of polaritons and their fundamental excitations.

Our system is a one-dimensional quantum fluid of exciton-polaritons whose flow velocity goes from being sub- to super-sonic, thus forming a sonic horizon where the local flow velocity of the fluid equals the local speed of sound. For simplicity, we may consider that the horizon separates two spatial regions whose properties are independent of space—two homogenous regions, although we shall eventually depart from this simplified picture. For now, we begin with the theoretical description of a homogeneous quantum fluid and of the propagation of quantum fluctuations of its phase and density.

In the majority of cases the interaction energy does not match the effective detuning and the dispersion curve is thus modified. Furthermore, writing the density of the fluid as a function of the intensity of the laser yields several solutions [43]. This degeneracy of fluid densities is due to optical bistability [39], which, as we will show, has tremendous influence on the emission and propagation of excitations of the fluid, including Bogoliubov excitations. Here we describe the influence of the bistability on the Bogoliubov dispersion.

We begin by describing the relationship between the density of polaritons, n, and the intensity of the pump laser, $|F_p|^2$ in the case where the energy of the laser is above that of the lower polaritons, $\Delta _p > \gamma \sqrt{3}/2$: we square Eq. (7) and find

$$\begin{aligned} \left( (gn - \Delta _p)^2 + \frac{\gamma ^2}{4}\right) n = |F_p|^2 \end{aligned}$$

(A1)

or, equivalently,

$$\begin{aligned} \left( \left( \frac{m^* c_B^2}{\hbar } - \Delta _p\right) ^2 + \frac{\gamma ^2}{4}\right) \frac{m^* c_B^2}{g\hbar } = |F_p|^2. \end{aligned}$$

(A2)

The physics at play may be investigated equivalently in terms of the relationship between the speed of excitations $c_B$ and the strength of the pump, as shown in Fig. 4. At first, $c_B$ increases slowly with $|F_p|$ (arrow (1)), until $|F_p|=F_1$ where it increases abruptly (arrow (2)). For $|F_p|>F_1$, $c_B$ increases slowly again. If the pump’s strength is decreased from $|F_p|\ge F_1$, $c_B$ decreases slowly until $|F_p|=F_2$ (arrow (3)), where it falls abruptly (arrow (4)). Since $F_1>F_2$, the $c_B$ to $|F_p|$ relationship presents a hysteresis cycle with two regimes of speed of sound: the low density regime when $|F_p|<F_1$ and $c_B$ is low, and the high density regime when $|F_p|>F_2$ and $c_B$ is high. This hysteresis cycle is the manifestation of optical bistability [39], so we will henceforth refer to it as the ‘bistability loop’. Note that the dashed line in Fig. 4 is unstable and the speed of sound will actually follow the hysteresis cycle schematised by arrows (1)–(4).

Fig. 4

Bistability loop for an homogeneous polaritonic fluid. $\omega _p - \omega _0>0$ and $k_p = 0$. Black, stable points; dashed, unstable points. The system is bistable for $F_2<|F_p|<F_1$ and follows the hysteresis cycle (1)–(4)

Now, in order to explicitly show the dependence of the Bogoliubov dispersion on the density of the fluid as well as the influence of optical bistability thereon, we generalise Eq. (11): We diagonalise the Bogoliubov matrix $\mathcal {L}$ for a homogeneous system pumped with arbitrary strength and obtain

$$\begin{aligned} \omega _{\pm } (k) = \pm \sqrt{ \left( \frac{\hbar k^2}{2 m^*} + 2gn - \Delta _p\right) ^2 - (gn)^2} - i \gamma /2.\nonumber \\ \end{aligned}$$

(A3)

Eq. (13) is the Doppler shifted version of Eq. (A3). In Fig. 5, we show the dispersion curve for 5 different fluid densities along the bistability loop. As can be seen in Fig. 5a and b, the shape of the dispersion does not change much in the linear regime: the two branches of the dispersion curve cross. When the fluid is bistable ($F_1\le |F_p|\le F_2$), in Fig. 5b, we observe the appearance of plateaus characteristic of an unstable fluid at the crossing points. On the other hand, the shape of the dispersion curve changes significantly in the high density regime depending on the position along the bistability loop: at high pump strength (Fig. 5e), the two branches are split in energy by a gap that increases with the pump strength. The sonic dispersion relation (11) is recovered at point C (Fig. 5c), while for slightly lower pump strength (Fig. 5d), the plateau at low k is characteristic of an unstable fluid (similarly to Fig. 5d) . Note that the dispersion curve has a linear slope at low k (and thus a sonic interpretation) at point C only, which is thus sometimes referred to as the ‘sonic point’ of the bistability.

As Eq. (A3) is of order four in k, the dispersion has four complex roots. The real part of these roots is non-zero in the low density regime (Fig. 5a, b) as well as at points C and $C'$ (Fig. 5c, d), but not at point D (Fig. 5e).

In this appendix, we have seen that the mean-field of a polariton system behaves as a fluid. We have reviewed the dispersion relation of Bogoliubov excitations in this fluid and seen that optical bistability of the fluid strongly influences the properties of this dispersion relation. These considerations may be generalised to a fluid whose density is not homogeneous.

Fig. 5

Bogoliubov dispersion for various fluid densities. Top row: Eq. (A3) is plotted in the fluid frame. a and b regime of low density; c unstable fluid; d, sonic dispersion; e regime of high density. Bottom, f bistability curve for a homogeneous fluid. A low density; B low density and bistable; C’ unstable; C sonic point; D high density and bistable

Fig. 6

Emission with a repulsive defect. a Flow profile. b Correlated emission

Appendix B Numerical method and correlation function

In this appendix we present the numerical method used to compute the correlation maps of the main text.

Our interest is in vacuum emission, that is amplification of the quantum vacuum fluctuations at the horizon (the spontaneous Hawking effect). The quantum description of the Bogoliubov excitations relies on the dispersion relation of the classical field, Eq. (13). In order to encompass quantum effects, we use a quantum Monte-Carlo method called the truncated Wigner approximation (TWA). In this method, the equation of motion is truncated so as to map it to a stochastic partial differential equation for a classical field $\Psi$:

$$\begin{aligned} i{\hbox {d}}\Psi= & {} \left( \omega _0 - \frac{\hbar }{2m^*}\frac{{\hbox {d}}^2}{{{\hbox {d}}x}^{2}} + V + g(|\Psi |^2 - 1/\Delta x) - i\frac{\gamma }{2}\right) \nonumber \\{} & {} \times \Psi {\hbox {d}}t+ F_p {\hbox {d}}t + \sqrt{\frac{\gamma }{4 \Delta x}} {\hbox {d}}W, \end{aligned}$$

(B4)

where ${\hbox {d}}W$ is complex white noise. In numerical simulations, sampling of the realisations obtained with (B4) starts when the steady state is reached. One must ensure that enough time is spent between each sampling to ensure independence of the realisations. Quantum observables are computed with statistical averaging over the realisations obtained with the TWA: the general rule for N arbitrary observables is [43]

$$\begin{aligned} \langle O_1...O_N\rangle _W = \frac{1}{N!} \sum _{\tiny {\begin{array}{c} \text {All} \\ N-\text {permutations} \end{array}}} \langle P({\hat{O}}_1,...,{\hat{O}}_N)\rangle ,\qquad \end{aligned}$$

(B5)

where $\langle \rangle _W$ denotes the statistical averaging over the realisations. In all figures of this paper we computed 1 million realisations.

Emission at the horizon by the HE is best detected via nonlocal correlations in the fluid density [21]. These may be quantified via the normalised spatial correlation function

$$\begin{aligned} g^{(2)}(x,x') = \frac{G^{(2)}(x,x')}{G^{(1)}(x)G^{(1)}(x')}. \end{aligned}$$

(B6)

$G^{(2)}(x,x')$ is the diagonal two-points correlation function of the field, which is calculated from (B5) and normally ordered using Bose statistics

$$\begin{aligned} G^{(2)}(x,x')= & {} \langle \hat{\Psi }^\dagger (x) \hat{\Psi }^\dagger (x') \hat{\Psi }(x') \hat{\Psi }(x)\rangle \nonumber \\= & {} \langle \Psi ^* (x) \Psi ^* (x') \Psi (x') \Psi (x) \rangle _W -\frac{1}{2\Delta x} (1+\delta _{x,x'}) \nonumber \\{} & {} \times \bigg (\langle \Psi ^* (x) \Psi (x) \rangle _W +\langle \Psi ^* (x') \Psi (x') \rangle _W -\frac{1}{2 \Delta x} \bigg ),\nonumber \\ \end{aligned}$$

(B7)

while the diagonal one-point correlation function is

$$\begin{aligned} G^{(1)}(x) = \langle \hat{\Psi }^\dagger (x) \hat{\Psi }(x)\rangle = \langle \Psi ^*(x)\Psi (x)\rangle _W - \frac{1}{2 \Delta x}.\qquad \end{aligned}$$

(B8)

Appendix C Constraints on the calculations

All configurations in Fig. 3 have been realised with the cavity parameters of [34]. When exploring all possible configurations of fluid density on either side of the horizon, some constraints must be abode by. We present them in this appendix.

The first constraint is on the upstream pump wavevector $k_{p,u}$ for a fluid near the sonic point. The fluid is at the sonic point for

$$\begin{aligned} c_u = \sqrt{\frac{\omega _p-\omega _0-\hbar k_{p,u}^2/2m^*}{m^*}}, \end{aligned}$$

(C9)

together with the upstream condition $v_u<c_u$, this yields an upper bound for the upstream fluid flow velocity and thus for the wavevector of the pump:

$$\begin{aligned} \omega _p-\omega _0>\frac{3}{2}m^* v_u^2. \end{aligned}$$

(C10)

For the value of detuning used throughout this paper the upper bound is around $k_{p,u} = {0.28}\,\upmu {\hbox {m}}^{-1}$. In most simulations, we used $k_{p,u} = {0.25}\,\upmu {\hbox {m}}^{-1}$ in order to be close to the bound while leaving a small interval for easier simulations.

Exploring all regimes of density in the downstream region comes with some constraints as well: For instance, placing the fluid in the upper part of the bistable regime as in configurations 3b is easier for a large bistable interval $F_1-F_2$. The turning point in configuration 4f is obtained at

$$\begin{aligned} c_d = \sqrt{\frac{\omega _p - \omega _0 - \hbar k_{p,d}^2/2m^*}{2m^*}}. \end{aligned}$$

(C11)

and the width of the interval is then given by

$$\begin{aligned} \Delta I^{bistable}= & {} |F_{p,\textrm{max}}|^2 - |F_{p,\textrm{min}}|^2 \nonumber \\= & {} \left( \frac{4}{9}\left( \omega _p-\omega _0-\frac{\hbar k_{p,d}^2}{2m^*}\right) ^2 - \frac{1}{2}\gamma ^2\right) \nonumber \\{} & {} \times \frac{\omega _p-\omega _0-\hbar k_{p,d}^2/2m^*}{3g}. \end{aligned}$$

(C12)

A bistable regime exists only if $\omega _p-\omega _0>\hbar k_{p,d}^2/2m^*+\gamma \nicefrac {\sqrt{3}}{2}$, hence an upper bound on $k_{p,d}$.

Furthermore, the speed of excitations right after the defect, $c_{d}$ is fixed by the upstream parameters and the strength of the defect, $V_\textrm{def}$. Pumping in the upper branch of the bistability requires

$$\begin{aligned} c_d > \sqrt{\frac{\omega _p-\omega _0 - \hbar k_{p,d}^2/2m^*}{m^*}}. \end{aligned}$$

(C13)

This critical point needs to be below $c_{d}$ for the fluid density to be on the upper branch. (Note that it would also be possible to change the value of $c_{d}$, which can be achieved for a weaker defect potential—energy conservation before and just after the defect links $V_\textrm{def}$ and $c_{d}$). The choice of different $k_{p,d}$ in the simulations of Fig. 3 is a consequence of all these constraints.

Appendix D Repulsive defect

Here we consider the configuration of Fig. 3e but with a repulsive defect ($V_\textrm{ext}=0.85$ meV) instead.

This configuration is comparable to that of ref [32] but with a pumping scheme such that the fluid density is supported as close to point C of the bistability as possible. In Fig. 6, we observe that, although stronger and longer than in [32], correlated emission by the Hawking effect is weaker than in Fig. 3e. So the attractive potential of the defect aids correlated emission. Furthermore, the horizontal/vertical traces signalling correlations with the horizon region (vacuum quantum excitation of a QNM of the acoustic field [48]) are absent.