Analogues of gravity-induced instabilities in anisotropic metamaterials

In the context of field theory in curved spacetimes, it is known that suitable background spacetime geometries can trigger instabilities of fields, leading to exponential growth of their (quantum and classical) fluctuations --- a phenomenon called $\textit{vacuum awakening}$ in the quantum context, which in some classical scenarios seeds $\textit{spontaneous scalarization/vectorization}$. Despite its conceptual interest, an actual observation in nature of this effect is uncertain since it depends on the existence of fields with appropriate masses and couplings in strong-gravity regimes. Here, we propose analogues for this gravity-induced instability based on nonlinear optics of metamaterials which could, in principle, be observed in laboratory.

In the context of field theory in curved spacetimes, it is known that suitable background spacetime geometries can trigger instabilities of fields, leading to exponential growth of their (quantum and classical) fluctuations -a phenomenon called vacuum awakening in the quantum context, which in some classical scenarios seeds spontaneous scalarization/vectorization. Despite its conceptual interest, an actual observation in nature of this effect is uncertain since it depends on the existence of fields with appropriate masses and couplings in strong-gravity regimes. Here, we propose analogues for this gravity-induced instability based on nonlinear optics of metamaterials which could, in principle, be observed in laboratory.

I. INTRODUCTION
The influence of a background material medium on the propagation of mechanic and electromagnetic waves is well known to be formally analogous to that of an effective curved spacetime geometry. This idea was first presented, in the electromagnetic/optical context, by Gordon in 1923 [1] and it has since been developed in a number of different scenarios, particularly after Unruh's [2] and Visser's [3] works on acoustic analogues of black holes and their associated Hawking-like radiation. More recent applications of this formal analogy include mimicking in material media quantum lightcone fluctuations [4] and anisotropy in cosmological spacetimes [5]. The most appealing feature of these condensed-matter analogues of gravitational backgrounds is the possibility of observing in laboratory subtle but conceptually interesting effects which can be virtually unobservable in their original contexts -Hawking radiation being certainly the most emblematic among them, with claims of having already been observed in laboratory [6][7][8].
An interesting effect in the context of (quantum) fields in curved spacetimes is the triggering of field instabilities due to the background spacetime geometry -a phenomenon called vacuum awakening in the quantum context [9][10][11][12]. These gravity-induced instabilities exponentially amplify vacuum fluctuations to the point they decohere and seed classical perturbations [13], which, depending on field parameters, eventually evolve to a nonzero classical field configuration ("spontaneous scalarization" in the case of scalar fields [14][15][16][17]), stabilizing the whole system. More recently, this mechanism was also predicted to occur for massless spin-1 fields through appropriate nonminimal couplings [18] and, in analogy with the scalar case, the stabilization process was termed "spontaneous vectorization." To the best of our knowledge, condensed-matter and optical analogues of these gravity-induced instabilities have not been proposed to * Electronic address: caiocesarribeiro@ifsc.usp.br † Electronic address: vanzella@ifsc.usp.br this date. In this work, we propose and explore possible analogues of gravity-induced instabilities in the context of electromagnetism in polarizable/magnetizable anisotropic (meta)materials. Electromagnetic instabilities in flat spacetime are expected to occur in some materials. One celebrated example appeared in the context of plasma physics in the late 1950s and became known as Weibel instability [19]. The system, a neutral plasma whose components have anisotropic velocity distribution, possesses growing electromagnetic transverse waves. Related effects have been studied since then, with recent applications to solar plasma instability [20] and solid state devices [21]. Moreover, causal aspects of classical propagation in active materials were discussed in Ref. [22], where properties of the refractive index were established. Nevertheless, besides the fairly recurrence in the literature, usually quantization in such scenarios is not considered [23][24][25] or it is regarded as inconsistent [26,27].
It is noteworthy that instability of the electromagnetic field is always accompanied by evolution of the background, ending with the stabilization of the system as a whole. In the case of gravity-induced instability, the gravitational field changes with time, whereas electromagnetic instability in the presence of plasmas involves growing plasmons. In the case of electromagnetic fields in the presence of matter, for whatever form of the interaction with the background, the field's evolution is ruled by Maxwell's equations in the presence of polarizable/magnetizable media, and the interaction with the background is encapsulated in the functional dependence of the electric displacement (magnetic) vector field D (H) with the true (microscopic) fields E and B. If the magnitudes involved are small (e.g., in the beginning of the instability action), these functional relations become linear and one may find the form of the coefficients for such systems. For the case of Weibel instability, for instance, if the velocity anisotropy is taken in the z direction, the instability is modelled by a negative squared refractive index in the direction perpendicular to z.
We apply Gordon's method to propose a family of optical-based analogue models for electromagnetic fields presenting instabilities in curved spacetimes. We show how anisotropies of the background enter the effective equations in the form of nonminimal couplings, and in the case of strong anisotropy (just like for the Weibel instability), this coupling results in unstable solutions. We also discuss that for these systems the stabilization process occurs through the nonlinear nature of the background, which may seed spontaneous vectorization in analogy to the Einstein's field equations in the gravitational scenario.
The paper is organized as follows. In Sec. II, we present the covariant formalism of electromagnetism in anisotropic polarizable/magnetizable materials, establishing the formal analogy with nonminimally-coupled electromagnetism in curved spacetimes. In Subsec. II A, we consider a particular type of nonminimal coupling inspired by one-loop quantum electrodynamics (QED) corrections to electromagnetism in curved spacetimes. In Sec. III, we apply the formalism presented in the previous section to the scenario of a plane-symmetric anisotropic medium at rest in an inertial frame. Although plane-symmetric curved spacetimes (in four dimensions) are not really (physically) appealing, we consider this scenario for its simplicity and for its possible implications for the physics of the material medium. We construct the electromagnetic quantum-field operator A (in the generalized Coulomb gauge) in the standardvacuum representation, discuss the conditions for appearance of instabilities and their types (Subsec. III A), and present a concrete example (homogeneous medium; Subsec. III B) where calculations can be carried over to the end. In Sec. IV, we repeat the treatment of the previous section, but now for a more appealing scenario on the gravitational side: spherically-symmetric, stationary anisotropic media. Conditions for triggering instabilities and their types are shown to be very similar to those in the plane-symmetric case (Subsec. IV A). As a concrete application, in Subsec. IV B we show how to mimic QED-inspired nonminimally-coupled electromagnetism in the background spacetime of a Schwarzschild black hole. Then, Sec. V is dedicated to discuss possible stabilization mechanisms which might bear analogy to some curved-spacetime phenomena, such as spontaneous vectorization [18] and particle bursts due to tachyonic instability [28]. Finally, in Sec. VI we present some final remarks. We leave for an appendix tedious calculations related to the orthonormalization of modes of Sec. IV. We adopt the abstract-index notation to represent tensorial quantities (see, e.g., Ref. [29]) and, unless stated otherwise, we use natural units (in which ̵ h = c = 1).

