Formation dynamics of black- and white-hole horizons in an analogue gravity model

We investigate the formation dynamics of sonic horizons in a Bose gas confined in a (quasi) one-dimensional trap. This system is one of the most promising realizations of the analogue gravity paradigm and has already been successfully studied experimentally. Taking advantage of the exact solution of the one-dimensional, hard-core, Bose model (Tonks-Girardeau gas) we show that, by switching on a step potential, either a sonic (black-hole-like) horizon or a black/white hole pair may form, according to the initial velocity of the fluid. Our simulations never suggest the formation of an isolated white-hole horizon, although a stable stationary solution of the dynamical equations with those properties is analytically found. Moreover, we show that the semiclassical dynamics, based on the Gross-Pitaevskii equation, conforms to the exact solution only in the case of fully subsonic flows while a stationary solution exhibiting a supersonic transition is never reached dynamically.


Introduction
Analogue models of gravity are useful tools bridging gravitation to other branches of physics and they have been intensively investigated since their proposal [1] in order to study effects whose experimental achievement is hardly possible in the cosmological context. In the past years a great scientific effort was dedicated to the observation and characterization of the Hawking effect [2,3] but several different phenomena can be investigated within the analogue gravity paradigm [4]. In particular, the physics of sonic horizons [5] has been vastly studied, both from a theoretical and an experimental point of view and analogue black-hole and white-hole horizons have been experimentally achieved in optics [6][7][8], in Bose-Einstein condensates (BECs) [9][10][11][12][13], in water [14][15][16][17][18], in superfluid Helium [19] and a in few other setups. In the context of BECs, a few recent experiments [10][11][12][13] successfully recreated sonic horizons in order to probe the elusive Hawking radiation but the results are still under debate in the scientific community [20][21][22][23][24][25][26][27][28]. Nonetheless, these experiments represent a remarkable achievement and they have stimulated a variety of theoretical studies on sonic black holes in 1D BECs. One feature, though, which has not been deeply investigated in these kind of configurations is the formation dynamics of an analogue black/white hole; this may be due to the fact that analogue models can hope to reproduce only the kinematical aspects of the gravitational phenomena, as the analogy breaks down for the dynamical properties. For this reason, often theoretical studies start from an "eternal" (analogue) black-hole, i.e., from a configuration which already shows a sonic horizon, neglecting the question of how it may have formed.
In this work we will study the dynamics of formation of a sonic horizon for a simple, realistic model of BEC which mimics existing experimental configurations [9][10][11][12][13]. In particular, we will exploit a recently-proposed exact model [26,27] in order to study the dynamics of the analogue of a gravitational collapse in the absence of external confining traps and to investigate whether the system only forms black-hole-like horizons or if a dynamical process could, in principle, lead to a formation of a white-hole horizon too. Furthermore, we will compare the results with those obtained with a semiclassical approach -commonly adopted in the theoretical interpretations of the experiments -in order to highlight possible inadequacies of the semiclassical paradigm. This will give insights on the formation process of analogue horizons and, interestingly, will also lead to challenge the "eternal" black-hole configuration as a tool to study the Hawking effect.

