{'fraction_non_alphanumeric': 0.03049265341400173, 'fraction_numerical': 0.010475367329299913, 'mean_word_length': 4.321192052980132, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Simulation of Acoustic Black Hole in a Laval Nozzle

A numerical simulation of fluid flows in a Laval nozzle is performed to observe formations of acoustic black holes and the classical counterpart to Hawking radiation under a realistic setting of the laboratory experiment. We determined the Hawking temperature of the acoustic black hole from obtained numerical data. Some noteworthy points in analyzing the experimental data are clarified through our numerical simulation.


I. INTRODUCTION
The existence of the black hole horizon is essential to the Hawking radiation [1,2,3,4].
The horizon causes an extremely large redshift on an outgoing wave of a matter field during propagating from a vicinity of the horizon to the asymptotically flat region. This redshift results in the Planckian distribution of quantum mechanically created particles at the future null infinity. Although the Hawking radiation has not yet been confirmed observationally, it has been pointed out that the similar phenomenon appears on phonons in a transonic flow of a fluid. This theoretical phenomenon in the fluid is called the sonic analogue of Hawking radiation, and the transonic fluid flow on which the sonic analogue of Hawking radiation occurs is called the acoustic black hole [5,6,7,8,9,10,11,12]. In the transonic flow, the boundary between a subsonic and a supersonic region is the sonic point. The sonic point corresponds to the black hole horizon and called the sonic horizon. The region apart form the sonic point in the subsonic region corresponds to the asymptotically flat region in a black hole spacetime. The sound wave which propagates from a vicinity of the sonic point to the subsonic region receives the extremely large redshift and the sonic analogue of Hawking radiation appears on phonons. The Hawking temperature T H of the acoustic black hole is given by the gradient of the fluid velocity at the sonic point x = x s [5]: where c s and v are the sound velocity and the fluid velocity, respectively. For an ordinary system with typical size ∼ 1 m and c s ∼ 340 m/s, we obtain T H ∼ 4 × 10 −10 K. Thus, in a practical experiment in a laboratory, it is very difficult to detect the sonic analogue of Hawking radiation with such an extremely low temperature.
While it is very hard to detect the sonic analogue of Hawking radiation on phonons, which is a quantum origin, a classical sound wave with a sufficiently large amplitude can be available for detecting the effect of the acoustic black hole on a propagation of sound waves in practical experiments. Indeed, the thermal nature of Hawking radiation is not a quantum origin; it comes from a large redshift on the wave propagating from a vicinity of the horizon to the asymptotically flat region, and has nothing to do with the quantum effect (see Appendix). Therefore, if we can construct a classical quantity from sound waves which shows the thermal distribution, then such a quantity can be interpreted as a classical counterpart to Hawking radiation. A candidate for the classical counterpart has already been introduced in the references [13,14] by using a Fourier component of an outgoing wave from a vicinity of the horizon.
By the way, since a realistic fluid is composed of molecules, the effect of the molecule size discrepancy of the fluid becomes one of the important issues of the acoustic black hole. This effect appears in the modified dispersion relation at high energies of phonons. However, in the expansion by the molecule size discrepancy, the zeroth order form of the quantum expectation value of the number operator 0| N |0 of phonons results in the thermal spectrum at least not so high energy region. That is, the effect of the molecule size discrepancy is not dominant in the sonic analogue of Hawking radiation. Therefore, before proceeding to the issue of the molecule size discrepancy, it is necessary to detect the ordinary Hawking radiation (thermal spectrum) by a real experiment in laboratory, and the molecule size discrepancy is not in the scope of this paper. Hence we concentrate on the construction of the practical procedure to extract the evidence of the Hawking radiation (thermal spectrum) from the experimental data. To do so, we make use of the classical counterpart to Hawking radiation.
We should emphasize here that the Hawking temperature (1) In this paper, we aim to construct a practical procedure of the data analysis to obtain the classical counterpart to Hawking radiation in a realistic setting of the acoustic black hole. To do so, we perform a numerical simulation of a fluid flow in a practical experimental setting, and demonstrate that our procedure works well to detect the classical counterpart to Hawking radiation.
We organize the paper as follows. In the section II, we review the acoustic black hole in a Laval nozzle and define the classical counterpart to Hawking radiation. The section III is devoted to the introduction of our numerical method to simulate a transonic flow, and to the construction of the practical procedure of the data analysis. In the section IV, our numerical results are presented and we determine the surface gravity of the sonic horizon.
Finally the summary and conclusion are in the section V.

II. ACOUSTIC BLACK HOLE WITH A LAVAL NOZZLE
We consider an acoustic black hole with a fluid in a Laval nozzle as proposed in [14]. The Here we note that, if a stationary transonic flow is prepared from the beginning (which corresponds to the eternal black hole), we can not expect to obtain the classical counterpart to the Hawking radiation; for the eternal case, the tunneling of phonons across the sonic point can result in the Hawking radiation and this is the purely quantum effect. Hence, when we are interested in the "classical" counterpart to Hawking radiation, it is necessary to consider the situation of the sonic point formation in course of the dynamical evolution of the fluid flow . As seen below, the sound wave which is prepared on the fluid flow before the formation of the sonic point becomes the classical counterpart to the quantum fluctuation which causes the classical counterpart to Hawking radiation.

A. Basic equations
We consider a perfect fluid and treat the flow in a Laval nozzle as one dimensional for simplicity. The basic equations are the mass conservation equation and the Euler equation: where ρ(t, x) is the mass density, v(t, x) is the fluid velocity, p(t, x) is the pressure of the fluid and S(x) is the cross section of the Laval nozzle. We assume the adiabatic ideal gas type equation of state p ∝ ρ γ where γ is the adiabatic index. The sound velocity c s is given by For a stationary background flow, we can obtain the sound velocity c s and the cross section where quantities with subscript "in" represent the values at the inlet of the nozzle. The spatial derivative of the Mach number at the sonic point x = x s is and this value relates to the surface gravity of the sonic horizon by Eq. (2). Figure and the surface gravity given by Eq. (2) is Here it should be emphasized that, because the cross section of the nozzle is set S(x) = 2 − cos(πx) with −1 < x < 1, the spatial scale is normalized by L/2, where L is the length of the nozzle, and that the temporal scale is normalized by L/(2c in ), where we set c in = 1.
We have not considered the perturbation of the fluid flow so far. In the following sections, the perturbation (sound wave) is introduced on the fluid flow, and the classical counterpart to Hawking radiation is defined. Then in the section IV, we compare this surface gravity (8) with the surface gravity of the sonic horizon which is obtained through the observation of the classical counterpart to Hawking radiation in our numerical simulation. B. Classical counterpart to Hawking radiation As originally shown in the reference [5], the perturbation of the velocity potential for This corresponds to the Klein-Gordon equation φ = 0 with the acoustic metric where the location of the sonic horizon (sonic point) x = x s is given by Here we should note that, because we consider the one dimensional flow along x-axis as mentioned at the beginning of the previous section II A, the two dimensional section of (y, z) coordinate is omitted in obtaining Eq. (9) from φ = 0. In the following analysis, we omit the two dimensional section of (y, z) coordinate.
In order to describe the extremely large redshift which causes the classical counterpart to Hawking radiation, we introduce the following null coordinates for the first, and the acoustic metric becomes ds 2 = (S ρ/c s ) (c 2 s − v 2 ) du dw, which has a coordinate singularity at x = x s . In order to eliminate the coordinate singularity, we introduce the new null coordinates where ′ = ∂ x and quantities with subscript 0 denotes values at the sonic point x = x s . Then the acoustic metric and the evolution equation (9) of the perturbation φ become The form of this metric near the sonic horizon is x is used. This denotes explicitly that the coordinate (U, W ) does not have coordinate singularities and corresponds to the Kruskal-Szekeres coordinate of the Schwarzschild spacetime. Here it should be recalled that, in the spacetime of the gravitational collapse (the left panel in FIG. 1), the Kruskal-Szekeres coordinate is appropriate to describe the rest observer at the spacetime region before the formation of the black hole [2].
Hence the outgoing normal mode of the sound wave which corresponds to the zero point fluctuation of the quantized matter field before the formation of the black hole (i.e. the mode function on the flat spacetime) is given by where ω 0 is the initial frequency and A φ is the amplitude of this mode. In the original coordinate (t, x), this mode behaves as This yields the temporal wave form at the observation point x = x obs in the asymptotic region: Fourier component of this mode is defined by and the power spectrum is where the quantity ω H is defined by Here it should be noted that, because the existence of the sonic point is assumed in the above discussion, the form of P (ω) of Eq. (19) must be obtained after the formation of the sonic horizon in the fluid flow.
One may think it is strange that the initial frequency ω 0 does not appear in the power spectrum (19). It is the fact that the initial frequency ω 0 does appear in the power spectrum in the context of the quantum field theory, because the amplitude A φ is determined by the normalization condition for the mode functions to be A φ ∝ 1/ √ ω 0 . However, in the calculation for the classical counterpart to Hawking radiation, the initial frequency ω 0 appears only in the phase of the Fourier component φ out ω and the power P (ω) does not explicitly depend on ω 0 .
The power spectrum (19) has the Planckian distribution for ω < 0 and this is the classical counterpart to the quantum Hawking radiation [13,14]. For the quantum Hawking radiation, the magnitude of the power is determined by the zero point oscillation of the quantized field.
However, for the classical counterpart to Hawking radiation, the magnitude of the power depends on the amplitude A φ of the input signal and it is expected that we can detect the thermal distribution of the power spectrum for the sound wave of a sufficiently large amplitude in practical experiments.
As long as concerning one dimensional fluid flow, the amplitude of the sound wave does not decrease during propagating from the outlet to the inlet of the nozzle. The one dimensional transonic fluid flow corresponds to the two dimensional black hole spacetime in which no curvature scattering occurs. The theoretically expected power (19) is equivalent to the formula of the Hawking radiation for black holes derived ignoring the curvature scattering.
Our purpose is to observe this thermal spectrum via the numerical simulation of a transonic flow in the Laval nozzle, and to construct the practical procedure of the data analysis.
Because the observable sound wave is a real valued one, we prepare the following real valued input mode at the outlet of the nozzle before the formation of the sonic horizon: where θ is the initial phase of the wave, φ in (U) is given by Eq.
The effect of the redshift due to the formation of the acoustic black hole gradually appears on the observed sound wave. The temporal evolution of the observed sound wave is affected by the formation of the acoustic black hole. That is, the information of the classical counterpart to Hawking radiation is encoded in the observed sound wave. Therefore, using the output form (17), the resulting output signal after the formation of the sonic horizon is given by and referring Eq. (18), the Fourier component of this output signal is given by Hence, because the initial phase θ appears in φ , we have to superpose two output modes φ out(θ) and φ out(θ+π/2) , and obtain It is obvious that the power spectrum of this Fourier component gives the same form as Eq. (19). Therefore, in order to obtain the Planckian distribution in the power spectrum P (ω) of the observed sound wave from real valued input signals, we have to combine at least two output modes of which the input phases differ by π/2.

III. NUMERICAL SIMULATION OF TRANSONIC FLOWS IN A LAVAL NOZ-ZLE
Our numerical simulation is designed to generate a one dimensional transonic fluid flow from the initial configuration with no sonic point. We prepare the input signal at the outlet of the nozzle and observe the sound wave at the inlet of the nozzle. At the same time, the practical procedure of the data analysis is constructed.

A. Basic equations for numerical simulations
For the numerical calculation, we rewrite the fluid Eqs.
(3) using c s and v, This set of equations is not suitable for a numerical simulation of wave propagation in a transonic flow. We transform these equations to the advection form. We introduce the Riemann invariants J ± (t, x) as follows Then the basic equations become where The left hand side of Eqs. (28) are of the advection form. That is, if ∂ x S = 0 and ∂ t V ± = 0, J + propagates upward along t − dx/V + = const. and J − propagates downward along t + dx/V − = const. This form of the equations is suitable to treat the propagation of waves numerically.
In the sonic analogue of Hawking radiation, the perturbation of the velocity potential (sound wave) corresponds to the scalar field in a black hole spacetime as shown by Eq. (9).
Therefore, in a practical experiment of an acoustic black hole, we should calculate the velocity potential by integrating the observed velocity of the fluid: where x obs is the observation point of sound waves and x 0 defines the origin of the velocity potential. However, since our experimental setting is designed to observe the sound wave at a fixed spatial point x = x obs , the observational data is a temporal sequence of v(t, x obs ) and we must devise the other method to obtain the velocity potential at the observation point.
In the upstream subsonic region of the transonic flow where the effects of the sonic horizon is negligible and the flow is stationary, the sound wave propagates along the "null" direction as indicated by Eq. (11) and we can obtain the value of the velocity potential at x obs by integrating the velocity with respect to time where t 0 defines the origin of the potential. We use this formula to evaluate the velocity potential at the observation point.

B. Experimental setting and numerical method
We use the Laval nozzle with the cross section of the following form, We make the fluid flow against x-axis and assume v < 0. The inlet of the nozzle is at x = 1, the throat at x = 0 and the outlet at x = −1. The input wave is prepared at the outlet of the nozzle and the wave propagates against the flow from x = −1 to x = 1. After the flow settles down to a stationary transonic flow, since ∂ x S| x=1 = 0 at the inlet, the set of Eqs. (28) gives ∂ x v| x=1 = 0 at the inlet. This means that the inlet of the nozzle corresponds to the asymptotically flat region in the gravitational collapse spacetime. Hence we make the observation of the sound wave at the inlet and set x obs = 1.
We prepare a flow with a homogeneous velocity distribution as the initial configuration, This initial configuration of flow has no sonic point and corresponds to the flat spacetime region before the formation of a black hole. Since our evolution equations (28) are of the advection form, it is enough to set the boundary conditions for J + at the outlet and J − at the inlet: where A J is a constant that represents the amplitude of the perturbation of J + at x = −1, where c s1 and v 1 are respectively the sound velocity and the fluid velocity at x = 1, and the relation of the perturbations φ(t, x = −1) = δΦ(t, −1) = (c s1 + v 1 ) t dt δJ + (t, −1)/2 at the outlet is used (see Eqs. (27) and (32)).
In order to extract the perturbation part φ from the "full" velocity potential Φ of Eq. (32) which includes the "background" flow of the fluid, we should prepare three input modes with different initial phases θ 1 = −π/2, θ 2 = 0 and θ 3 = π/2. If the amplitude A J of the perturbation part of J + is small enough, then the backreaction effects in the velocity potentials Φ out(θ i ) (i = 1, 2, 3) at the observation point (x = 1) are negligible and we can subtract the common background contribution by taking the difference among them. Hence we obtain the Fourier components of the perturbation (sound wave) at the observation point as follows where we used cos(ω 0 t) − cos(ω 0 t ± π/2) = √ 2 cos(ω 0 t ± π/4). Hence we obtain the pure Fourier component φ out ω of the output sound wave by Eq. (25), and its power spectrum is given by where the factor 1/2 is introduced to eliminate the factor √ 2 which appears in the amplitude √ 2A J of the input signals (38). When we carry out a realistic experiment of the acoustic black hole, the power spectrum (39) gives the classical counterpart to Hawking radiation and it has to be compared with the theoretical form (19). Hence, it is recognized that we have to perform at least three independent experiments with three different input phases and observe three independent output signals to obtain the classical counterpart to Hawking radiation.
With our setting of the fluid flow in the Laval nozzle, we carry out the numerical simulation by the finite difference method. We use the Cubic-Interpolated Pseudoparticle (CIP) method [15] for the interpolation between neighboring spatial mesh points. This method enables us to treat the shock which appears after the formation of the sonic point. We note here that the shock arises in the supersonic region and it does not affect the subsonic region where we observe the output signal. The procedure of the numerical simulation to observe the thermal power spectrum (19) of the perturbation of the velocity potential is as follows: Step 1: Generate the transonic fluid flow three times with three different initial phase values, θ 1 = −π/2, θ 2 = 0 and θ 3 = π/2. Then, using Eq. (32), obtain the full velocity potentials Φ out(θ i ) (i = 1, 2, 3).
Step 3: Calculate the power spectrum This is the observationally obtained power spectrum for the perturbation of the velocity potential.
Step 4: Then compare this power spectrum with the theoretically expected one and determine the surface gravity κ H of the sonic horizon. Here note that the quantitỹ ω 0 in the amplitude of P theory is the effective input frequency given by the following discussion.
We must note here that there arises an extra redshift effect on the sound wave. In our numerical simulation, the initial fluid flow has homogeneous velocity distribution. Then the system starts to evolve dynamically and settles down to the stationary transonic fluid flow finally. Therefore, the observed frequency shifts due to the evolution of the "background" flow. According to the null coordinate (11), the outgoing wave with frequency ω 0 is and the effective frequency of the wave at x = 1 is given bỹ The second term in Eq. (43) gives a negative contribution in our non-stationary setting of numerical experiment and the observed frequencyω 0 at x = 1 becomes smaller than ω 0 .
The theoretically expected form of the power spectrum in our setting should be given by replacing ω 0 toω 0 in the formula (19).

B. black hole formation case
For the initial fluid velocity V i = −0.2, we observe the formation of the acoustic black hole. Figure 10 shows the evolution of the Mach number.  The spacetime distribution of the time derivative of the fluid velocity and the perturbation is shown in FIG. 13 and FIG. 14, respectively. The time evolution of the power spectrum of the observed signal is shown in FIG. 16,   FIG. 17 and FIG. 18. The effective frequencyω 0 of the input perturbation observed at the inlet of the nozzle is 74.6. Figure 16 shows that, after t ∼ 10, the power spectrum of the observed perturbation spreads out toward the low ω range by the redshift effect due to the formation of the sonic horizon, and the divergence of the power P obs (ω) appears at ω = 0.
This indicates that the observed spectrum approaches the theoretically expected form (41) which diverges at ω = 0 as ∝ 1/ω. To remove this divergence, we plot ω 2 P obs (ω) in FIG. 17   and FIG. 18. C. Surface gravity of the sonic horizon Figure 19 shows the observed power spectrum ω 2 P obs (ω) at t = 18 given by Eq. (40) and the theoretically expected power ω 2 P theory (ω) given by Eq. (41).
This value is consistent with the theoretically expected value 0.355 given by Eq. (8).
Here we discuss about two points. The first is about the oscillation of the numerical power spectrum. This oscillation is due to the finite size Fourier transformation. This oscillation will disappear if the numerical calculation is carried out with larger value of ω 0 .
The second point is the deviation of the numerical power from the theoretical one in the high ω range. This deviation is due to the input sound wave given in the boundary condition Through our numerical simulation, we have obtained two noteworthy points for data analysis of the experiments of the acoustic black hole: the first one is that a single input wave can not give us necessary information of the classical counterpart to Hawking radiation.
We need to carry out several independent observations of the sound waves with different initial phases to retrieve the thermal distribution of the power spectrum. The second one is that we can evaluate the velocity potential at the observation point by integrating the fluid velocity at that point with respect to time (see Eq. (32)). Therefore, it is not necessary to observe the sound wave at every spatial points in the fluid.
In our calculation presented here, we have considered the sound wave with small ampli-tude. However our numerical code is applicable beyond perturbation, in which the non-linear effect of the sound wave becomes important. By analyzing such a situation, we may be able to discuss the backreaction effect on the classical counterpart to Hawking radiation. As the next step, we are planning to extend our simulation of transonic flows in a Laval nozzle to include quantum effects. Then, we expect to obtain implications for quantum aspects of Hawking radiation using the sonic analogue model of black hole.

APPENDIX: CLASSICAL AND QUANTUM EFFECTS IN HAWKING RADIA-TION
We briefly review the Hawking radiation in a spacetime of a gravitational collapse forming a Schwarzschild black hole, and clarify the distinction between the quantum effects and the classical effects in the occurrence of the Hawking radiation. For simplicity, we consider a massless free scalar field Φ as a representative of the matter field, and set c = k B = 1.
For the first, we start with the classical effects. Because the spacetime is dynamical, the positive frequency modes and the negative frequency modes of the scalar field are mixed as the system evolves in time. This mixing is represented by the Bogoliubov transformation between the positive frequency mode φ ω of Φ at the past null infinity and that modeφ ω at the future null infinity whereφ * ω =φ −ω is the negative frequency mode at the future null infinity. The Bogoliubov coefficients A and B are obtained by solving the classical wave equation Φ = 0 with appropriate boundary conditions at the past and future null infinities and with the following relation, which comes from the normalization condition (φ ω 1 , φ ω 2 ) = δ ω 1 ω 2 with respect to the Klein-Gordon inner product of Φ [2]. This means that the derivation of A and B is purely classical.
Next, in the classical framework, we consider the wave mode which propagates from the past null infinity to the future null infinity via a vicinity of the black hole horizon (see the left panel in FIG. 1). This wave mode is ingoing at the past null infinity and becomes outgoing at the future null infinity. The ingoing positive frequency mode φ ω at the past null infinity is given by where w is the ingoing null coordinate appropriate for the rest observer at the past null infinity. When this wave mode φ ω (w) propagates to the future null infinity via a vicinity of the horizon, it evolves to be φ ω (w(u)) at the future null infinity under the geometrical optics approximation, where u is the outgoing null coordinate appropriate for the rest observer at the future null infinity. The function w(u) expresses the extremely large redshift which the mode φ ω (w) receives during propagating from the past null infinity to the future null infinity. This redshift effect can be decomposed into two parts: one is the redshift during the propagation from the past null infinity to the vicinity of the horizon, and the other part is the redshift after passing through the vicinity of the horizon. The first contribution is not so large and the mixing of positive and negative frequency modes does not occur. However, the second contribution is large enough to make the Bogoliubov coefficient B ωω non zero.
By matching the null coordinates w and u along a null geodesic which connects the past and future null infinities via the vicinity of the horizon, the function w(u) is obtained [1,2]: where κ is the surface gravity of the black hole horizon, δ is an arbitrary constant denoting the freedom of choosing the origin of w, and the constant α is determined by the first part of the redshift. The exponential form exp(−κu) in the Eq. (A.4) comes from the second part of the redshift and implies that the wave length of the outgoing wave is exponentially stretched during propagating from a vicinity of the horizon to the future null infinity. That is, the wave φ ω (w(u)) is no longer a pure positive frequency mode but becomes a superposition of positive and negative frequency modes at the future null infinity. In order to calculate the superposition, we need the outgoing positive frequency modeφ ω at the future null infinity: φ ω (u) = 1 √ 4πω exp (−iω u) , ω > 0. (A.5) Then, using the definition of the Bogoliubov transformation (A.1), the wave φ ω (w(u)) at the future null infinity is decomposed as The Bogoliubov coefficients are obtained using the inner product (φ ω 1 ,φ ω 2 ) = A ω 1 ω 2 and (φ ω 1 ,φ * ω 2 ) = −B ω 1 ω 2 , and we find |A ω 1 ω 2 | 2 = e ω 2 /ω H |B ω 1 ω 2 | 2 , |B ω 1 ω 2 | 2 = 1 2πκω 1 1 e ω 2 /ω H − 1 , ω H = κ 2π . (A.7) The square of the Bogoliubov coefficients do not depend on the constants α and δ, and the following relation holds: |B ω 1 ,−ω 2 | 2 = |A ω 1 ω 2 | 2 . (A.8) All the above phenomena are the classical effects.
Finally we proceed to the quantum effects. When Φ is quantized, the harmonic operators a ω and a † ω with respect to the past mode φ ω are related to thoseã ω andã † ω with respect to the future modeφ ω as follows [1,2]: Then the number of particles at the future null infinity is obtained where |0 is the vacuum state at the past null infinity, a ω |0 = 0. Therefore, for the mode φ ω which passes a vicinity of the horizon and propagates to the future null infinity, using Eqs. (A.7) and an appropriate regularization method [1], we obtain N ω = 0|ã † ωã ω |0 = e ω/ ω H − 1 . (A.12) It is concluded that a black hole emits a thermal radiation of Φ with the Hawking temperature T H ≡ ω H .
It should be emphasized that, although the creation of particles is just the quantum effect, however the Planckian distribution (A.12) of the emitted particles is purely the classical effect due to Eqs. (A.7). The thermal nature of the spectrum comes from the Bogoliubov coefficient |B ω 1 ω 2 | 2 which has the Planckian distribution with respect to ω 2 . That is, the thermal nature of the Hawking radiation comes from the extremely large redshift for the wave mode φ ω which passes the vicinity of the horizon. If we take the classical limit → 0, the number of created particles and the Hawking temperature become zero but |B| 2 is not affected.