II. COVARIANT ELECTROMAGNETISM IN ANISOTROPIC MATERIAL MEDIA
Electromagnetism in material media, in flat spacetime and in the absence of free charges, is described by two antisymmetric (observer-independent) tensors, F ab and G ab , satisfying the macroscopic covariant Maxwell's equations, where ∂ a is the derivative operator compatible with the flat metric η ab (but in arbitrary coordinates) and the square brackets denote antisymmetrization over the indices enclosed by them. These equations must be supplemented by medium-dependent constitutive relations between F ab and G ab , as well as initial and boundary conditions, in order to provide a well-posed problem. These constitutive relations are usually set at the level of (observer-dependent) fields E a , B a , D a , and H a , related to F ab and G ab through where u a is the four-velocity of the observer measuring these fields and abcd is the Levi-Civita pseudo-tensor (with 0123 = √ −η, η ∶= det(η µν )). Moreover, the constitutive relations usually take a simpler form in the reference frame in which the medium is (locally and instantaneously) at rest. Here, we consider a polarizable and magnetizable medium whose constitutive relations in its instantaneous rest frame take the form where the tensors ε ab and µ ab may depend on spacetime coordinates, and the system is assumed dispersionless. We return to this point later. The fact that Eqs. (7,8) are valid in the medium's instantaneous rest frame means that the fields E a , B a , D a , and H a appearing in them are related to F ab and G ab through Eqs. (3)(4)(5)(6) with u a = v a , the medium's four-velocity field. We proceed by splitting the "spatial" [30] tensors ε ab and µ ab into isotropic and traceless anisotropic parts, where h a b ∶= δ a b + v a v b is the projection operator orthogonal to v a . Inverting Eqs. (4,6) (with u a = v a ), and substituting Eqs. (7)(8)(9)(10) and (3,5), we obtain where we have defined the tensors and the squared refractive index n 2 = µε. The idea, then, is to consider the symmetric tensor g ab , defined through g ab g bc = δ c a , as an effective metric of a curved background spacetime perceived by the electromagnetic field F ab . Note that the components of g ab and η ab satisfy det g αβ = det η αβ (15) and, thus, One can easily check that g ab is explicitly given by Therefore, in an arbitrary coordinate system, Eq. (1) reads Up to this point, it was understood that the physical background metric η ab and its inverse η ab were responsible for lowering and raising tensorial indices. Now, with the introduction of an effective metric g ab , we should be careful when performing these isomorphisms. In order to minimize chances of confusion, we shall avoid lowering and raising tensorial indices using the effective metric, making explicit most appearances of g ab and g ab in the equations below, with few exceptions which will be clearly stated. One obvious exception is the definition of g ab as the inverse of g ab . Another such exception is the use of ∇ a to denote covariant derivative compatible with g ab . With this in mind, from Eqs. (2) and (17), the electromagnetic tensor F ab satisfies Notice that Eqs. (18) and (19) applied to homogeneous ab ) materials, with arbitrary 4-velocity field v a , lead to the same equations which rule minimally-coupled vacuum electromagnetism in a curved spacetime with metric µ n g ab . Optical analogue models in these configurations with µ = 1 were studied in [31,32]. Here, we shall focus on electromagnetism in anisotropic materials, more specifically, materials with only "shear-like" anisotropies: χ (ε) = 0 = χ (µ) [ab] . In this case, the tensor χ abcd defined in Eq. (14) has the same algebraic symmetries as the Riemann curvature tensor, namely, χ abcd = χ cdab and χ a[bcd] = 0 -in addition to χ abcd = χ [ab] [cd] , which is always true. The Eqs. (18) and (19) can be seen as analogous to some nonminimally-coupled electromagnetic field equations in curved spacetime. Although in general χ abcd is independent of the Riemann tensor associated with the effective metric g ab , one can construct cases where they are related. This is interesting because some one-loop QED corrections to Maxwell's field equations in curved spacetime [33,34] can be emulated by such nonminimal coupling, as we shall discuss below, in Subsec. II A.
Before considering particular applications of the equations above, let us define a sesquilinear form on the space of complexified solutions, which will be relevant when applying the canonical quantization procedure. As usual, let us solve Eq. (18) by introducing the 4-potential A a such that F ab = ∇ a A b − ∇ b A a . Then, let F ab and F ′ ab be two complex solutions of Eq. (19), associated to A a and A ′ a , respectively. With overbars representing complex conjugation, we contractĀ b (resp., A ′ b ) with Eq. (19) applied to F ′ cd (resp.,F cd ) and subtract one from the other, arriving at This continuity-like equation ensures that the quantity is independent of the space-like hypersurface Σ where the integration is performed -provided we restrict attention to solutions satisfying "appropriate" boundary condition -, where dΣ is the physical volume element on Σ and N a = η ab N b , with N a being a unit, future-pointing vector orthogonal to Σ (according to η ab ). More specifically, considering that the system of interest is contained in the spacetime region M ≅ T × Σ, where T ⊆ R is a real open interval, then the appropriate boundary condition amounts to imposing that the flux of the (sesquilinear) current appearing in Eq. (20) vanishes through T ×Σ (whereṠ denotes the boundary of the space S). In particular, in stationary situations which we shall treat here, this condition translates to where dS is the physical area element onΣ and s a is the unit vector field normal to T ×Σ (according to η ab ). Thus, these conditions being satisfied, Eq. (21) provides a legitimate sesquilinear form on the space S C of complexvalued solutions of Eqs. (18) and (19). Notice that for pure-gauge solutions -i.e., A a = ∇ a ψ, for some scalar function ψ -, (A, A) = 0. (The converse, however, is not true.) The relevance of this sesquilinear form is that it provides a legitimate inner product on a (non-unique choice of) subspace S + C ⊊ S C of "positive-norm solutions," which, together with its complex conjugate S − C ⊊ S C , generates all solutions: S C : S + C ⊕ S − C = S C . Loosely speaking, upon completion, S + C yields a Hilbert space H from which the (symmetrized) Fock space F s (H) is canonically constructed to represent states of the electromagnetic field. In particular, choosing S + C to be generated by positive-frequency solutions (those proportional to e −iωt , with ω > 0), the vacuum state of this Fock representation corresponds to the usual physical vacuum state of the field.
A. QED-inspired nonminimal couplings As mentioned earlier, Eqs. (18) and (19) can be interpreted as ruling electromagnetism in curved spacetimes with some QED-inspired nonminimal coupling χ abcd with the background geometry. In fact, in the one-loop-QED approximation [33,34], where α is the fine-structure constant, m e is the electron's mass, and R abcd , R ab , and R are, respectively, the Riemann, Ricci and Ricci-scalar curvature tensors associated with the (effective) metric g ab . By leaving α 1 , α 2 , α 3 unconstrained, Eq. (23) represents a three-parameter family of couplings of the electromagnetic field with the background effective geometry -see Ref. [35] for some interesting particular cases.
For a generic medium, χ abcd is not related to the geometry associated with g ab . However, we can simulate couplings given by Eq. (23) by conveniently relating n and v a (which determine g ab ) with µ and the anisotropic tensors χ ab (ε) and χ (µ) ab (which appear in χ abcd ). From Eqs. (14) and (23), and their contractions with g ab , g ac g bd χ abcd = 6 n µ − 1 = (α 1 + 3α 2 2 + 6α 3 )R, (25) we can solve for µ and the anisotropic tensors, obtaining: where V a = n 3 4 v a is the 4-velocity of the medium normalized according to the effective metric g ab and H ab ∶= g ab + V a V b . In Eqs. (27) and (28) indices are lowered and raised by the effective metric and its inverse. Notice that, ab v b , only geometries associated with g ab which can be put in the form given by Eq. (16) and satisfying for some timelike 4-vector v a , can be emulated by these anisotropic media -with v a then set as the medium's 4-velocity. Using Einstein's equations to map this constraint to the stress-energy-momentum tensor T ab of the corresponding gravitational source, we have that where, again, the effective metric and its inverse are used to lower and raise indices (and T ∶= T a a ). One can easily check that in case of perfect fluids -characterized by a proper energy density ρ and (isotropic) pressure p -, Eq. (30) is only satisfied for p = −ρ; i.e., for a cosmological-constant-type "fluid." However, if one allows for sources with anisotropic pressures (p 1 , p 2 , p 3 ), described by the stress-energy-momentum tensor -with {u a , e a 1 , e a 2 , e a 3 } being a tetrad and u a timelike -, then and (V a e a j ) (ρ + p j ) = 0, j = 1, 2, 3.
In particular, if V a = u a , then Eq. (32) is the only additional constraint to be enforced. Returning attention to the background effective geometry and recalling that all the geometric tensors are obtained from g ab given in Eq. (16), we see that Eq. (29) actually comprises a system of four differential equations which n and v a must satisfy. Electromagnetism with nonminimal coupling described by Eq. (23) can only be simulated in these anisotropic media if the background spacetime geometry is associated to solutions of this system [via Eq. (16)]. We shall treat a particular solution to these differential equations later.