Tonks-Girardeau gas
We consider a one-dimensional Bose fluid of hard-core, point-like particles in an external potential V (x). The model hamiltonian is written, in first quantization, as where the hard-core constraint corresponds to the strong interaction limit c → ∞. In the absence of external potential V (x) = 0, eigenvalues and eigenfunctions of this hamiltonian have been found for any value of c by Lieb and Liniger [29] using Bethe Ansatz techniques. In the hard-core limit, the Lieb-Liniger solution can be more easily obtained by performing a Jordan-Wigner transformation which maps the hard core Bose fluid into a non-interacting, spinless, fermion model [30]. A simple argument allows to understand this exact result: in one dimension the spatial ordering of hard-core bosons cannot be modified by a local hamiltonian. Therefore, the different phases acquired by bosons and fermions upon particle exchange do not come into play and the matrix element of the hamiltonian between two arbitrary bosonic configurations coincides with the corresponding matrix element of the free spinless Fermi gas. Interestingly, this result remains valid also in the presence of an arbitrary external potential, but clearly the argument fails in more than one dimension, where a natural ordering cannot be defined. Due to this mapping, the full spectrum of the hamiltonian (1) and the many-particles dynamical density correlations coincide with those of a non-interacting Fermi fluid in the same external potential V (x). Note that, instead, the momentum distribution of the Bose gas differs from that of the corresponding Fermi system due to the non-locality of the Bose-Fermi mapping defined by the Jordan-Wigner transformation.
Exactly solvable models represent a unique tool as a test ground for approximate theories, which are the only possible option in most cases. Among the approximate descriptions of a weakly interacting Bose fluid, the celebrated Gross-Pitaevskii equation (GPE) [31] stands out as a simple and accurate method to portray the dynamics of the condensate wavefunction. The inclusion of a cubic non-linearity in the Schrödinger equation suffices to account for weak interparticle interactions, providing a satisfactory description of many experiments in cold atom physics [31]. The Tonks-Girardeau (TG) gas is a strongly interacting Bose fluid and therefore the usual semiclassical approaches for the approximate description of the dynamics are not expected to provide accurate results. However, it is known that a simple generalization of the Gross-Pitaevskii equation does indeed faithfully reproduce the physics of a Tonks-Girardeau gas [32,33], at least in homogeneous systems. Such a "strong coupling limit" of the usual (cubic) GPE displays a quintic non-linearity with a specific value of the coupling constant g = 2 π 2 2m : In the Gross-Pitaevskii approximation, the global properties of the flowing Bose gas are identified with those derived from the single particle wavefunction ψ(x, t). In particular, in the homogeneous case (V (x) = 0), the phonon excitation spectrum derived from Eq. (2) coincides with the exact one, defined by a sound velocity c, related to the fluid density n by the same expression valid for a one dimensional Fermi gas: c = π n m . The dynamics described by the Gross-Pitaevskii equation will be compared with the exact results in order to verify the accuracy of the widely adopted semiclassical approximation in the context of analogue gravity experiments in BECs.

The model
The Tonks-Girardeau gas model has been recently studied in the framework of analogue gravity [26,27] because it allows to describe, from a microscopic point of view, the dynamics of formation of a sonic horizon and to investigate the physical conditions leading to the expected Hawking-like phonon flux emerging from the horizon. The analytical and numerical studies of Ref. [26,27] proved that a precise correspondence between the physics of the analogue model and that of the Hawking effect can be obtained only when the TG gas flows against an extremely smooth potential barrier. The rationale for this result is that the identification of the quasiparticles (phonons) as a free quantum field is possible only for very low excitation energies and (almost) homogeneous states. Instead, sharp potentials, like those usually adopted in the experiments, give rise to a steady state characterized by a rapidly varying density profile, whose elementary excitations cannot be expressed as a free quantum field with a locally varying sound velocity, as requested by the analogy to the gravitational problem. Remarkably, a sharp external potential does not seem to suppress the Hawking effect: in the steady state a phonon flux may indeed be present but it is not described by a thermal distribution and therefore it cannot be associated to a well-defined Hawking temperature. We stress that these results have been obtained [26,27] for the Tonks-Girardeau gas in one dimension, where the occurrence of thermalization after a quantum quench is hampered by the integrability of the model [34]. On the other hand, at variance with the case of a potential barrier, an extremely smooth potential step does not lead to the formation of a sonic horizon, making such a configuration unsuitable for the investigation of the Hawking effect [26,27].
In this work we will examine the formation dynamics of black-and white-hole horizons when a flowing TG gas is perturbed by the sudden switching-on of a sharp step potential. In the experiments, this procedure has been adopted in quasi one-dimensional BECs in order to create either a sonic black-hole horizon [11][12][13] or a black/white hole pair in order to trigger the so-called "black-hole laser" effect [10]. According to the experimental configurations and to our previous theoretical studies [26,27], we choose a left-moving gas flowing over a "waterfall" potential as a preferred setup but we will generalize our investigation to also describe the opposite flow.
Thus, we start from a flowing TG gas in a stationary state which is described, via the Bose-Fermi mapping, by a Slater determinant of plane waves 1 √ 2π e ipx ; the wavevectors p belong to the interval −k F − k 0 ≤ p ≤ k F − k 0 , where k F > 0 is the Fermi momentum in the co-moving reference frame and k 0 is the momentum shift due to the fluid flow; these two parameters, characterizing the initial state of the system, are related to the uniform particle density n = k F π , sound speed c = k F m and fluid velocity v = − k0 m . We therefore adopt the convention that a positive value of k 0 corresponds to a left-moving flow. Then, at t = 0, we suddenly turn on the step potential where V 0 is the potential height which is conveniently parametrized through the wavevector Q: V 0 = 2 Q 2 2m . Immediately after the quench the fluid is not in a stationary state any more and its many particle wavefunction evolves in time. Adopting the fermionic representation, at any time t > 0 the wavefunction is given by the Slater determinant of the time evolution of the initial plane waves: where φ k (x) are the exact single particle eigenfunctions in the presence of the external potential, k = ω k the corresponding energy, and η → 0 + is the usual convergence factor. For a step potential k = 2 k 2 2m and the eigenfunctions φ k (x) are easily obtained by elementary methods: • for k > Q The explicit expressions of the reflection and transmission coefficients are reported in [27].