III. PLANE-SYMMETRIC ANISOTROPIC MEDIUM AT REST
In this section, we consider the simplest case of an anisotropic medium: a plane-symmetric medium at rest in the inertial lab frame. The purpose of this section is not yet to establish an analogy with some interesting gravitational system, but to present the analysis in a simple context. In Sec. IV we apply the analysis to a more appealing scenario.
Let us consider a medium at rest in an inertial laboratory, such that in inertial Cartesian coordinates and with where convenience, we shall work in the generalized Coulomb gauge [36] in which . In this gauge, the t component of Eq. (19) is automatically satisfied, while the spatial components lead to First, let us consider solutions A such that A z = 0, which describe electric fields which are perpendicular to the z directiontransverse electric modes, A (TE) , for short [37]. In this case, our gauge condition ensures that there exists a scalar field ψ such that A with ζ being a spatial coordinate such that dζ = µ dz. The Eq. (42) must be supplemented by boundary conditions for f ωk . Imposing Eq. (22) to these modes leads where [] İ denotes the flux of the quantity in square brackets throughİ. This condition restricts the possible values of ω 2 . Let E (TE) k be the (k -dependent) set of ω values for which Eqs. (42) and (43) are satisfied for f we can orthonormalize these modes according to , where the sesquilinear form given in Eq. (21), applied to the current scenario, takes the form We obtain (up to a global phase) and n z ∶= (0, 0, 1). The second set of solutions of Eqs. (40) and (41), which describe magnetic fields which are perpendicular to the z directiontransverse magnetic modes, A (TM) , for short [37] -, is obtained by conveniently setting A where φ is an auxiliary function.
Our gauge condition then leads to A Using, again, staticity and planar symmetry, with ξ being a spatial coordinate such that dξ = ε dz. The boundary condition imposed by Eq. (22) now leads to k be the (k -dependent) set of ω values for which Eqs. (49) and (50) are satisfied for f k ∩R * + , we can normalize these modes according to Moreover, modes A ωk . The solutions expressed in Eqs. (47) and (53), dubbed positive-frequency normal modes, play a central role in the construction of the Fock (Hilbert) space of the quantized theory, as described at the end of the previous section. With these solutions, the quantum-field operator A is represented bŷ where "H.c." stands for "Hermitian conjugate" of the preceding term andâ As an application of our quantization scheme one can use the above formulas to obtain, for instance, the Carniglia-Mandel quantization [38] in a straightforward way. The system in this case is composed by a dielectricvacuum interface at z = 0 and a non-magnetizable (µ = µ = 1) homogeneous isotropic non-dispersive dielectric (ε = ε = ε ≡ n 2 ) filling the half-space z < 0. These data enter Eqs. (42) and (49), thus describing the background in terms of effective potentials of one-dimensional Schrödinger-like problems.

A. Instability analysis
In the analysis presented above, it was implicitly assumed that all constitutive functions ε , ε ∥ , µ , and µ ∥ are positive functions of z ∈ I. This condition ensures that the field modes presented in Eqs. (47) and (53), together with their complex conjugates, constitute a complete set of (complexified) solutions of Maxwell equations in R 3 × I; in other words, the boundary-value problems defined by Eqs. (42,43) and Eqs. (49,50) admit solutions only for (a subset of) ω 2 > 0. This is easily seen by interpreting them as null-eigenvalue problems for the linear operators defined in the square brackets of Eqs. (42) and (49). Experience with Schrödinger-like equations teaches us that these equations have solutions provided the associated effective potentials (terms in parentheses) become sufficiently negative in a given region -which implies ω 2 > 0 and, typically, the larger the k 2 , the larger the ω 2 .
Here, however, we shall consider a more interesting situation. It has been known for almost two decades that materials can be engineered so that some of their constitutive functions can assume negative values [39][40][41][42][43]. These exotic materials have been termed metamaterials. In this case, the effective potentials appearing in Eqs. (42) and (49) may become sufficiently negativegranting solutions to these boundary-value problemswithout demanding ω 2 > 0. For instance, if µ ∥ < 0 (with µ , ε > 0), then the larger the value of k , the more negatively it contributes to the effective potential of Eq. (42), favoring the appearance of solutions with smaller (possibly negative) values of ω 2 . The same is true for Eq. (49) if ε ∥ < 0 and similar analysis can be done if any other constitutive function becomes negative.
At this point, we must introduce an element of reality concerning the constitutive functions. We have been treating these quantities as given functions of z alone -neglecting dispersion effects, since we are, here, interested in gravity analogues. However, these material properties generally depend on characteristics of the electromagnetic field itself, particularly on its time variation (i.e., on ω), in which case Eqs. (7) and (8) would be valid mode by mode, with the constitutive tensors ε ab and µ ab possibly being different for different modes. When translated to spacetime-dependent quantities, Eqs. (7) and (8) would be substituted by sums over the set of allowed field modes [44]. Therefore, the precise key assumption about our metamaterial media is that some of their anisotropic constitutive functions ε , ε ∥ , µ , µ ∥ can become negative for some ω on the positive imaginary axis, ω 2 < 0. Notwithstanding, the less restrictive condition Im(ω) > 0 would suffice for our purposes. However, dealing with the case Im(ω)Re(ω) ≠ 0 would involve quantization in active media, which we shall treat elsewhere [44]. Moreover, our focus here is to show that the electromagnetic field itself can exhibit interesting behavior without need to exchange energy with the medium (which occurs in dispersive/active media). This justifies our focus on ω 2 < 0 in what follows. The possibility of having this type of material will be discussed later.
Let ω 2 = −Ω 2 (with Ω > 0) be such value for which at least one of the constitutive functions is negative for z ∈ I. Thus, both the effective potentials of Eqs. (42) and (49) take the general form with C 1 and C 2 being functions of z. Two interesting possibilities arise: • (i) C 1 < 0: In this case, the larger the value of k , the more negative the effective potential gets. Therefore, it is quite reasonable to expect that, for a given size of the interval I, one can always find "large enough" values of k -certainly satisfying k 2 > C 2 Ω 2 C 1 -such that the Schrödinger-like equation with effective potential V eff admits null-eingenvalue solutions. We shall refer to this situation as large-k instability; • (ii) C 1 > 0 and C 2 < 0: Under these conditions, the effective potential V eff , as a function of k , is bounded from below: V eff ≥ − C 2 Ω 2 . Therefore, a Schrödingerlike equation with effective potential V eff only admits null-eigenvalue solutions provided k is "sufficiently small" -certainly satisfying k 2 < C 2 Ω 2 C 1 -and the size of the interval where V eff is negative is "sufficiently large." We shall refer to this situation as minimum-width instability.
Let us call g (J) Ωk the null-eigenvalue solutions mentioned in either case above, with J ∈ {TE, TM} depending on whether it refers to Eq. (42) or (49) with ω 2 = −Ω 2 (without loss of generality, Ω > 0). These solutions are associated with unstable electromagnetic modes whose temporal behavior is proportional to e ±Ωt . Although it might be tempting not to consider these "runaway" solutions, [25,26], they are essential, if they exist, to expand an arbitrary initial field configuration satisfying the boundaryvalue problems set by Eqs. (42,43) and (49,50); in other words, the stationary modes alone do not constitute a complete set of solutions of Maxwell's equations with the given boundary conditions. And even if, on the classical level, one might want to restrict attention to initial field configurations which have no contribution coming from these unstable modes -which is certainly unnatural, for causality forbids the system to constrain its initial configuration based on its future behavior -, inevitable quantum fluctuations of these modes would grow, making them dominant some time e-foldings (t ∼ N Ω −1 , N ≫ 1) after the proper material conditions having been engineered. Therefore, these modes are as physical as the oscillatory ones. In fact, artificial inconsistencies have been reported in the literature, regarding field quantization in active media [25,26], which are completely cured when unstable modes are included in the analysis [44].
It is interesting to note that depending on which constitutive function is negative, Eqs. (42) and (49) may incur in different types of instabilities. For instance, if µ < 0 for a given ω 2 = −Ω 2 < 0, with all other constitutive functions being positive, then Eq. (42) exhibits case-(i) instability, while Eq. (49) incur in case-(ii) instability. This means that unstable TE modes -with some k > √ µ ∥ ε Ω -would certainly be present, while unstable TM modes -with some k < µ ε ∥ Ω -would only appear if the width of the material (size of the interval I) is larger than some critical value. We shall illustrate these facts in a simple example below. But first, let us analyze some features of these unstable modes. In order not to rely on particular initial field configurations, let us focus on the inevitable quantum fluctuations of these modes.