Steady states
While the full dynamics after the quench requires a numerical analysis starting from Eqs. (4)(5), the long time behavior can be analytically evaluated. As shown in Ref. [27], at any given fixed position x, the local physical properties of the model approach those given by the stationary state built upon the wavefunctions ψ p (x, t) given by: As a consequence, for each choice of the two parameters (k F , k 0 ) a steady state exists, characterized by the number density profile (10) n and by the mass current which is constant in space and time. By direct substitution of the scattering eigenfunctions in Eq. (9) we get the explicit results: for x < 0 and for x > 0, where q = p 2 − Q 2 and k = p 2 + Q 2 . The occupation number f (p) identifies the region in momentum space occupied by fermions: f (p) = 1 for −k F − k 0 ≤ p ≤ k F − k 0 and f (p) = 0 elsewhere. These expressions allow to evaluate several important quantities in the framework of the analogue gravity: the local fluid and sound velocity, v(x) and c(x) respectively, together with their ratio, i.e. the Mach number: The analytical expressions show that in the limiting cases k 0 ≥ k F the density profile is rigorously flat in the region x < 0, leading to a uniform fluid velocity. The same can be said for the region x > 0 when −k 0 ≥ Q + k F . At first, tough, we will consider only the cases of |k 0 | ≤ k F as the condition ensures that the flow is initially subsonic.
The integrals can be easily evaluated for both the particle current and the asymptotic densities at x → ±∞ and allow to evaluate the range of parameters leading to a steady state characterized by a single supersonic transition. Figure 1 summarizes the results: the blue curve refers to a left moving fluid (k 0 > 0) with subsonic velocity at x → +∞ which gets accelerated by "falling into the waterfall" and becomes supersonic at x → −∞ if k F − k 0 lies below the blue line. Instead, the red curve shows that a right moving fluid (k 0 < 0), which starts subsonic at x → −∞, can get supersonic flowing against the step potential if k F + k 0 belongs to the region below the red line. In this case, however, k F must exceed the minimum value Q 2 . In both cases, a black-hole horizon forms close to x = 0, only for sufficiently high initial velocities: |k 0 | k F . The steady states corresponding to an initially still Fermi gas, i.e. a vanishing k 0 , display a finite particle current due to the sudden appearance of the waterfall potential, leading to a leftwards motion (i.e. with negative velocity) of the fluid. For this case, two velocity profiles are shown in Figure 2 for k F = π 2 Q and k F = π 4 Q: while in the former case the steady state is always subsonic, in the latter a supersonic flow sets in a finite region around x ∼ 0, giving rise to a pair of black-hole/white-hole horizons.
When we start from a fluid which is already moving (i.e. k 0 = 0), the long time characteristics of the flow velocity and of the sound speed show different possible behaviors. Illustrative plots of the stationary velocity profiles for a few representative choices of k 0 at fixed k F = π 2 Q are shown in Figure 3 for both cases of a fluid falling into the waterfall (k 0 > 0, left panels) and flowing against the step potential (k 0 < 0, right panels). While for k 0 = 0 Figure 2 shows that the fluid motion is subsonic, when k 0 is increased (i.e. when the fluid falls from the potential step) a small supersonic region appears close to the step (see the case k 0 = 3π 10 Q). Further increasing k 0 (see k 0 = 9π 20 Q) the flows becomes supersonic in the whole downstream (i.e. left) region. At k 0 = k F = π 2 Q the velocity becomes constant for all x < 0, as previously discussed. The case of a fluid flowing against the step potential (right panels) is qualitatively similar, showing a fully subsonic flow for k 0 = − 3π 10 Q, a black-hole/white-hole pair at k 0 = − 9π 20 Q and a single supersonic transition at k 0 = −k F = − π 2 Q. In the latter case, however, the fluid velocity is not constant in the downstream (i.e. right) region. In all the cases, the undulations present both in the upstream and downstream regions are due to the presence of a sharp Fermi surface.
Note that the steady state solutions discussed here are not invariant by time reversal. Being the hamiltonian real, the complex conjugate of the single particle wavefunctions (9) are still exact eigenfunctions for the free Fermi gas; the time reversal of a scattering state, however, is not a scattering state but is written as a particular linear combination of the two scattering states. Therefore, the Slater determinant of the time-reversed states describes a different stationary state of the system with the same density profile and opposite current j, as shown by Eqs. (10,11): wherever the original solution displays a black-hole (BH) horizon, the time-reversed state predicts a white-hole (WH) horizon. However, we stress that, according to the analysis leading to Eq. (9), the time reversed state does not describe the stationary solution spontaneously reached by the system after a quench.
It is interesting to compare these exact results with the stationary states of the corresponding Gross-Pitaevskii equation, to check whether the semiclassical approach is able to capture the variety of behaviors displayed by the TG gas model. Looking for solutions to Eq. (2) of the form ψ(x, t) = Φ(x) e − i µt we obtain a non-linear differential equation whose solutions can be found analytically in the case of a step potential V (x) = V 0 Θ(x), where Θ(x) is the Heaviside function. The generic stationary solution depends on two parameters: the chemical potential µ and the current j. The equation has real coefficients and then each stationary solution has a time reversal partner with opposite current: the direction of the particle flow is reversed by keeping the same density profile. Restricting our analysis to solutions leading to density profiles asymptotically flat 1 for x → ±∞, we generally find only monotonic decreasing density profiles and fully subsonic flows. Only for a special relationship between µ and |j| a different solution appears, displaying a supersonic transition, which might represent either a BH or a WH horizon according to the sign of j. The density (and velocity) profile is flat for x < 0 while is monotonically increasing for x > 0: this solution corresponds to the "half soliton" profile in full analogy with the known results valid for the standard (cubic) GPE [35]. The analytical form of the half soliton for the quintic GPE is Φ(x) = n(x) e i 2 ϕ(x) where the modulus and the phase are expressed in terms of the unique parameter δ Along the positive semi-axis x > 0 we have 1 A different class of solutions of the cubic GPE, characterized by an oscillating density profile, was discussed in Ref. [35].
Although this class of solutions is also present in the quintitc GPE, we do not discuss their properties in detail.
where α = (4 − 3 δ 2 ) 1 2 (2 − 3 δ 2 ) −1 and the asymptotic density n ∞ is expressed in terms of the parameter δ by the algebraic equation: For x < 0, as previously stated, the density is constant and the phase linear in x. This analysis shows that the GPE is not able to reproduce, even qualitatively, the exact solution when a supersonic transition is present and the supersonic flows have non-monotonic profiles with a well defined limit for x → ±∞ (see Figure 3). Only the special case k 0 = k F appears to be related to the "half soliton" solution of the GPE equation, although the agreement is not quantitative (see Figure 4): if the parameter δ is chosen so that the x → +∞ limit of n(x) coincides with the exact value, the uniform density in the supersonic region is underestimated and the profile close to the transition point is not correctly reproduced.