Unstable TE modes
Repeating the procedure which led us from Eq. (42) to Eq. (47) (and orthogonal to all other modes) read (up to a time translation) with 0 < κ < π, g Ωk normalized according to and s ε being the sign of the integral above. Calculating the electric E Ωk fields associated to these modes, we have: (and orthogonal to all other modes) read (up to a time translation) where, again, 0 < κ < π, g Ωk is normalized according to and s µ is the sign of the integral above. Calculating the electric E Ωk fields associated to these modes, we have: The modes given by Eqs. (61) and (67), if present, must be added to the expansion of the field operatorÂ given in Eq. (55), along with their complex conjugateswith corresponding annihilationâ (uJ) Ωk and creationâ (uJ) † Ωk operators, J ∈ {TE, TM}. The resulting operator expansion can then be used to calculate electromagnetic-field fluctuations and correlations. In the presence of unstable modes, it is easy to see that the field's vacuum fluctuations are eventually (t ≫ Ω −1 ) dominated by these exponentially-growing modes. Obviously, this instability cannot persist indefinitely as these wild fluctuations will affect the medium's properties, supposedly leading the whole system to a final stable state. In some gravitational contexts, stabilization occurs by decoherence of these growing vacuum fluctuations [13], giving rise to a nonzero classical field configuration -a phenomenon called spontaneous scalarization (for spin-0) [14][15][16][17] or vectorization (for spin-1 fields) [18]. It is possible that something similar might occur in the analogous system. We shall discuss this point further in Sec. V.

B. Example
Let us consider a very simple system just to illustrate the results above in a concrete scenario: a slab of width L (in the region −L 2 < z < L 2), made of a homogeneous material with, say, µ < 0 for a given ω 2 = −Ω 2 (Ω > 0) and all other constitutive functions positive. For concreteness sake, here we assume that this value ω 2 = −Ω 2 is isolated and that it is the most negative value of ω 2 for which µ < 0. This latter assumption is merely a matter of choice, while the former only affects the measure on the set of quantum numbers k : where L is the legth scale associated with the area of the "infinite" slab (L ≫ L).
According to the discussion presented earlier, in this scenario, TE modes incur in case-(i) (large-k ) instability, while TM modes undergo case-(ii) (minimum-width) instability. The solutions g (J) Ωk of Eqs. (42) and (49) with ω 2 = −Ω 2 are given by the normalizable -according to Eqs. (62) and (68) -solutions of the null-eigenvalue, Schrödinger-like equation with V eff being the well potential represented in Fig. 1. The depth of the potential is given by Although here we focus only on unstable modes, associated with g (J) Ωk , note that in this example there would also appear stationary bound solutions associated with f (TE) ω0k -if µ < 0 for some ω 0 ∈ R -, for some ∶= √ µ ∥ ε is the transverse refractive index for the TE modes. For such a hypothetical mode, the slab would act as a waveguide, keeping the mode confined due to total internal reflections at its boundaries. The only peculiar feature here is that k would assume arbitrarily large values (in practice, limited only by the inverse length scale below which the continuous-medium idealization breaks down) for a given ω 0 . Back to the unstable modes, a straightforward calculation leads to the familiar even and odd solutions to the square-well potential, with g while for the TM modes, 0 ≤ a m ≤ ΩL n ∥ 2 and The transverse momentum k is given in terms of a m by The explicit form of N (J) m is not particularly important, so we only present its asymptotic behavior for k → ∞ for the TE modes, and for k → 0 for both TE and TM modes, In Fig. 2 variable a -, and the functions − tan a and cot a (blue dashed lines and red dotted lines, respectively). Crossing of the blue dashed lines (respectively, red dotted lines) with a fixed solid black curve determines values a = a m for even (resp., odd) solutions g

(J)
Ωk , for the corresponding value of ΩL. The figure clearly corroborates our preliminary analysis, showing that unstable TE modes appear with arbitrarily large values of a m (and, therefore, of k ) and that unstable TM modes only appear if L is larger than some minimum width L 0 , given by The unstable TE and TM modes inside the slab can then be put in the form with (⌈x⌉ represents the smallest integer larger than, or equal to, x, while ⌊x⌋ represents the largest integer smaller than, or equal to, x). The corresponding electric and magnetic field modes are Let us recall that these modes give information about fluctuations and correlations of the electromagnetic field; as long as decoherence does not come into play, the expectation values of the field are null, ⟨Â⟩ = ⟨Ê⟩ = ⟨B⟩ = 0. We shall use these modes later, when discussing possible consequences of these analogue instabilities. But first, let us explore more interesting analogies.

IV. SPHERICALLY-SYMMETRIC, STATIONARY ANISOTROPIC MEDIUM
In the previous section, we presented with great amount of detail the canonical quantization scheme for the electromagnetic field in flat spacetime in the presence of arbitrary plane-symmetric anisotropic polarizable/magnetizable media at linear order. The vacuum of such system was then identified with the vacuum of some nonminimally-coupled spin-1 field in a true curved spacetime described by the effective metric g αβ = √ n diag(−n −2 , 1, 1, 1). The analysis had the advantage of generalizing in a unified language the quantization of various interesting models coming from quantum optics in terms of simple equations (e.g., the Carniglia-Mandel modes [38]). However, the analogue spacetime for these configurations is of mathematical interest only and does not capture the symmetry of physical spacetimes. In order to study more appealing analogues, in this section we turn to spherically symmetric configurations, presenting them in a more concise way -for the nuances of the quantization were already explained previously. In this context, we may obtain interesting analogues by also assuming that the medium is able to flow. If the refractive index in a flowing material is high enough, such that the velocity of light becomes smaller than the medium's velocity, then it is clear that a sort of event horizon will form (restricted only to some frequency band which may contain unstable modes). This kind of phenomenon enable us to study analogues of unstable black holes, for instance.
We start working in standard spherical coordinates (t, r, θ, ϕ), such that η µν = diag(−1, 1, r 2 , r 2 sin 2 θ). Let the medium's four-velocity field be v µ = γ (1, v, 0, 0), where v = v(r) and γ = (1 − v 2 ) −1 2 . The effective-metric components then take the form where the isotropic parts of the constitutive tensors (in the local, instantaneous rest frame of the medium) are functions of r -ε = ε(r), µ = µ(r) -and, as usual, n 2 = µε. As for the traceless anisotropic tensors χ ab (ε) and ab , their components read and Similarly to the plane-symmetric case, these anisotropic tensors simply mean that in the instantaneous local rest frame of the medium, its electric permitivity and magnetic permeability in the radial direction (ε ∥ and µ ∥ ) and in the angular directions (ε and µ ) satisfy the same relations given below Eqs. (37)(38)(39): Not surprisingly, the lab coordinates (t, r, θ, ϕ) are not the most convenient ones to express Eqs. (18) and (19) in the case of a moving medium. One might initially think that coordinates (τ, r, θ, ϕ) which diagonalize the components of the effective metric, obtained by defining τ ∶= t − p(r), with p(r) satisfying would lead to the simplest form of the field equations. In these coordinates, the effective line element ds 2 eff becomes where F = γ 2 (1 − n 2 v 2 ). It is noteworthy that for n = constant > 0 (such that the factors of n in ds 2 eff can be absorbed via τ ↦ n 3 4 τ and r ↦ n −1 4 r), then the line element above can be made to represent Schwarzschild spacetime by tuning v so that F ≡ (1 − r s r), where r s is some positive constant. This is achieved by a velocity field satisfying v 2 = 1 + (n 2 − 1)r r s −1 (n ≠ 1).
Despite this apparent simplification, the coordinate τ = t − p(r) with p satisfying Eq. (91) is not convenient to express Maxwell's equations in anisotropic media. This is due to the kinematic polarization (resp., magnetization) caused by the magnetic (resp., electric) field. In the case of small velocities and isotropic materials, this effect is modeled by Minkowski's equations [45]. The coordinates (τ, r, θ, ϕ) defined using Eq. (91) "diagonalizes" only the isotropic part of the theory and do not take into account the anisotropies. It turns out that a much better choice is obtained by setting τ ∶= t−p(r) and replacing condition given in Eq. (91) by where, again, n 2 ∥ ∶= µ ε . This choice fully decouples the electromagnetic field modes in the anisotropic, moving material medium, as we shall see below.
Introducing again the 4-potential A µ via F µν = ∂ µ A ν − ∂ ν A µ , in these new coordinates (τ, r, θ, ϕ), the convenient (generalized Coulomb) gauge conditions read A τ = 0 and where is merely an auxiliary variable such that dr d ≡ , ∂ is the derivative operator on the unit sphere compatible with its metric, and it is understood that r is a function of the auxiliary variable . In this gauge, Maxwell's equations lead to where ρ appearing in Eq. (96) is another auxiliary variable defined through dr dρ ≡ γ 2 (1 − n 2 ∥ v 2 ) µ and ∆ S are the Laplacian operators defined on the unit sphere, acting on scalar and covector fields, respectively.
In order to solve these equations, we proceed in close analogy to the plane-symmetric case. First, let us find solutions with A r = 0 -the transverse electric modes, A (TE) . The gauge conditions imply that these solutions can be written as A (TE) = (0, ∂ ϕ ψ sin θ, − sin θ∂ θ ψ), where ψ is an auxiliary function to be determined. Making use of the stationarity and spherical symmetry of the present scenario, we can look for field modes of the form ψ = e −iωτ Y m (θ, ϕ)f where it is understood that r is a function of the auxiliary variable ρ. Notice the similarity between this equation and Eq. (42). In fact, the boundary condition given by Eq. (22) assumes the same form here as it does in the plane-symmetric case: This boundary condition ensures that these modes can be orthonormalized according to the sesquilinear form given in Eq. (21), which in this spherically-symmetric scenario assumes the form with Σ t being a spacelike surface t = constant. After some tedious but straightforward manipulations (presented in the appendix), we obtain the final form of normalized, positive-frequency TE modes: with f (TE) ω satisfying Eqs. (97) and (98), and normalized according to Note that the integration variable is [instead of ρ appearing in Eq. (97)] and I stands for the domain of integration in this variable corresponding to I in coordinate r. Now, let us look for solutions with A r ≢ 0the transverse magnetic modes, A (TM) . Let φ be such that ∆ Using again stationarity and spherical symmetry, Notice, again, the similarity between this equation and Eq. (49). And, again, the boundary condition imposed by Eq. (22) to these modes take the same form as in the plane-symmetric case: Properly orthonormalizing these modes using Eq. (A1) -see appendix -, leads to the positive-frequency TM normal modes with f (TM) ω satisfying Eqs. (102) and (103), and normalized according to Similarly to the TE case, note that the integration variable is not the same which appears in the differential equation, Eq.
requires that the canonical commutation relations hold.