Black hole formation
We now examine the formation dynamics of an event horizon in this analogue model. We first note that all the steady state solutions we have found, together with their time reversal, are dynamically stable because the exact dynamics is linear, being that a free Fermi gas. A direct numerical integration of the Schrödinger equation starting from the analytical steady state solution confirms this expectation: a small random noise superimposed to the initial condition does not grow in time. This result might appear puzzling because we have seen that, for suitable choices of the parameters (k F , k 0 ) a black-hole/white-hole pair appears (see, for instance, the case k 0 = 3π 10 Q in Figure 3). Under these circumstances a dynamical instability is generally expected due to the black-hole laser mechanism [36,37], although detailed studies [38][39][40] found that additional requirements must be met in order to trigger the instability. The stability of the exact stationary states is shared by the half soliton solution of the GPE: a Bogoliubov-de Gennes analysis shows that, analogously to the known result valid for the dark soliton solution in the cubic GPE [31,35], no complex eigenfrequency is present. The stability of the half soliton solution can be also easily checked numerically by direct integration of the GPE equation.

Exact dynamics.
Given the previous considerations, we might expect that a uniformly flowing TG gas perturbed by the sudden switching on of a step potential should evolve spontaneously towards one of the steady state solutions investigated in Section 4. Starting with a left moving fluid, this is indeed the case for the exact dynamics in all the different regimes previously identified as reported in Figure 5, where a few snapshots of the evolving density profile are reported for k F = π 2 Q and k 0 ≥ 0. As argued in our previous studies [26,27] and shown in Figure 5, the exact dynamics approaches the stationary state which we have determined by analytical means in the previous Section: i.e., the Slater determinant of the scattering states identified by the same parameters (k F , k 0 ) of the initial condition. The time scale necessary to reach the steady state depends on the size of the region we are interested in; we see that a few tens of the natural time unit (τ = m Q 2 ) are enough to converge to the analytic solution in the range |xQ| 10 and no shock waves appear during the exact time evolution. Thus, the dynamics drives a homogeneous, left moving TG gas to a known steady state which may exhibit the presence of sonic horizons, depending on how the initial density and velocity of the gas are chosen. Note again that the particular case of k 0 = k F shows a constant density downstream and the half-soliton shape discussed in the previous Section.

Semiclassical dynamics.
Remarkably, the case is quite different for the semiclassical dynamics, as shown in Figure 6. We know that only the fully subsonic stationary states and the half soliton solution have a direct analogue in the GPE, while there is no semiclassical solution corresponding to non monotonic density profiles and a well defined asymptotic limit at |x| → ∞. For an initial condition leading to a subsonic steady state (see for instance the k 0 = 0 case), the GPE evolution does indeed proceed smoothly, converging to a stationary state as expected. Furthermore, if we compare Figure 6 with Figure 2 and 5, we see that the agreement between the semiclassical solution and the exact one is also quantitative, as the values of the density, of the fluid velocity and of the speed of sound show good agreement.
Instead, when the initial velocity increases, the semiclassical evolution does not approach a steady state even in the case k 0 = k F (which is shown in Figure 6), where we have shown that the half soliton solution is stable and closely mimics the exact result; here, instead, a sonic transition is formed during the evolution, as expected, but the density profile shows oscillations in space and time in the supersonic region and the flattening of the profile does not occur. Moreover, by carefully examining the dynamics, it can be noticed that the envelopes constituting this pattern are emitted one at the time, they have equal width, they move downstream with a constant velocity and their shape does not change during the evolution: they thus resemble the soliton trains present in the cubic GPE [41,42]. This may not be evident from Figure 6, where the downstream behavior looks more like a regular oscillation than a collection of solitons, but the structure becomes much clearer if we take lower values of k 0 : in these cases the single envelopes are well separated and the succession of solitons is well defined. On the other hand, when we increase k 0 the interval between the emissions decreases but we recover the same qualitative behavior for all the intermediate cases k 0 < k F which feature a sonic horizon: there, in fact, the semiclassical equation does not allow for stationary solutions and indeed the soliton emission prevents to attain a steady state up to the k 0 = k F limiting case, where a stationary solution of the GPE equation (the half soliton) is present but is not reached by the GPE dynamics.