A. Instability analysis
The close similarity between Eqs. (42) and (97) and between Eqs. (49) and (102) make the instability analysis in this spherically-symmetric scenario essentially identical to the one performed in the plane-symmetric case, with ( + 1) playing the role k 2 did in Eq. (58). So, putting the effective potentials of Eqs. (97) and (102), with ω 2 = −Ω 2 , in the form we again have two types of instabilities: (i) large-instability, when C 1 < 0 somewhere, and (ii) minimumthickness instability, when C 1 > 0 but C 2 < 0 in a sufficiently thick spherical shell -see discussion below Eq. (58). The only additional feature is that, by allowing the medium to flow, type-(i) (large-) instability for both TE and TM modes can arise when the medium's velocity v(r) exceeds the radial light velocity n −1 ∥ .
e Ωτ −is ε κ 2 + e −Ωτ +is ε κ 2 2 Ω ( + 1) sin κ g with κ being a constant (0 < κ < π), g (TE) Ω normalized according to and s ε being the sign of the integral above. Calculating the electric E vector fields associated to these modes in the lab frame, we have: (115)

Unstable TM modes
Finally, the unstable TM modes orthonormalized according to the analogous of Eqs. (107) and (108) read (up to global phase and time translation) 2 Ω 3 ( + 1) sin κ g vector fields associated to these modes in the lab frame, we have: As argued in the previous case, when instability is triggered and modes A (uJ) Ω m appear, they must be included in the field expansion given by Eq. (106), along with their complex conjugates. Eventually (t ≫ Ω −1 ), these modes dominate the field fluctuations.