In search of the white hole
A white hole is the time reversal of the black hole [43] and we already noted that the corresponding stationary state is indeed a stable solution of the exact dynamical equations, which is in agreement with previous studies on analogue white holes in BECs [44,45]. So one might expect that for each initial condition leading to the formation of a black-hole horizon, starting from a time reversed initial state, the dynamical evolution of a free uniform TG gas would give rise to a white hole. In other words, if an initial condition defined by the two parameters (k F , k 0 ) leads to the formation of a black hole, does the time evolution of the state (k F , −k 0 ) give rise to a white-hole horizon? 6.1. Exact dynamics. In order to test this claim, we have to investigate whether the dynamics drives a right moving fluid flowing against a step potential towards a white-hole configuration. This is not the case, as already shown in Figure 3, where the properties of the steady states identified by the pairs (k F , ±k 0 ) are compared. To better illustrate this point, in Figure 7 we plot the exact dynamical evolution of a TG gas starting from different initial conditions defined by k F = π 2 Q and k 0 < 0, i.e. a right moving uniform fluid. When the potential step is switched on the density profile smoothly converges towards the analytical steady state shown in the right panels of Figure 3. The stationary flow may then be fully subsonic (see the case of k 0 = − 3π 10 Q), may display a finite supersonic region (see k 0 = − 9π 20 Q) or give rise to a black-hole horizon (see k 0 = −k F = − π 2 Q). We never observe the spontaneous formation of a white-hole horizon, i.e. a flow going from supersonic to subsonic when passing through the potential step: the external potential always impresses an acceleration on the flow, irrespective of the direction of motion. The observed dynamics can be understood in terms of scattering states: a white hole is formed by taking the time reversal of a black hole, which, due to the dynamics induced by the step, is a particular superposition of different states. Thus, a initially homogeneous fluid (i.e. a single "incoming" state) can never lead to a white-hole configuration.
6.2. Semiclassical dynamics. Also in this case, the GPE dynamics of a right moving fluid differs significantly from the exact result. In Section 4 we showed that, while the GPE admits a class of stable subsonic stationary solutions, only the "half soliton" wavefunction displays a sonic horizon, which may correspond either to a BH or a WH horizon according to whether the mass current is negative or positive. Figure 8 shows the main features of the dynamics according to the GPE: when the initial condition corresponds to a generic velocity in the range −k F < k 0 < 0, the evolution proceeds smoothly. Density waves are ejected towards x → −∞ and a subsonic stationary flow is reached at long times (see left panels). In the limiting case k 0 = −k F (right panels) the stationary flow becomes asymptotically sonic, with v = c at x → +∞. It is interesting to notice how, also in this case, the semiclassical picture agrees with the exact one also from a quantitative point of view only in the subsonic regimes: indeed, by comparing Figure 8 with Figure 3 and 7 it is clear how the values of the density, the gas velocity and the sound speed are comparable whenever the flow is subsonic. Yet, when a supersonic transition sets in, the quantitative agreement is spoiled.
Finally, we want to address an important point regarding the semiclassical treatment: recalling Eq. (2), the question may arise whether the dynamical behavior of the gas in presence of a sonic transition may be due to the form of the non-linearity, while an approach based on the canonical GPE would have given a different picture. Figure 9 shows that this is not the case, as the main features (horizon formations, oscillations etc.) observed in the semiclassical dynamics are not a peculiarity due to the precise form of Eq. (2). Indeed, the same qualitative behavior is retrieved when the above analysis is repeated with the standard (cubic) Gross-Pitaevskii equation and the quantitative differences between the two cases are irrelevant to the conclusions drawn in this study. This is important as theoretical studies on sonic horizons in BECs usually rely on the standard GPE and a description of the dynamical behavior of a condensate in these regimes has been deemed useful several times in this context [35,46,47]. To the best of our knowledge, only two other works investigated the formation dynamics of the sonic horizons [48,49] but the results only partially agree with ours due to the particular configurations chosen 2 . Figure 9. The two figures compare the long time behavior of the density when the usual (cubic) GPE is used (blue curves) in contrast with the modified (quintic) form of Eq. (2) (red curves). In both cases, the same dimensionless coupling constantg = π 2 has been used. The profiles are qualitatively the same for both plots. The top plot shows a subsonic flow with k 0 = 0 and k F = π 2 Q; the bottom plot, instead, is the case in which a sonic horizon form, i.e. k 0 = k F = π 2 Q. Similar results are obtained in the case of a right moving fluid. All the curves refer to a time t = 20 m Q 2 .

Further possible configurations
The analysis carried forward in the previous Sections starts from the assumption that |k 0 | ≤ k F , which means that initially, when the potential is switched off, the fluid has a subsonic (at most sonic) velocity. For completeness, we will now relax this hypothesis and investigate what happens if the gas initially flows at supersonic speeds. Eqs. (12,13,14) provide the expected density profile and the particle current at stationarity also for |k 0 | > k F .
Following the previous discussion, we can study the dynamics of the gas after the quench numerically and we can investigate the possible formation of a sonic horizon both within the exact dynamics and according to the semiclassical (GPE) approximation. One would expect that a supersonic homogeneous gas subject to an external waterfall potential would not give rise to a sonic transition, but, instead, the switching-on of the step potential introduces various interesting configurations as we will now describe in some detail.
Let us start with a left moving gas flowing over the potential step. If we choose a value of the initial gas velocity which is slightly bigger than the sound speed (k 0 k F ), the exact evolution drives the system towards Figure 10. Absolute value of the fluid velocity (red) and of the sound speed (blue) according to the GPE dynamics. The initial state is a homogeneous fluid with k F = π 2 Q and different values of k 0 . Left panels correspond to a left moving fluid (k 0 > 0), right panels to right moving fluid (k 0 < 0, i.e. against the step). Velocities are in unit of Q m . The snapshots are taken at a time t = 60 m Q 2 for the right moving cases, where the GPE slowly drives the system towards a stationary state. For a left moving fluid, the snapshots correspond to t = 15 m Q 2 . In these cases, the system reaches a steady state in the upstream region, while the velocity profiles oscillate in time (and space) downstream. Snapshots taken at a time t = 14 m Q 2 are also shown (magenta and cyan lines) to emphasize the temporal oscillations. The black dashed line represents the stationary state reached by the exact evolution in each case. a stationary state displaying a sonic horizon. By further increasing k 0 this region shrinks until the fluid is supersonic on the whole axis. The dashed line of Figure 10 illustrates this behavior: indeed, the upper left panel shows a value of the initial velocity for which the subsonic region consists in one point, while in the central and bottom left panel the flow is always supersonic. Instead, the long time behavior of the semiclassical (GPE) solution shows an oscillation pattern in the upstream region and a train of solitons in the downstream region, thus never leading to the exact stationary state. This is found to be true for any value of k 0 > k F , as Figure 10 shows. Also, as noted in Section 4, under this circumstances the exact solution is always flat in the downstream region. The right panels of Figure 10, instead, report the long time evolution of a right moving fluid. The exact results (dashed lines) show that in a range of velocities |k 0 | > k F a BH horizon is formed: the flow upstream is subsonic and a supersonic transition occurs near the position of the step. By further increasing the value of the initial velocity the subsonic region becomes finite until it disappears, leaving a supersonic flow everywhere. Interestingly, it can be seen that both the exact and the semiclassical solution tend to a finite value in the downstream region while the two are significantly different in the upstream region: while the exact solution has smooth behavior, the semiclassical one develops spatial oscillations. Also, as already noted in Section 4, the exact solution is always flat in the downstream region if −k 0 > Q + k F (bottom panel).