B. Example
Now, let us consider a concrete scenario where electromagnetism in a gravitationally interesting system, nonminimally coupled to the background geometry via χ abcd given by Eq. (23) (but with arbitrary α 1 , α 2 , α 3 ), can be mimicked by an anisotropic, stationary moving medium. We have already seen that setting n = constant and v 2 = [1 + (n 2 − 1)r r s ] −1 , leads to an effective line element which describes the vacuum Schwarzschild spacetime. In this case, Eq. (29) is trivially satisfied and Eqs. (26)(27)(28) give which lead to the material properties We promptly see that n ∥ ∶= √ µ ε = n, which shows that the analogue horizon for these nonminimally-coupled modes, located where v 2 = n −2 ∥ , coincides with the analogue Schwarzschild radius r s . [Note, however, that this system is analogous to a physical black hole with Schwarzschild radius R s = n 1 4 r s , due to absorption of √ n in Eq. (92).] As for the other refractive indices, Fig. 3 shows their squared values (in black and red, respectively) for positive (solid lines) and negative (dashed lines) values of α 1 . Note that, depending on the values of α 1 (n 1 2 r 2 s ), some kind of metamaterial (possibly with some negative squared refractive indices) may be needed in order to mimic this nonminimal coupling of the electromagnetic field with the Riemann curvature tensor in the exterior region of a Schwarzschild black hole. Conversely, regardless how difficult it may be to set up such an experimental configuration in the lab, it is interesting in its own that QED-inspired nonminimally-coupled electromagnetism in the background of a black hole behaves as in such an exotic metamaterial in flat spacetime.
Turning to the question of possible instabilities, in Fig. 4 we show the behavior of the terms C 1 and C 2 appearing in Eq. (111) for the TE (in blue) and TM (in red) modes -extracted, respectively, from Eqs. (97) and (102): where the first and second lines in the expressions above refer to the TE and TM modes, respectively. The Fig. 4(a) is representative of the behavior of C 1 for −r 2 s √ n 2 < α 1 < r 2 s √ n, while Fig. 4(b) gives the correct qualitative behavior of C 1 for α 1 < −r 2 s √ n 2 or It is clear, from the expressions above, that C 2 is everywhere non-negative, while C 1 assumes negative values in the region with radial coordinate r between (α 1 r s √ n) 1 3 and r s (if α 1 > 0) or between [ α 1 r s (2 √ n)] 1 3 and r s (if α 1 < 0). Therefore, according to the discussion of Subsec. IV A, this nonminimally-coupled electromagnetic theory in Schwarzschild spacetime exhibits large-instability. In particular, if α 1 > r 2 s √ n or α 1 < −2r 2 s √ n, then the unstable modes influence the exterior region of the back hole. We now turn our attention to discussing what can possibly happen to the analogous system when the vacuum instability is triggered. In the gravitational scenario, it has been shown that in some cases (for instance, According to the instability discussion, only large-instability can appear in this case, since C2 ≥ 0 everywhere. Moreover, for α1 < −r 2 s √ n 2 or α1 > r 2 s √ n, the unstable modes can be mostly supported outside the analogous event horizon, r > rs. depending on the field-background coupling), stabilization occurs due to the appearance of a nonzero value for the field (spontaneous scalarization/vectorization) [14][15][16][17][18], seeded by decoherence of the growing initial-vacuum fluctuations [13]. In this process, field particles/waves are produced [14,28] and carry away the energy excess of the initial vacuum state in comparison to the stabilized configuration.
If we transpose these conclusions, mutatis mutandis, to our analogous systems, then an electromagnetic field should spontaneously appear in the material, bringing the whole system to a new equilibrium configurationthrough nonlinear effects brought in by field-dependent constitutive tensors ε ab and µ ab [see Eqs. (7,8)] -, with photons being emitted, carrying away the energy excess. Although the detailed dynamics of the stabilization processes in the gravitational and in the analogous sys-tems are quite different -ruled by Einstein equations in the gravitational case and by the macroscopic Maxwell's equations with field-dependent ε ab and µ ab in the analogous systems -, the qualitative features of the whole process, described above, seem quite reasonable to occur in generic field stabilization processes.
It is important to mention that the time scale set by the instability, Ω −1 , is typically of the order of the time light takes to travel the typical size of the system, L. Therefore, in the analogous lab scenarios, the stabilization process would occur almost instantaneously (∼ L (1 cm) × 10 −10 s) once the instability conditions are met -which, for a given system, may depend on external parameters such as temperature, external fields, etc., through their influence on the constitutive functions ε , ε ∥ , µ , µ ∥ . The whole process would most likely be interpreted as a kind of phase transition, where the "long-range" emergent correlations in the material would come from interaction of its constituents with a common (initially-unstable vacuum) fluctuating mode and/or the stabilized field configuration.
For concreteness sake, let us consider the explict form of the unstable modes found in the example of Sec. III, where instability occurs due to a negative value of µ -for some (isolated) ω 2 = −Ω 2 < 0 -in a homogeneous slab of width L. Although this system is not analogous to vacuum nonminimally-coupled electromagnetism in any realistic spacetime, it serves to illustrate general features of the mechanism itself, in addition to being much simpler to setup in the lab. This is no different than looking for fingerprints of analogue Hawking radiation in systems whose only similarity with realistic black holes is the pres-ence of an effective event horizon -which is the common approach in condensed-matter and optical experimental analogues.
As argued before, once instability sets in, the unstable modes must be added to the expansion of the field oper-atorÂ, along with their complex conjugates, with corresponding annihilationâ (uJ) Ωk and creationâ (uJ) † Ωk operators. It is easy to see that the field's vacuum fluctuations and correlations are eventually (t, t ′ ≫ Ω −1 ) dominated by these unstable modes -at least as long as decoherence does not come into play. The dominant contribution to the vacuum correlations in the example of Subsec. III B reads (the reader should refer to Subsec. III B for the definition of all quantities appearing in these expressions): terials with appropriate constitutive functions. This follows from the formal analogy between electromagnetism in anisotropic media and nonminimally-coupled electromagnetism in curved spacetimes, presented in Sec. II. We explored two concrete scenarios: (i) a plane-symmetric, static slab -whose main interest is its simplicity regarding experimental setup (see Sec. III) -and (ii) a spherically-symmetric, moving media -whose main feature is its analogy with QED-inspired nonminimally-coupled electromagnetism in Schwarzschild spacetime [33,34] for given velocity and constitutive-functions profiles (see Subsec. II A and Sec. IV).
Once instability is triggered in the analogous systems, some stabilization process must take place, leading the system to a new stable configuration. The details of this stabilization process and of the final configuration will most likely depend on specific nonlinear properties of the metamaterial, but it seems reasonable that they might involve the appearance of nonzero electromagnetic fields in the material (analogous to spontaneous vectorization in curved spacetimes) and photo production which carries away the energy excess with respect to the stable configuration. As discussed earlier, the time scale involved in the stabilization process can be very short (∼ 10 −10 s), which would make it very difficult to even identify the unstable phase. This is similar to what might occur with negative conductivity, which has never been directly measured but which is predicted to lead to zero-dc-resistance states [46] which were observed in laboratory [47,48] -although an alternative explanation has been proposed [49].
Clearly, the feasibility of such analogues is bound to the existence of material configurations with the required constitutive functions. As briefly pointed out in the introduction, this can be achieved at least for anisotropic neutral plasmas, and the recent advances in metamate-rial science offer a plethora of possible candidates, specially the hyperbolic metamaterials [42,43], that possess precisely the form given in Eqs. (9,10) with the required "negativeness." In particular, we call attention to the increase in the spontaneous light emission in such configurations, which may be related to the process of stabilization in active scenarios.
It is also important to mention that the QED-inspired analogues (Subsec. II A) are not restricted to the study of vacuum instability. For instance, they can be used to study light ray propagation in the corresponding spacetimes and one possible application is the QED-induced birefringence in the Schwarzschild spacetime [33]. For this particular experiment, one can work far from the effective horizon, where the constitutive coefficients (123)-(126) are positive.
Our main purpose here was to lay down a novel class of analogue models of curved-spacetime phenomena, with main interest on the gravitational side of the analogy. Notwithstanding, the consequences of the analogue gravity-induced instability to the metamaterial side may be interesting on its own. The electromagnetic field instability may mark, lead or mediate some kind of phase transition in the metamaterial, where the spontaneously created field and/or its amplified "long-range" correlations may play some important role (see discussion in Sec. V). Investigation in these lines are currently in course and will be presented elsewhere.