Conclusions
We have studied, both analytically and numerically, an exactly solvable model describing a flowing onedimensional Bose gas, in order to shed light on the physics underlying the analogue gravity paradigm. This model of hard-core bosons enables to follow the formation dynamics of black-and white-hole horizons after an external step potential is switched on: the analogue of a gravitational collapse event. Indeed, the fluid starts from a situation of constant density and velocity and, after a quench, it reaches another stationary configuration. From a curved spacetime point of view, this process is the exact analogue of Hawking famous gedankenexperiment [3]. Furthermore, the shape of the external potential represents a realistic reproduction of existing experimental setups [9][10][11][12][13] and a model used in various theoretical works (for example [46][47][48][49]).
We analytically found the stable stationary states of this model which allows for the existence of both blackand white-hole horizons, as well as black/white hole pairs. By numerical integration of the exact evolution equations we found that if a horizon is formed after the quench, either an isolated black-hole horizon or a tight black/white hole pair is formed but a white-hole horizon can never be reached dynamically. This somewhat conforms with the analogue gravity paradigm: the positive black-hole entropy allows for the spontaneous formation of an horizon, while the negative white-hole entropy does not. Moreover this study suggests that a sufficiently tight black/white hole pair is not necessarily unstable, at least in analogue models, while our solutions confirm that both a single black/white hole is a stable solution of the dynamical equations.
Furthermore, we compared these findings with the solutions of a semiclassical equation usually adopted to describe this class of models: i.e., the Gross-Pitaevskii equation, here suitably generalized to deal with a strongly interacting Bose gas. We found that only a subset of all possible stationary states are described by the GPE, including the subsonic solutions (with no event horizons) and the single black hole case, for a specific value of the mass current. However, the dynamics described by the GPE does not lead to the formation of a stationary solution whenever a sonic transition is present: indeed, starting from a uniformly flowing Bose gas, the semiclassical evolution gives rise to density waves propagating outwards which are not damped in time, preventing the convergence towards the stable stationary solution of the GPE. This behavior, which is a general feature also for the usual cubic GPE, suggests that results based on the integration of non-linear Schrödinger equations should be critically considered, at least as far as the long time dynamics is concerned.
Finally, to complete the description, we show other possible configurations leading to the formation of one (or more) sonic horizon which can be achieved by further increasing the initial velocity of the gas.
On a final note, the results obtained in this work provide an important step for future possible investigations as they give insights for new possible experimental configurations, they represent a key ingredient for the study of the hydrodynamical instability (and gravitational effect) known as black-hole laser and they help understanding the bounds on validity of semiclassical methods. To this extent, a question which arises from these results is what the optimal black-hole configuration should be in order to study the analogue of the Hawking effect. In fact, on one hand, as Hawking already pointed out [2], the dynamical history of the gravitational collapse should not affect the emitted radiation, but, on the other hand, his line of thought heavily rests on the comparison between a precise initial and final stationary state, which we are able to identify through the exact solution of the model but not by means of semiclassical approaches.