Evolving black hole horizons in General Relativity and alternative gravity

From the microscopic point of view, realistic black holes are time-dependent and the teleological concept of event horizon fails. At present, the apparent or the trapping horizon seem its best replacements in various areas of black hole physics. We discuss the known phenomenology of apparent and trapping horizons for analytical solutions of General Relativity and alternative theories of gravity. These specific examples (we focus on spherically symmetric inhomogeneities in a background cosmological spacetime) may be useful as toy models for research on various aspects of black hole physics.


Introduction
In the literature on gravitation and quantum field theory in curved space, several kinds of horizon are studied, including Rindler horizons for accelerated observers in Minkowski space, black hole horizons, and cosmological horizons. Research on classical and semiclassical black hole physics has unveiled inner, outer, Cauchy, and extremal horizons. The early literature on black holes and the 1970s development of black hole thermodynamics focused on stationary (or even static) black holes and on event horizons (see, e.g., Refs. [1,2,3]). However, highly dynamical situations are of great interest for theorists, including gravitational collapse, the merger of a black hole with a compact object, evaporation of a small black hole due to Hawking radiation, and black holes interacting with non-trivial environments. Examples of nontrivial environments occur in the case of black holes accreting or emitting gravitating matter such as Vaidya spacetimes (which recur in mathematical studies of horizon dynamics); black holes immersed in a cosmological background other than de Sitter space; black holes emitting (and possibly accreting) Hawking radiation in the final stages of their evolution when backreaction is significant; or black holes with variable mass due to other conceivable processes. In all these situations the concept of event horizon must be generalized. Moreover, if the black hole is located in a non-Minkowskian background, its mass-energy (which is usually the internal energy appearing in the first law of black hole thermodynamics) needs to be defined carefully through some notion of quasi-local energy, and is related to the notion of horizon.
In Rindler's words, an horizon is "a frontier between things observable and things unobservable" [4]. The horizon concept, which is the product of strong gravity, is perhaps the most impressive feature of a black hole spacetime and, traditionally, the one that best characterizes the black hole concept itself. Various notions of black hole horizon studied in the technical literature include event, Killing, inner, outer, Cauchy, apparent, trapping, quasi-local, isolated, dynamical, and slowly evolving horizons (for reviews see Refs. [1,3,5,6,7,8]). The notions of event, apparent, trapping, and dynamical horizon usually coincide for stationary black holes but they are quite different from each other in the case of dynamical black holes with masses evolving in time.
The usual definition of black hole event horizon turns out to be pretty much useless for practical purposes in highly dynamical situations because it requires the knowledge of the entire causal structure of spacetime (including future null infinity I + ), which is physically impossible [9,10,6]. Time-evolving black holes are not just a theoretician's playground but they are also important for astrophysics. The remarkable improvements in astronomical techniques in recent years and their projected developments in the near future have stressed the important roles that stellar mass and supermassive black holes play in the modelling of astrophysical systems. The improvement of ground-based gravitational wave detectors (presumably nearing their first detection in the next few years) and the development of space-based detectors prompt enormous theoretical efforts to predict in detail the gravitational waves emitted by astrophysical black holes and build template banks for interferometric detectors. Progress from attempts (see [25,26,27] for reviews) to explain the current acceleration of the universe discovered with type Ia supernovae without invoking an ad hoc dark energy [28]. One should also mention attempts to replace the concept of dark matter at galactic and cluster scales by modifying not just relativistic, but even Newtonian gravity, given that dark matter particles seem to elude direct detection (there is yet no agreement from the various groups on reports of possible signals).
From the astrophysical point of view, if primordial black holes formed in the early universe, they would have had a scale comparable to the Hubble scale and very dynamical horizons evolving on the Hubble time scale, and an important question is how fast these black holes could accrete and grow.
To summarize, there are currently many avenues of research in theoretical gravitational physics in which evolving horizons play some role. Here we review the main definitions and properties of horizons and we focus on apparent and trapping black hole horizons. While it is questionable that apparent horizons are the "correct" notion of horizon in the dynamical case (and there are indications that they may not satisfy a quantum generalized second law [29]), they are the ones that are used in practice in numerical relativity and there seems to be no better candidate at the moment for the concept of "horizon" when time-dependence and interaction with the environment are allowed. Since only a few exact solutions of GR (and even less of other theories of gravity) are known for which the horizons are explicitly time-dependent, here we focus on solutions describing black holes embedded in cosmological backgrounds, which have been studied in some detail. Fig. 1 shows how a conformal diagram of an hypothetical cosmological black hole may look like.
There are many additional reasons for pursuing the study of analytical solutions of GR and of gravitational theories representing a central inhomogeneity in a cosmological background. The non-linearity of the field equations prevents the splitting of solutions into a "background" and a "deviation" from it in general, but we drop the quotation marks and from now on the term "background" refers to the asymptotic structure of spacetime. We have already mentioned the use of alternative theories of gravity, and of f (R c c ) gravity in particular (where R c c is the Ricci curvature of spacetime and f is the Lagrangian density of the gravitational field), to explain the present acceleration of the universe without dark energy: since these theories are designed to produce a time-varying effective cosmological "constant", black hole spacetimes in these theories are naturally asymptotically FLRW, not asymptotically flat, and are dynamical. It is also of interest to study the spatial variation of fundamental constants throughout the universe, and scalar-tensor theories of gravity [30,31,32,33] embody the variation of the gravitational "constant", hence the search for analytical solutions describing condensations in cosmological spacetimes in these theories. Overall, very few such solutions are known in alternative gravity. But then one realizes that also analytical solutions of GR interpretable as central objects in cosmological backgrounds are quite interesting. The first solution of this kind is the McVittie spacetime [34], which was invented to address the problem of whether, or how, the expansion of the universe affects local systems (see Ref. [35] for a recent review on is a spacelike black hole singularity, and an apparent horizon (marked AH) can change from timelike, to null, to spacelike, and it can be located inside or outside the event horizon (fortyfive degrees line marked EH) according to whether the energy conditions are satisfied or not. If R ab l a l b ≥ 0 for all null vectors l a , then the apparent horizon lies inside the event horizon ( [3], p. 311).
this subject). The McVittie solution of GR has a complex structure and is not yet completely understood [36,37,38,39,40,41,161]. Relatively few other solutions of GR with similar features have been discovered, including Swiss-cheese and other models [44].
Recent interest in cosmological condensations in GR arises also from a different attempt to explain the present cosmic acceleration without dark energy and without modifying gravity. The idea is that the backreaction of inhomogeneities on the expansion of the universe could be sufficient to produce the observed acceleration [45,46,47,48,49,50,51,52,53,54,55]. However, the implementation of this idea has several formal problems and its proponents have not yet shown convincingly that this idea explains the magnitude or even the sign of the cosmic acceleration [56,57] (more mathematically oriented work puts this proposed solution to the cosmic acceleration problem in jeopardy [58]). The study of exact inhomogeneous universes has also been pursued in yet another attempt to explain the current acceleration of the universe, the dominant idea being that we live inside a giant void which mimics an accelerated expansion; some of the analytical GR solutions considered are related to black holes in expanding universes (see [59] for a review).
Independent motivation for the study of evolving horizons, and one not insignificant for astrophysics, comes from the renewed interest in exact models of spherical accretion by black holes, in particular the accretion of dark or phantom energy [60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78]. This issue may again be relevant for primordial black holes which need to grow fast if they are to survive until the present era.
The plan of this review paper is the following: first, we review basic material in the next section. The following section discusses analytical solutions of GR, and is followed by a section on spacetimes with similar features in other theories of gravity. Due to space limitations it is not possible to discuss all the known solutions and their details, or to review all the works on cosmological black holes, but we do provide more detail for a few solutions to illustrate the techniques used, and the selection made is no doubt biased. Moreover, we do not discuss here the more mathematical approaches to time-evolving horizons and the various existence and uniqueness theorems for horizons. The metric signature used is − + ++ and we follow the conventions of Ref. [3] (the speed of light c and Newton's constant G are set to unity except where, occasionally, we restore them explicitly).

Various notions of horizon
Let us review briefly the geometry of the congruences of null geodesics crossing a horizon, which are used in the definition of non-stationary horizons.

Null geodesic congruences and trapped surfaces
Consider a congruence of null geodesics with tangent l a = dx a /dλ, where λ is an affine parameter along each geodesic and l a l a = 0, l c ∇ c l a = 0. The metric h ab in the 2-space orthogonal to l a is determined by the following [1]: pick another null vector field n a such that n c n c = 0 and l c n c = −1, then we have h ab ≡ g ab + l a n b + l b n a .
h ab is purely spatial and h a b is a projection operator on the 2-space orthogonal to l a , i.e., The choice of n a is not unique but the geometric quantities of interest to us do not depend on it once l a is fixed. Let η a be the geodesic deviation vector (it corresponds to a one parameter subfamily of the congruence since its choice is not unique) and define the tensor field [3,1] and is orthogonal to the null geodesics, B ab l a = B ab l b = 0. The transverse part of the deviation vector is and the orthogonal component of l c ∇ c η a , denoted by a tilde, is [1] ( The transverse tensor B ab is decomposed into its symmetric and antisymmetric parts, and the symmetric part is further decomposed into its trace and trace-free parts as [3,1] is the expansion of the affinely parametrized congruence, is the expansion tensor, is the shear tensor, and ω ab =B [ab] (2.9) is the vorticity tensor. The expansion, shear, and vorticity tensors are purely transversal (i.e., orthogonal to l a ) and the shear and vorticity are trace-free. The shear scalar and vorticity scalar σ 2 ≡ σ ab σ ab , ω 2 ≡ ω ab ω ab (2.10) are non-negative. The expansion propagates along a null geodesic according to the celebrated Raychaudhuri equation [3,1], which was the main tool in the proof of the singularity theorems of Hawking and Penrose [79,3], similar propagation equations hold for σ ab and ω ab [3]. If the congruence of null geodesics with tangent l a is not affinely parametrized, the geodesic equation assumes the form where the quantity κ which measures the failure of l a to be affinely parametrized is sometimes used, on a horizon, as a possible definition of surface gravity [80,82] (there are various inequivalent definitions of surface gravity in the literature). The expansion is now (2.14) Eq. (2.14) is independent of the field equations of the theory and can be applied when l c and n c are not normalized to satisfy l c n c = −1 as usual [80]. With non-affine parametrization, the Raychaudhuri equation picks up an extra term [1], A compact and orientable 2-surface embedded in 4-space has two independent directions orthogonal to it, corresponding to ingoing and outgoing null rays. One is naturally led to study congruences of ingoing and outgoing null geodesics with tangent fields l a and n a , respectively, and the way they propagate in strong gravity. Let us provide now some basic definitions for closed 2-surfaces [10,5,6,83,8] (usually it is assumed that these 2-surfaces are spacelike [10, 84, 5] but we will not impose this requirement here): • A normal surface corresponds to θ l > 0 and θ n < 0 (for example, a 2-sphere in Minkowski space satisfies this property).
• A trapped surface [85] corresponds to θ l < 0 and θ n < 0. The outgoing, in addition to the ingoing, future-directed null rays converge here instead of diverging and outwardpropagating light is dragged back by strong gravity.
• A marginally outer trapped (or marginal) surface (MOTS) corresponds to θ l = 0 (where l a is the outgoing null normal to the surface) and θ n < 0.
• An untrapped surface is one with θ l θ n < 0.
• An antitrapped surface corresponds to θ l > 0 and θ n > 0 (both outgoing and ingoing future-directed null rays are diverging).
• A marginally outer trapped tube (MOTT) is a 3-dimensional surface which can be foliated entirely by marginally outer trapped (2-dimensional) surfaces.
In GR Penrose has proved that, if a spacetime contains a trapped surface, the null energy condition holds, and there is a non-compact Cauchy surface for the spacetime, then this spacetime contains a singularity [85]. Trapped surfaces are probably essential features in the concept of black hole and notions of "horizon" of practical utility will be identified with boundaries of spacetime regions containing trapped surfaces. At present, the mathematical conditions for the existence and uniqueness of MOTSs are not completely clear. It is known that, in general, a MOTT can be distorted smoothly, hence MOTTs are non-unique [86,87,5].
Let us review the various kinds of horizons appearing in the literature on black holes, cosmology, quantum field theory in curved spaces, and the corresponding thermodynamics.

Event horizons
The event horizon is the traditional notion of horizon for stationary black holes in GR. An event horizon is a connected component of the boundary ∂ (J − (I + )) of the causal past J − (I + ) of future null infinity I + [88,79,3,1]. This is the most peculiar feature of a black hole: the horizon is a causal boundary separating a region from which nothing can come out to reach a distant observer from a region in which signals can be sent out and eventually arrive to this observer. An event horizon is generated by the null geodesics which fail to reach future null infinity and, therefore (provided that it is smooth) is always a null hypersurface.
In astrophysics the concept of event horizon is implicitly taken as a synonym of black hole. However, since to define and locate an event horizon one must know all the future history of spacetime (one must know all the geodesics which do reach future null infinity and, tracing them back, the boundary of the region from which they originate), an event horizon is a globally defined concept. For an observer to state that a black hole event horizon has formed requires knowledge of the spacetime outside his or her future light cone, which is impossible to achieve unless the spacetime is stationary, the black hole has existed forever, and nothing changes (a common expression is that the event horizon has a teleological nature). It has been shown [7,84] that, in a collapsing Vaidya spacetime, an event horizon forms and grows starting from the centre and an observer can cross it and be unaware of it even though his or her causal past consists entirely of a portion of Minkowski space: the event horizon cannot be detected by this observer with a physical experiment. In other words, the event horizon "knows" about events belonging to a spacetime region very far away and in its future but not causally connected to it ("clarvoyance" [5,89,90]).
Due to its global nature, an event horizon is not a practical concept and it is nearly impossible to locate precisely an event horizon in a dynamical situation. Realistic astrophysical black holes have not existed forever but are formed by gravitational collapse. In numerical relativity, codes designed to follow a collapse situation, a binary system merger, or other dynamical situations generating black holes, eventually crash and it is impossible to follow the evolution of a system to future null infinity. It is routine in numerical relativity to employ marginally trapped surfaces as proxies for event horizons (e.g., [12,13,14]).
The event horizon H is a tube in spacetime; a very common abuse of terminology consists of referring to the intersections of H with surfaces of constant time (which produce 2-surfaces) as "event horizons" (this abuse of terminology extends to all the notions of horizon that we define below).

Killing horizons
When present, a Killing vector field k a satisfying the Killing equation ∇ a k b + ∇ b k a = 0 defines a Killing horizon H of the spacetime (M, g ab ), which is a null hypersurface which is everywhere tangent to a Killing vector field k a which becomes null, k c k c = 0, on H. This Killing vector field is timelike, k c k c < 0, in a spacetime region which has H as boundary. Stationary event horizons in GR are Killing horizons [91], for example in the Schwarzschild geometry the event horizon R = 2M is also a Killing horizon and the timelike Killing vector k a = (∂/∂t) a in the R > 2M region outside the event horizon becomes null at R = 2M and spacelike for R < 2M . An event horizon in a locally static spacetime is also a Killing horizon for the Killing vector k a = (∂/∂t) a associated with the time symmetry [1]. If the spacetime is stationary and asymptotically flat (but not necessarily static), it must be axisymmetric and an event horizon is a Killing horizon for the Killing vector which is a linear combination of the vectors associated with time and rotational symmetries, and where Ω H is the angular velocity at the horizon (this statement requires the assumption that the Einstein-Maxwell equations hold and some assumption on the matter stress-energy tensor [79,92]). When present, a Killing horizon defines a notion of surface gravity κ Killing , as we will see below.
Of course, the concept of Killing horizon is useless in spacetimes which do not admit timelike Killing vectors. Attempts to use conformal Killing horizons in spacetimes conformal to the Schwarzschild one ( [93,94], see also [95,96,98,99,100]) have not been fruitful. However, the introduction of the Kodama vector, which is defined in spacetimes without Killing vectors, in place of a Killing field is much more useful to introduce a surface gravity.

Apparent horizons
A future apparent horizon is the closure of a surface (usually a 3-surface) which is foliated by marginal surfaces; it is defined by the conditions on the time slicings [22] θ l = 0 , (2.17) where θ l and θ n are the expansions of the future-directed outgoing and ingoing null geodesic congruences, respectively (this more practical definition differs from that of Hawking and Ellis [79], which is rather impractical [5]). Eq. (2.17) expresses the fact that the congruence of future-pointing outgoing null rays momentarily stops expanding and, presumably, these rays turn around at the horizon, while the condition (2.18) originally served the purpose of distinguishing between black holes and white holes. Apparent horizons are defined quasi-locally and are independent of the global causal structure of spacetime, contrary to event horizons. However, apparent horizons (and also trapping horizons, see below) depend on the choice of the foliation of the 3-surface with marginal surfaces [101,102] and, of course, also the ingoing and outgoing null geodesics orthogonal to these surfaces do, as well as their expansions θ l and θ n [103]. While the expansions are scalars, and are therefore independent of the coordinate system chosen, sometimes a coordinate choice makes it easier to specify locally the foliation (for example by choosing spacelike surfaces of constant time coordinate -different time coordinates identify different families of hypersurfaces of constant time), which is a geometric object and is coordinate-independent. Congruences of outgoing and ingoing null geodesics orthogonal to these surfaces will, of course, change when changing the foliation. The dependence of apparent horizons on the spacetime slicing is illustrated by the fact that non-symmetric slicings of the Schwarzschild spacetime can be found for which no apparent horizons exist [101,102].
Apparent horizons are, in general, quite distinct from event horizons: for example, event and apparent horizons do not coincide in the Reissner-Nordström black hole (inner horizon) and in the Vaidya spacetime [1]. Also in static black holes which are perturbed, the apparent and the event horizons do not coincide [136]. During the spherical collapse of uncharged matter an event horizon forms before the apparent horizon does and the two come closer and closer until they eventually coincide asymptotically as the final static state is reached [79].
In GR, a black hole apparent horizon lies inside the event horizon provided that the null curvature condition R ab l a l b ≥ 0 ∀ null vector l a is satisfied [79]. This requirement coincides with the null energy condition T ab l a l b ≥ 0 ∀ null vector l a if the Einstein equations are imposed, and in this case it is believed to be a reasonable condition on classical matter. However, Hawking radiation itself violates the weak and the null energy conditions [105], as does quantum matter. A simple scalar field non-minimally coupled to the curvature can also violate all of the energy conditions. The null curvature condition is easily violated also in alternative theories of gravity (for example, Brans-Dicke [30] and scalar-tensor [31,32,33] gravity) and the black hole apparent horizon has been observed to lie outside of the event horizon during spherical collapse in Brans-Dicke theory, although it eventually settles inside of it when the static Schwarzschild state is achieved [106] (note that the GR black hole is the endpoint of collapse in general scalartensor gravity for asymptotically flat black holes without matter other than the Brans-Dicke scalar field outside the horizon [107,108]).
To summarize, the cherished notion of event horizon seems rather useless in general dynamical situations and apparent horizons appear to be more practical in spite of their fundamental limitations of depending on the spacetime slicing and of possibly being timelike surfaces (this last drawback is probably the most puzzling one [83]).

Trapping horizons
A future outer trapping horizon (FOTH) is the closure of a surface (usually a 3-surface) foliated by marginal surfaces such that on its 2-dimensional "time slicings" ([22], see also [109] and references therein) where θ l and θ n are the expansions of the future-directed outgoing and ingoing null geodesic congruences, respectively. The inequality (2.21) serves the purpose of distinguishing between inner and outer horizons, e.g., in the non-extremal Reissner-Nordström solution, and also distinguishes between apparent horizons and trapping horizons (it is not imposed for apparent horizons but it is required for trapping ones), and its sign distinguishes between future and past horizons. The definition of a past inner trapping horizon (PITH) is obtained by exchanging l a with n a and reversing the signs of the inequalities, The past inner trapping horizon identifies a white hole or a cosmological horizon. As one moves just inside an outer trapping horizon, one encounters trapped surfaces, while trapped surfaces are encountered as as one moves just outside an inner trapping horizon.
As an example, consider the static Reissner-Nordström black hole with the natural spherically symmetric foliation: the event horizon r = r + is a future outer trapping horizon (FOTH), the inner (Cauchy) horizon r = r − is a future inner trapping horizon (FITH), while the white hole horizons are past trapping horizons (PTHs).
Black hole trapping horizons have been associated with thermodynamics, and it is claimed that it is the trapping horizon area and not the area of the event horizon which should be associated with entropy in black hole thermodynamics [110,111,21,6]. This claim is controversial [112,113,114]. The Parikh-Wilczek "tunneling" approach [115] is in principle applicable also to apparent and trapping horizons, not only to event horizons [116,117,118,119,120,121,122,123] but also this aspect is not entirely free of controversy [124].
In general, trapping horizons do not coincide with event horizons. Dramatic examples are spacetimes which possess trapping horizons but not event horizons [125,126]. The difference between the areas of the trapping and the event horizon in particular spacetimes have been studied in Ref. [127].

Isolated, dynamical, and slowly evolving horizons
Isolated horizons correspond to isolated systems in thermal equilibrium not interacting with their surroundings, which are described by a stress-energy tensor T ab . The concept of isolated horizon has been introduced in relation with loop quantum gravity [128,130,131,129,132,133,134,135] and, in a general perspective, it is too restrictive when one wants to allow mass-energy to cross the "horizon" (whichever way the latter is defined) in one direction or the other.
A weakly isolated horizon is a null surface H with null normal l a such that θ l = 0, −T ab l a is a future-oriented and causal vector, and L l n b ∇ a l b = 0. In this context l a is a Killing vector for the intrinsic geometry on H, without reference to the surroundings, and can therefore be used to define a "completely local Killing horizon" when there are no energy flows across H. The vector field l a generates a congruence of null geodesics on H, which can be used to define a surface gravity κ via the (non-affinely parametrized) geodesic equation This surface gravity κ is constant on the weakly isolated horizon H, which corresponds to the zeroth law of thermodynamics. Since the vector field n a is not unique also this surface gravity is not unique. A Hamiltonian analysis of the phase space of isolated horizons, identifying boundary terms with the energies of these boundaries, leads to a first law of thermodynamics for isolated horizons with rotational symmetry [131], where J is the angular momentum, H H the Hamiltonian, A the area of the 2-dimensional cross-sections of H, and Ω H the angular velocity of the horizon. A dynamical horizon [7] is a spacelike marginally trapped tube foliated by marginally trapped 2-surfaces (MTT). This definition allows for energy fluxes across the dynamical horizon. A set of flux laws describing the related changes in the area of the dynamical horizons have been formulated [7]. An apparent horizon which is everywhere spacelike coincides with a dynamical horizon, but an apparent horizon is not required to be spacelike. Being spacelike, dynamical horizons can be crossed only in one direction by causal curves, while this is not the case for apparent horizons which can be partially or entirely timelike.
Finally, slowly evolving horizons have also been introduced and studied [137,136,138,5]: these are "almost isolated" FOTHs and they are intended to represent black hole horizons which evolve slowly in time, as is expected in many astrophysical processes but not, for example, in the final stages of black hole evaporation. They are analogous to thermodynamical systems in quasi-equilibrium.

Kodama vector
In the literature one finds several notions of surface gravity associated with horizons. In stationary situations, in which a timelike Killing vector field outside the horizon becomes null on it, these notions of surface gravity coincide. In dynamical situations there is no timelike Killing vector and these surface gravities turn out to be inequivalent. In spherical symmetry, the Kodama vector mimics the properties of a Killing vector and originates a (miracolously) conserved current and a surface gravity. The Kodama vector [139] is a generalization of the notion of Killing vector field to spacetimes which do not admit one, and has been used in place of a Killing vector in the thermodynamics of dynamically evolving horizons. The Kodama vector is defined only for spherically symmetric spacetimes (see [140] for an attempt to introduce a Kodama-like vector in nonspherical spacetimes). Let the metric be where (a, b) = (t, R), R is the areal radius, and dΩ 2 (2) = dθ 2 + sin 2 θ dϕ 2 is the line element on the unit 2-sphere. Let ǫ ab be the volume form associated with the 2-metric h ab [3]; then the Kodama vector is [139] K In a static spacetime the Kodama vector is parallel (in general, not equal) to the timelike Killing vector.
In the region in which it is timelike, the Kodama vector defines a class of preferred observers with four-velocity u a ≡ K a / |K c K c | (the Kodama vector is timelike in asymptotically flat regions). The Kodama vector is divergence-free [139,141], which has the consequence that the Kodama energy current [139] even if there is no timelike Killing vector (a property referred to as the "Kodama miracle" [141]). By writing the spherical metric in Schwarzschild-like coordinates, the Kodama vector assumes the simple form (e.g., [139,142]) (2.33) The Noether charge associated with the Kodama conserved current is the Misner-Sharp-Hernandez energy [143,144] of spacetime (which, again, is defined only for spherically symmetric spacetimes) [145].

Surface gravities
Traditionally, surface gravity is defined in terms of geometric properties of the metric tensor and it also shows up in black hole thermodynamics as the proportionality factor between the variation of the black hole mass (which plays the role of internal energy) dM and the variation of the event horizon area (proportional to the entropy) dA. Since it is unclear which definition of black hole mass is appropriate in non-trivial backgrounds (see the review [146]), also the definition of surface gravity suffers from the same ambiguities. Surface gravity is also a semiclassical quantity since for stationary black holes it coincides, up to a constant, with the Hawking temperature of a black hole. The textbook definition of surface gravity is given on a Killing horizon [3]. Given that Killing fields are not available in non-stationary situations, a different concept of surface gravity is necessary there. The recurrent definitions are reviewed in Ref. [80] and are briefly recalled here.
Killing horizon surface gravity: a Killing horizon defines the surface gravity κ Killing as follows [3]: on the Killing horizon the Killing vector k a satisfies (e.g., [80]) In static spacetimes, κ Killing is interpreted as the limiting force required at spatial infinity to hold in place a unit test mass just above the event horizon by means of an infinitely long massless string [3] (which shows the non-local nature of the notion of Killing surface gravity). Since the Killing equation ∇ a k b + ∇ b k a = 0 determines the Killing vector k a only up to an overall normalization, there is freedom to rescale k a and the value of the surface gravity depends on the non-affine parametrization chosen for k a . However, in static/stationary situations one has the freedom of imposing that k c k c = −1 at spatial infinity. The Killing surface gravity can be generalized to any event horizon that is not a Killing horizon by replacing the Killing vector k a with the null generator of the event horizon [80].
Surface gravity for marginally trapped surfaces: let l a and n a be the outgoing and ingoing null normals to a marginally trapped (spacelike compact 2-dimensional) surface, with the expansion of l a vanishing, and assume that l a and n a are normalized so that l c n c = −1. In general, l a is not a horizon generator but is a non-affinely parametrized geodesic vector on the trapping horizon, which allows one to define a surface gravity κ as The value of κ depends on the parametrization of l a and there are various proposals for this. In general, writing l a as the tangent to a null curve x µ (λ) with parameter λ, a parameter change (dependent on the spacetime point) λ → λ ′ means that the components of l a change according to and finally The Hayward proposal (for spherical symmetry) [81] is based on the Kodama vector K a [139].
In spherical symmetry the Kodama vector satisfies and it is taken to be future-directed. The ensuing surface gravity for a trapping horizon is given by This definition is unique because the Kodama vector is unique. κ Kodama agrees with the surface gravity on the horizon of a Reissner-Nordström black hole but not with other definitions of dynamical surface gravity. An expression equivalent to (2.44) is [81] where R is the areal radius and h is the determinant of the metric h ab in the 2-space orthogonal to ∇ a R. The Hamilton-Jacobi approach, a variant of the Parikh-Wilczek method [115], employs the Kodama-Hayward definition of surface gravity [147] (for a review of tunneling methods see Ref. [148]).
The Fodor et al. surface gravity: this proposal for spherically symmetric asymptotically flat spacetimes [149] is based on the ingoing null normal n a being normalized so that n a t a = −1, where t a is the asymptotic time-translational Killing vector at spatial infinity. n a is affinely parametrized everywhere and at spatial infinity is parametrized by the proper time of static observers. Requiring that l c n c = −1 fixes the parametrization of l a , yielding The isolated horizon surface gravity: this proposal of Ashtekar, Beetle, and Fairhurst [129] applies to an isolated horizon. n a is normalized so that its expansion agrees with that of the Reissner-Nordström case and with l a n a = −1. This choice fixes a unique surface gravity as a function of the horizon parameters. However, this concept appears to be limited, for example it cannot be extended to the Einstein-Yang-Mills case [150,80].
The Booth and Fairhurst proposal for slowly evolving horizons: this definition [137] extends the previous proposal. On the isolated horizon the normal is τ a = Bl a + Cn a , with B and C scalar fields defined there, which weight the contributions of l a and n a (for an isolated horizon it is B = 1 , C = 0). The surface gravity is Other proposals: other proposals for surface gravity include Hayward's trapping gravity [22] κ trapping ≡ 1 2 −n a ∇ a θ l (2.48) and the Mukohyama and Hayward proposal [151]. The surface gravities listed here are computed in Ref. [80] for a general spherically symmetric metric in Eddington-Finkelstein coordinates and in terms of the Misner-Sharp-Hernandez mass [143,144]. Ref. [82] compares these definitions for spherical black holes in Painlevé-Gullstrand coordinates.

Spherical symmetry
Assuming spherical symmetry greatly simplifies the solution of the field equations and the study of horizons. Although this is an unrealistic assumption for rotating astrophysical black holes and for universes with realistic inhomogeneities, it is important for the fundamental theory.
In spherically symmetric spacetimes a useful tool is the Misner-Sharp-Hernandez mass M [143,144], which here coincides with the Hawking-Hayward quasi-local mass [152,153]. The Misner-Sharp-Hernandez mass is defined in GR and for spherical symmetry. Using the areal radius R and angular coordinates (θ, ϕ), a spherical line element can be written as The Misner-Sharp-Hernandez mass M is defined by [143,144] 1 Horizons in spherical symmetry are discussed in a clear way in the formalism of Nielsen and Visser [109,119]. These authors consider the most general spherically symmetric metric (not necessarily stationary or asymptotically flat) with a spherically symmetric spacetime slicing, which assumes the form in Schwarzschild-like coordinates, where M (t, R) a posteriori turns out to be the Misner-Sharp-Hernandez mass. The line element (2.52) can be recast in Painlevé-Gullstrand coordinates as where φ(τ, R) and M (τ, R) are implicit functions of (τ, R) and the spacelike hypersurfaces τ =constant are flat. Using the implicit functions of (τ, R) [109] c the line element becomes A number of practical results are then obtained [109]. The outgoing radial null geodesic congruence has tangent field with components (in Painlevé-Gullstrand coordinates (τ, R, θ, ϕ)) while the ingoing radial null geodesics have tangent field where the normalization g ab l a n b = −2 (2.59) is adopted [109]. The expansions of these radial null geodesic congruences are The apparent horizon corresponding to θ l = 0 and θ n < 0 is given by where the last equation holds in both Painlevé-Gullstrand coordinates and in the gauge (2.52) and is obtained by using the fact that the inverse of the metric (2.56) has components (2.63) In practice, the condition g RR = 0 is a very convenient recipe to locate the apparent horizons in the presence of spherical symmetry when the areal radius R is used as a coordinate and it is often convenient to perform a coordinate transformation to this radial coordinate and to rewrite the line element using R.
The gradient of the areal radius R and the normal n a = ∇ a R to the surfaces R =const. become null at the apparent horizon; this recipe is reminiscent of the change in the causal character of the Schwarzschild radial coordinate on the Schwarzschild event horizon. However the apparent horizon is not, in general, a null surface. We have also [109] where a prime and an overdot denote partial differentiation with respect to R and τ , respectively, and the subscript AH identifies quantities evaluated on the apparent horizon. 1 − 2M ′ AH > 0 is required for the horizon to be outer in a spacetime with regular asymptotic region, hence the condition for the apparent horizon to be also a trapping horizon is [109,119] If matter satisfies the null energy condition, and assuming the Einstein equations, the area of the apparent horizon cannot decrease. Various energy fluxes across the apparent horizon are also discussed and computed in [109]. The Nielsen-Visser surface gravity at the horizon is computed from l b ∇ b l a = κ l l a , which gives [109] κ An extremal horizon will be one with vanishing surface gravity, The fact that the Misner-Sharp-Hernandez mass can be used to define and locate apparent horizons in spherically symmetric spacetimes shows that the apparent horizon is a quasi-local concept and is independent of the global causal structure. However, it does not appear to be a completely local notion (it depends on a surface, not only on the spacetime point).

Evolving horizons, cosmological black holes, and naked singularities in GR
Let us review briefly some dynamical and spherically symmetric solutions of the Einstein equations of particular significance, paying attention to the structure and dynamics of their apparent horizons.

The Schwarzschild-de Sitter-Kottler spacetime
The Schwarzschild-de Sitter-Kottler [154] spacetime is locally static but it is useful to review it in order to understand the apparent horizons of more complicated dynamical solutions. It has line element where the constant H = Λ/3 is the Hubble parameter of the de Sitter background, Λ > 0 is the cosmological constant, and m > 0 is a second parameter related to the mass of the central inhomogeneity (e.g., [79,155]). The static coordinates (t, R, θ, ϕ) cover the region R 1 < R < R 2 . The apparent horizons are located by g RR = 0, which is equivalent to the cubic equation with roots with sin(3ψ) = 3 √ 3 mH. m and H are both necessarily positive in an expanding universe, then R 3 is negative and there are at most two apparent horizons. When R 1 and R 2 are real, R 1 is a black hole apparent horizon which reduces to the R = 2m Schwarzschild horizon in the limit H → 0, while R 2 is a cosmological apparent horizon which reduces to the R = 1/H de Sitter horizon in the limit m → 0. The metric (3.68) is static in the region between these two horizons.
Both apparent horizons exist only if 0 < sin(3ψ) < 1 and, since the metric is locally static, the apparent black hole and cosmological horizons are also event horizons. If sin(3ψ) = 1 these horizons coincide (extremal Nariai black hole). For sin(3ψ) > 1 the roots are complex-valued and there is a naked singularity. To summarize there are no apparent horizons. The interpretation seems to be that the would-be black hole horizon would become larger than the cosmological one but, strictly speaking, the roots corresponding to the apparent horizons are complex in this case.
The black hole horizon has area A = 4πR 2 1 which is, of course, time-independent. In the non-extremal case the central singularity is eternal and spacelike ( [156], see this reference also for a conformal diagram) and is surrounded by the black hole event horizon at all times for the parameter values for which this horizon exists.
A sphere of radius R has the Misner-Sharp-Hernandez mass where ρ = Λ 8π . The Schwarzschild-de Sitter-Kottler black hole has been studied extensively in relation to its thermodynamics. Here we do not discuss anti-de Sitter black holes corresponding to Λ < 0, which are the subject of much recent interest due to the fluid-gravity duality [103,104].

The McVittie solution
The 1933 McVittie solution of the Einstein equations [34] is a generalization of the Schwarzschildde Sitter-Kottler solution and represents a central object embedded in a FLRW (not necessarily a locally static de Sitter) background. Even after many works [157,44,158,159,160,36,39,40,41,37,38,161], this solution is not completely understood. The McVittie solution with negative cosmological constant was analyzed in Ref. [163] and an electrically charged version of the McVittie spacetime was found in Ref. [162]. In this subsection we restrict to a spatially flat FLRW background and to zero electric charge.
A simplifying assumption of McVittie consists of the no-accretion condition Gr t = 0 (in spherical coordinates, where G µν is the Einstein tensor) which forbids any mass-energy flow (which, in spherical symmetry, could only be radial), T tr = 0. Generalizations of the McVittie solution allowing radial energy fluxes are more complicated and will be considered later.
McVittie was led to his solution [34] by the problem of the effect of the cosmological expansion on local systems. Different approaches to this problem generated other solutions, such as the Swiss-cheese model [164,165] (this problem has seen an extensive literature devoted to it but is not completely solved [35]). Unlike the Schwarzschild-de Sitter-Kottler spacetime, black holes in more general FLRW backgrounds are dynamical.
The McVittie line element in isotropic coordinates is where the function m(t) is required to satisfy the McVittie no-accretion condition Tr t = 0 on the stress-energy tensor T ab , which becomeṡ where m 0 is a constant, therefore, The McVittie metric reduces to the Schwarzschild one in isotropic coordinates if a ≡ 1 and to the FLRW metric if m 0 = 0 and is singular on the 2-spherer = m 0 /2 (which reduces to the Schwarzschild horizon if a ≡ 1) [166,158,159,160,157]. This singularity is spacelike [158,159,160,37,161] (and is represented as a horizontal line in conformal diagrams). There is another spacetime singularity atr = 0. McVittie's original interpretation of the line element (3.72) as describing a point mass atr = 0 is made untenable by the fact that this point mass would be surrounded by ther = m 0 /2 singularity [157,166,158,159,160]. We will only consider the regionr > 2m 0 here. The energy density of the source fluid is finite but its pressure diverges atr = m 0 /2 with the Ricci scalar R a a = 8π (3P − ρ) [157,166,158,159,160,96,97], with the exception of a de Sitter background withḢ = 0 [158,159,160,167,38].
The apparent horizons were studied in Refs. [159,42,41] and interpreted in [41], which we follow here. We rewrite the line element (3.72) in terms of the areal radius the differentials dr and dR are related by where the relation (3.73), which gives where m/r = ma/R = m 0 /R (here ma is constant because of eq. (3.73)) leads to where H ≡ȧ/a. The cross-term in dtdR is then eliminated by defining a new time T (t, R) such that with the integrating factor F (t, R) and function β(t, R) to be determined. dT is an exact differential if the 1-form (3.83) is closed, or Now replace dt with F dT − βdR in eq. (3.82), obtaining (3.85) Imposing now that the line element is diagonalized, The singularityr = m/2 corresponds to the proper radius R = 2m a(t) = 2m 0 and does not expand. Let us study now the apparent horizons of the McVittie spacetime [158,159,41]. For simplicity, we restrict ourselves to a spatially flat FLRW background. The Einstein equations provide the density of the fluid The McVittie metric admits arbitrary FLRW backgrounds generated by cosmic fluids satisfying any constant equation of state. For brevity, we restrict to a cosmic fluid which reduces to dust at spatial infinity and corresponding to an equation of state parameter w = 0. Then the pressure is [39,40] The apparent horizons are located by g RR = 0 or, using eq. (3.87) 1. for t < t * it is m 0 > 1 3 √ 3 H(t) and both R 1 (t) and R 2 (t) are complex. There are no apparent horizons.
2. The critical time t = t * corresponds to m 0 = 1 3 √ 3 H(t) . R 1 (t) and R 2 (t) coincide at a real value and there is a single apparent horizon at 3 H(t) and there are two apparent horizons of real positive radii R 1 (t) and R 2 (t).
The behaviour of the apparent horizons is described in fig. 2. At times t < t * there is a naked singularity at R = 2m 0 : while the Hubble parameter H(t) diverges near the Big Bang, the mass coefficient m 0 stays supercritical at m 0 > 1 3 √ 3 H(t) . As the Schwarzschild-de Sitter-Kottler experience teaches us, a black hole horizon cannot be accommodated in this small universe and the singularity is naked: the putative black hole is too large to fit in the observable universe (varying speed of light cosmologies have a related phenomenology -the radii of primordial black holes and the Compton wavelengths of massive particle states can become larger than the Hubble radius [43]). At the critical time t * a black hole apparent horizon and a cosmological apparent horizon appear together at radius R 1 (t * ) = R 2 (t * ) = 1 √ 3 H(t * ) , in analogy with the Nariai black hole of the Schwarzschild-de Sitter-Kottler solution. This critical black hole is instantaneous. As time progresses to t > t * , this single horizon splits into an evolving black hole apparent horizon surrounded by an evolving cosmological horizon. The black hole apparent horizon shrinks, asymptoting to the spacetime singularity at 2m 0 from above as t → +∞, while the cosmological apparent horizon expands monotonically, tending to 1/H(t) in the same limit.
The well known singularity R = 2m 0 [158,159,160,39,40,37], where the Ricci scalar diverges, separates the two disconnected spacetime regions R < 2m 0 and R > 2m 0 [158,159,160] and is spacelike [41]. One can compare the rate of change of the apparent horizon radii with respect to that of the cosmic substratum, obtaining [41] R AH the apparent horizons are not comoving except for trivial cases. The sum of the areas of the two apparent horizons of the McVittie spacetime is a non-decreasing function of time but undergoes a discontinuous jump from zero at the critical time t * [41].

A phantom background
A background FLRW universe dominated by a phantom fluid with equation of state parameter w ≡ P/ρ < −1 which violates the weak energy condition can be considered. Phantom fluids, studied in conjunction with the present cosmic acceleration [28], cause a Big Rip singularity at a finite future t rip [168]. A phantom background universe for the McVittie solution was studied in Ref. [41]. The scale factor of a phantom-dominated spatially flat FLRW universe is where A is a constant. The Hubble parameter  . 3). An idealized interior solution for the McVittie metric describing a relativistic star of uniform density in a FLRW background was found by Nolan [169] and it generalizes the Schwarzschild interior solution with a Minkowski background [3], to which it reduces when a =const. It belongs to the Kustaanheimo family of shear-free solutions. The star surface is comoving with the cosmic substratum [167]. The generalization of the Tolman-Oppenheimer-Volkoff equation [3] for this crude star model was written down in [167].
Recent works on the McVittie spacetime study its conformal structure [36,37,161], which means integrating numerically the null geodesics or deriving general analytical results upon assuming something on the expansion. Lake and Abdelqader [37] find that null geodesics asymptote to the singularity without entering it. Depending on the form of the scale factor, a bifurcation surface may appear which splits the spacetime boundary into a black hole horizon in the future and a white hole horizon in the past. This behaviour seems a reflection of the McVittie no-accretion condition which applies to a timelike dust and, in the limit, also to a null dust. da Silva et al. [161] find that the presence or absence of this white hole horizon depends crucially on the expansion history of the universe a(t), and prove a theorem in this regard for McVittie spacetimes for which the background is non-superaccelerating (i.e.,Ḣ ≤ 0) and de Sitter at late times. It would be desirable to extend the result to backgrounds which at late times asymptote to any FLRW space, not just de Sitter. See Refs. [37,161] for the corresponding conformal diagrams.

Area quantization and McVittie solutions as toy models
As an example of the use of cosmological black holes as toy models to exemplify unintuitive physics, we quote the current issue of the quantization of black hole areas. Inspired by certain stringy black holes, there has been excitement in the string community about the fact that the areas A ± of black hole inner (−) and outer (+) horizons satisfy the relation or where l pl is the Planck length [170,171]. These area-quantizing relations have somehow come to be seen as universal [172,173,174,175]. While certain stringy black holes remarkably do satisfy these relations, this property is certainly not universal, as shown by Visser [176,177]

Generalized McVittie solutions
Generalized McVittie solutions with spacetime metric of the form (3.72) but without the noaccretion restriction (3.73) were introduced in Ref. [167]. In principle such metrics could be meaningless: in the "Synge approach" one can always impose that an invented metric solves the Einstein equations and run these equations from left to right to compute the corresponding formal stress-energy tensor T ab . This T ab is usually found to be completely unphysical and violates all reasonable energy conditions, beginning with the positivity of the energy density. Rather surprisingly, generalized McVittie solutions with reasonable matter sources exist. In isotropic coordinates, generalized McVittie solutions can be presented as where m(t) ≥ 0 and The only non-vanishing components of the mixed Einstein tensor are

100)
(the unusual feature that Gr r = G θ θ is named "spatial Ricci isotropy" in Ref. [179]). The quantity diverges on the surfacer = m/2 unless m is a constant. Imperfect fluids can be contemplated as matter sources for this metric.

Single perfect fluid
If the matter source of the generalized McVittie metric is a single perfect fluid with stress energy tensor T ab = (P + ρ) u a u b + P g ab (3.106) and a radial fluid flow described by the fluid four-velocity u µ = u 0 , u, 0, 0 , is allowed, then the only possible solution of the Einstein equations is the Schwarzschild-de Sitter-Kottler black hole [167,65,37]. This is easily seen, since the normalization u c u c = −1 yields and, using eqs. (3.99)-(3.101), the Einstein equations imply thaṫ is the area of a sphere of isotropic radiusr and

Imperfect fluid and no radial mass flow
Consider now the imperfect fluid stress-energy tensor as a source for the generalized McVittie solutions, where the purely spatial vector q c describes a radial energy flow, and u c u c = −1 (in principle one could take q c to be spacelike instead of purely spatial [179,38]). 3.116) and the area of a sphere Σ of constant time and constant isotropic radiusr is

The (t,r) component of the Einstein equations yieldṡ
then energy flow, area A, and accretion rate are related bẏ In the case of inflow (q < 0), this condition can be written on a sphere of radiusr ≫ m aṡ M ≃ aA |q|; for a 2-sphere, M increases due to the inflow of matter alone (but it receives another contribution from the evolution of the cosmological fluid contained in it). The energy density and pressure obtained from the Einstein equations are clearly the energy density is always non-negative. In terms of the quantity C, eq. (3.120) becomes the generalization of the Raychaudhuri equation of FLRW spacė It reduces to the usual Raychaudhuri equation of FLRW cosmology in the limit m → 0,

Imperfect fluid and radial mass flow
Let us consider now an imperfect fluid with stress-energy tensor of the form (3.113) with both radial mass flow and energy current present and of the form By using the components (3.99)-(3.101) of the Einstein tensor, the field equations becomė Adding the last two equations yields (equivalently, this equation can be seen as a consequence of the spatial Ricci isotropy Gr r = G θ θ ), i.e., to an ingoing radial mass flow there corresponds an outgoing radial heat current if P > −ρ. By substituting eq. (3.130) into eq. (3.126), one obtains the accretion ratė where (P + ρ) Au can be seen as the flux of gravitating energy through the surface of area A (remember that u < 0). The energy density is given by (3.132)

The "comoving mass" solution
In the class of generalized McVittie solutions of GR, the choice M (t) = m 0 a(t) where m 0 is a constant, selects a special one which is a late-time attractor within this class. The corresponding line element in isotropic coordinates is The apparent horizons of this metric were studied in Ref. [65] by transforming to areal radius. A nice feature of this solution is that the apparent horizons are given analytically by R c is a cosmological and R b is a black hole apparent horizon. The surface r = m 0 /2 (or r = 2m 0 , or R = 2m 0 a = 2M (t)) is a spacetime singularity contained inside the black hole apparent horizon when the latter exists, since R c,b > 2m 0 a = 2M [65]. The black hole and cosmological apparent horizons have the qualitative behaviour already discussed for all the McVittie and generalized McVittie solutions [65].

The general class of solutions
For the wider class of generalized McVittie solutions with arbitrary dependence m(t) ≥ 0, an analysis of the apparent horizons using the areal radius [65] identifies them as the roots of the equation wherer ≡ R/a. Since M (t) = m(t)a(t), the left hand side can be written as where the factor M 1 + m 2r quantifies the deviation of the radius from 2M (r > m/2 corresponds to R > 2M and to M 1 + 2m r > 2M ), while the factor Ṁ M − H is the difference between the percent rate of change of M and that of the scale factor of the substratum. The vanishing of this factor corresponds to an analog of stationary accretion for a time-dependent background. Then, the special solution with M (t) = m 0 a(t) corresponds to stationary accretion relative to the FLRW background.
Eq. (3.136), which becomes is not a quadratic algebraic equation, but it can be treated formally as such, providing the formal roots Since r = r(R), this is really an implicit equation for the radii R c,b of the cosmological and black hole apparent horizons. When the argument of the square root is positive there are a cosmological apparent horizon at R c and a black hole apparent horizon at R b . When this argument vanishes, these two apparent horizons coincide at 2M H . If this argument becomes negative, the apparent horizons disappear leaving behind a naked singularity [65].
3.4.6 Attractor behaviour of the "comoving mass" solution "Comoving mass" solutions are generic under certain assumptions, in the sense that all other generalized McVittie solutions approach them at late times [180]. In fact, assume that the universe always expands, that m(t) ≥ 0, and that the function m(t) is continuous with its first derivative. Then, usingr ≡ R/a, one obtains Since m ≥ 0 the left-hand side is always non-negative andṁ 1 + m 2r < 1 a . Then, given that 1 + m 2r > 0, in an expanding universe in which a → +∞, one hasṁ ∞ ≡ lim t→+∞ṁ (t) ≤ 0. Iḟ m ∞ = 0, the quantity m(t) becomes asymptotically comoving.
In the first case one has, at t = t * , HR = |ṁ * | a + 1, which yields the radius of the black hole apparent horizon at t * Late in the history of the universe we have a black hole of zero mass M (t * ) = a(t * )m(t * ) but finite radius r * . As time evolution continues, one would have negative mass M and finite radius of the black hole apparent horizon. This unphysical situation for m(t * ) = 0 with m(t > t * ) < 0 is discarded. In the second caseṁ(t) → 0 at late times and t → +∞ if the cosmic expansion continues forever or t → t rip if a Big Rip occurs at t rip . The physical meaning ofṁ → 0 is that, at late times, the rate of increase of the black hole mass is at most the Hubble rate and becomes comoving.

The Sultana-Dyer solution
The Sultana-Dyer solution of GR [95] is a Petrov type D metric interpreted as a black hole embedded in a spatially flat FLRW universe. This solution was generated by extending a metric resulting from the conformal transformation of the Schwarzschild metric g ab with conformal factor Ω = a(t) = η 2 equal to the scale factor of a dust-filled k = 0 FLRW universe in conformal time η. That is, this spacetime is conformally static and admits a conformal Killing vector ξ a (fig. 4).
The authors of [95] aimed at changing the Schwarzschild timelike Killing field ξ c into a conformal Killing field defined for ξ c ∇ c Ω = 0, thus generating a conformal Killing horizon (which, however, seems of little relevance in modern studies of time-evolving horizons).
where m 0 is a constant and a(η) = η 2 . The coordinate transformation turns the line element into the form which is explicitly conformal to the Schwarzschild metric with conformal factor The matter source of the Sultana-Dyer spacetime is a mixture of two non-interacting perfect fluids with stress-energy tensor T ab = T = ρ n k a k b describes a null dust with density ρ n and k c k c = 0 [95]. A problem of this solution is that the cosmological fluid becomes tachyonic with negative energy density at late times nearr = m 0 /2 [95].
Let us use, in the rest of this subsection, the quantity which is not constant in the Sultana-Dyer solution. The locus r = 2m 0 is not a singularity, but the conformal factor Ω vanishes there. The metric (3.142), however, is not singular there. The Ricci curvature is and is not singular at r = 2m 0 (where η → −∞) but is singular at r = 0 (central singularity) and for η = 0 (Big Bang singularity). The problem of Hawking emission from the Sultana-Dyer black hole was approached in Ref. [181]. These authors considered quantum radiation from a massless conformally coupled scalar field φ and computed the renormalized stress-energy tensor T ab of φ taking advantage of the simplifications introduced by the fact that the Sultana-Dyer spacetime is conformal to the Schwarzschild one and taking into account the conformal anomaly and particle creation by the FLRW background. Discarding complicated corrections which are small if the black hole is evolving slowly, the effective Hawking temperature from the Sultana-Dyer black hole was computed as [181] T ef f = 1 8πm 0 a(t) = T Schw a(t) , (3.149) where T Schw = (8πm 0 ) −1 is the Hawking temperature of the Schwarzschild black hole which seeds the Sultana-Dyer metric. The more general relation for spacetimes conformally related to the Schwarzschild spacetime by a trasformation with conformal factor Ω is conjectured in [181]. Independent support for eq. (3.150) comes from dimensional considerations related to the use of conformal transformations [182].

The Husain-Martinez-Nuñez solution
In the 1994 Husain-Martinez-Nuñez solution of GR [183] a new phenomenology of the apparent horizons appears. This spacetime describes an inhomogeneous universe with a spatially flat FLRW background sourced by a free, minimally coupled, scalar field. The coupled Einstein-Klein-Gordon equations reduce to 152) and the Husain-Martinez-Nuñez solution to them is [183] where A 0 , B 0 , C, and D are non-negative constants, α = ± √ 3/2, and η > 0. The additive constant B 0 becomes irrelevant and can be dropped whenever A 0 = 0. When A 0 = 0, the Husain-Martinez-Nuñez metric degenerates into the static Fisher spacetime [184] where V (r) = 1 − 2µ/r, µ and ν are parameters, and the Fisher scalar field is The Fisher solution of the coupled Einstein-Klein-Gordon equations, also referred to as the Janis-Newman-Winicour-Wyman solution, has been rediscovered many times [185,186,187,188,189,190]. Its features are a naked singularity at r = 2C and its asymptotic flatness. It is claimed that this solution is the most general static and spherically symmetric solution of the Einstein equations with zero cosmological constant and a massless, minimally coupled, scalar field [191], but it is unstable [192]. The general Husain-Martinez-Nuñez metric is conformal to the Fisher metric with conformal factor Ω = √ A 0 η + B 0 equal to the scale factor of the background FLRW space and with only two possible values of the parameter ν. From now on, we set the constant B 0 to zero by labelling the Big Bang by η = 0. The sign in eq. (3.154) is not associated with the sign of α. The full metric is asymptotically FLRW for r → +∞ and is FLRW if C = 0 (in which case the constant A 0 can be eliminated by rescaling the time coordinate η). The Ricci scalar immediately identifies a spacetime singularity at r = 2C (for both values of the parameter α).
The scalar φ also diverges there, and a Big Bang singularity is present at η = 0. Only the coordinate range 2C < r < +∞ is physical and the lower limit r = 2C corresponds to zero areal radius Let us introduce the comoving time t defined by dt = adη (where a(η) = √ A 0 η is the FLRW scale factor) in place of the conformal time η, then it is by choosing η = 0 at t = 0, or and This power law for the scale factor is consistent with the stiff equation of state P = ρ/3 of a free massless scalar field in a FLRW universe and with the general solution a(t) = const. t The areal radius (3.158) increases with r for r > 2C. It is useful to rewrite the line element in terms of the areal radius R. By setting we have R(r) = a(t)rA 1−α 2 (r) and .
The metric is then The time-radius cross-term is eliminated by introducing a new time T with differential where β(t, R) is a function to be determined and F (t, R) is an integrating factor which must satisfy ∂ ∂R in order for dT to be an exact differential. Using dt = F dT − βdR in eq. (3.166) and choosing the line element becomes (3.170) The apparent horizons, located by g RR = 0, must satisfy where now r = r(t, R), or using the original coordinates (η, r) [183]. For r → +∞ (corresponding to R → +∞), this equation reduces to R ≃ H −1 , the radius of the cosmological apparent horizon in spatially flat FLRW space. Eq. (3.171) must be solved numerically. Let x ≡ C/r, then the equation locating the apparent horizons is The left hand side can be written as This is a parametric representation of the function R(t) and can be used to plot this function. The result is illustrated in figs. 5 and 6. If α = √ 3/2, between the Big Bang and a critical time t * there is only one expanding apparent horizon, then two other apparent horizons are created at t * . One is a cosmological apparent horizon which expands forever and the other is a black hole horizon which contracts until it meets the first (expanding) black hole apparent horizon [183]. When they meet, these two annihilate and a naked singularity appears at R = 0 in a FLRW universe. This phenomenology of apparent horizons differs from that of the McVittie Figure 6: The radius of the Husain-Martinez-Nuñez apparent horizon (vertical axis) versus comoving time (horizontal axis) for α = − √ 3/2. There is always only one, expanding, cosmological apparent horizon and there is a naked singularity at R = 0. and generalized McVittie solutions. The "S-curve" phenomenology of fig. 5 appears also in Lemaître-Tolman-Bondi spacetimes already for a dust fluid much simpler than a scalar field [83] (multiple "S"s are possible, for example see fig. 9 of Ref. [83]) and in analytical solutions of Brans-Dicke and f (R c c ) gravity. The scalar field is regular on the apparent horizons.
For α = − √ 3/2 there is only one, forever expanding, cosmological apparent horizon and the universe contains a naked singularity at R = 0 ( fig. 6), with the usual Big Bang singularity at t = 0.
The apparent horizons are spacelike [183], as can be seen by studying the normal vector to these surfaces and checking that it always lies inside the light cone in an (η, r) diagram. Eq. (3.172) yields along the apparent horizons. Differentiate this relation with respect to R to obtain Along radial null geodesics it is which follows from ds 2 = 0 with dθ = dϕ = 0. Therefore, it is [183] η ,r AH η ,r light cone 3.180) and the normal to the apparent horizons is always pointing inside the light cone, except at the spacetime points at which this vector becomes tangent to the light cone and is null, which occurs when a pair of apparent horizons is created or destroyed [183]. This occurrence is in agreement with a general result of Ref. [83] stating that a trapping horizon created by a massless scalar field must be spacelike (however, even simple potentials V (φ) can make the trapping horizon be non-spacelike). The nature of the singularity at r = 2C (or R = 0) is easily established: all surfaces R =const. have equation Φ(R) ≡ R − const. = 0 and gradient N µ ≡ ∇ µ Φ = δ µR in coordinates (t, R, θ, ϕ). The norm squared is and, because B(r) → 1−α 2 and A(r) → 0 + as r → 2C + , it is N c N c > 0 and N c N c → +∞ as r → 2C + . The singularity at R = 0 is timelike for both values of the parameter α.
For α = + √ 3/2, the Husain-Martinez-Nuñez spacetime is interpreted as describing the creation and annihilation of pairs of black hole apparent horizons. The central singularity at R = 0 is created with the universe in the Big Bang and does not result from a collapse process (this is also the case for α = − √ 3/2).

The Fonarev and generalized Fonarev solutions
The Fonarev solution of the Einstein equations with a minimally coupled scalar field in an exponential potential as the matter source [193] generalizes the Husain-Martinez-Nuñez solution.
It describes a central inhomogeneity embedded in a scalar field FLRW universe. The action is where κ ≡ 8πG and and V 0 and λ are two positive constants (this potential has been investigated in great detail in cosmology). The coupled Einstein-Klein-Gordon equations simplify to The spherically symmetric Fonarev line element and scalar field are with w and A 0 constants and η is the conformal time. For simplicity we choose A 0 = 1. When w = 0 the metric (3.186) reduces to a spatially flat FLRW one while, when a ≡ 1 and α = 1, it degenerates into the Schwarzschild solution (however, the value α = 1 is not possible if α = λ √ λ 2 +2 ). The line element becomes asymptotically that of spatially flat FLRW space as r → +∞. The Husain-Martinez-Nuñez class of solutions (3.153) is recovered by setting λ = ± √ 6 and V 0 = 0. See Refs. [193,257] for the corresponding conformal diagrams.

A generalized Fonarev solution
A generalized Fonarev solution corresponding to a dynamical phantom scalar field solution of GR is known [65]. It is obtained from the Fonarev solution via the transformation The corresponding action is (3.191) and it contains a phantom field endowed with the "wrong" sign of the kinetic term. The generalized Fonarev line element representing a dynamical black hole immersed in a phantom FLRW background is 193) where Assuming that λ > √ 2, it is of interest to understand the physical meaning of the constant w. When λ ≫ √ 2 it is a ≈ 1 and α ≈ −1 and the metric approximates to The coordinate transformation [65] y = r 1 − 2w r , (3.197) transforms the line element (3.196) into (3.198) this is the Schwarzschild spacetime with mass −w. The parameter w corresponds to the negative of the mass in this limit and from now on we will use −M instead of w. Let us locate the apparent horizons as the parameters M and α vary. This phantom black hole solution can be cast in the form the replacement of the conformal time η with the comoving time t leads to where the integration constant t 0 marks the time of the Big Rip and it is α < −1 since λ > √ 2. The exponent α is determined by the slope of the potential according to eq. (3.194). When M = 0 the spacetime (3.200) reduces to a phantom-dominated FLRW cosmos. By setting, for simplicity, α = −3 or λ = 3/2, the line element (3.200) reduces to In terms of the areal radius R = ar (1 + 2M/r) −1 , the equation locating the apparent horizons is where H ≡ȧ/a is the Hubble parameter of the background. Further setting This quartic equation has only two real positive roots corresponding to a cosmological apparent horizon R c and a black hole apparent horizon R b [65]. The qualitative behaviour of the apparent horizons is the same as that of the McVittie and generalized McVittie classes of solutions with a phantom FLRW substratum: a black hole apparent horizon inflates while a cosmological apparent horizon shrinks. At a critical time these two apparent horizons meet and disappear leaving behind a naked singularity [65].

The Swiss-Cheese model
In 1945, apparently unaware of McVittie's work from a decade earlier, Einstein and Straus [164,165] constructed the solution of GR now called "Einstein-Straus vacuole" or "Swisscheese model" by pasting a Schwarzschild-like region of spacetime onto a dust-dominated FLRW universe across a timelike hypersurface (for reviews see [44,35,199], which we partially follow here). There is a black hole event horizon in this spacetime and the usual energy conditions are satisfied. Let the interior Schwarzschild region be denoted with M − and the exterior FLRW region with M + and let Σ be a spacelike 2-sphere of constant comoving radius r Σ . The coordinate charts covering Σ are the FLRW {t, θ, ϕ} and the Schwarzschild chart {T (t), θ, ϕ}. The metric in the two regions is given by . (3.205) Choose on Σ the triad of orthonormal vectors e α (t) , e α (θ) , e α (ϕ) = δ α t , δ α θ ar , δ α ϕ ar sin θ , (3.206) where α, β = r, θ, ϕ. The equation of Σ is Φ(r) ≡ r − r Σ = 0 and the gradient of Φ is N a ≡ ∇ a Φ = δ ar , with norm squared N a N a = g rr = 1−kr 2 a 2 . The unit normal to Σ, therefore, has components The extrinsic curvature of Σ is given by the usual formula β ∇ a n b (3.208) which is used to compute K αβ in M − and M + . The continuity of the first and second fundamental forms on Σ requires the combination of which yields Using the Hamiltonian constraint of FLRW space H 2 = 8π 3 ρ − k a 2 , one obtains In the absence of surface distributions of mass-energy on Σ, the stress-energy tensor of the matter source of this solution of the Einstein equations must also be continuous across Σ.
Since the interior is vacuum, the pressure on the outside is forced to vanish, P (+) = P (−) = 0, which implies that only a dust-dominated FLRW background can match the Schwarzschild solution. Moreover, the energy density must be continuous at Σ, implying that 214) which means that the mass of the black hole inside the vacuole must equal the mass that a sphere of volume 4πR 3 Σ /3 would have in the FLRW background (note that this volume is not the proper volume of such a sphere unless the FLRW curvature index k vanishes). This condition yields 8π 3 ρR 2 Σ = 2m ar Σ . Eq. (3.213) then reduces to The continuity of the matter distribution across Σ can be seen as the continuity of the Misner-Sharp-Hernandez mass M (+) M SH [194] (see [35,44] for a detailed discussion). The Einstein-Straus model has no accretion onto the central inhomogeneity. The interior Schwarzschild region is shielded from the cosmological expansion (and also the exterior FLRW region sees no effect from the central hole) and is static and, because of this fact, the Swiss-cheese model is often used as supporting evidence that the cosmological expansion does not affect local systems. However, the boundary of the vacuole is expanding and perfectly comoving; if the vacuole is regarded as the "local object" (instead of the black hole in it, which is insulated by a vacuum region), this argument fails. The Einstein-Straus vacuole has few drawbacks: it is unable to describe the Solar System [195,44] and is unstable with respect to non-spherical perturbations [196,197,158,198] and to perturbations of the matching condition M  [44]. The Einstein-Straus vacuole was generalized to include a cosmological constant, obtaining a Schwarschild-(anti-)de Sitter instead of Schwarzschild interior [200], or to include a fluid with pressure in the interior region [201]. Also the generalization obtained by matching a Schwarschild interior with an inhomogeneous Lemaître-Tolman-Bondi exterior has been studied [202]. The Hawking radiation emitted by the Einstein-Straus black hole has been studied in [203,181]. It is found that such a black hole in an expanding universe is excited to a nonequilibrium state and emits with stronger intensity than a thermal one.

Other GR solutions
There are several other analytical solutions of the Einstein equations describing central inhomogeneities in FLRW backgrounds. While one has to be careful as many of them do not have reasonable matter sources, they are of some interest. They cannot be included here due to space limitations (for a more general and rigorous treatment of inhomogeneous cosmologies see the book by Krasiński [44]). They include, among others, the well known Lemaître-Tolman-Bondi and Szekeres solutions (e.g., [44,206,83,207,208,209,210,211,212,213]), the Oppenheimer-Snyder solution [214,9,83], members of the large Barnes family [215], the solutions of Dyer, McClure, and collaborators [96,97,98,99,100], the Roberts solution with a scalar field [204,205], Patel and Trivedi's [216] and Vaidya's Kerr-FLRW solutions [217], Balbinot's evaporating black hole [218], and other solutions can be obtained from the previous ones with cut-and-paste techniques [219,220], possibly to excise regions in which the energy conditions are violated. Also asymptotically flat metrics describing transient and timedependent horizons have been studied [221,222,223,224,225,226,227,228,229] and other solutions are of potential interest [230,231,232,233,234,235,236,237,238,239,240,241].

Some cosmological black holes and naked singularities in alternative gravity
Few solutions of theories of gravity alternative to GR and representing cosmological black holes at least part of the time are known, most of them in scalar-tensor gravity. Here we review a few. The simplest scalar-tensor theory, Brans-Dicke gravity, is described by the action [30] where L (m) is the matter Lagrangian, φ is the Brans-Dicke scalar field (roughly speaking, the inverse of the effective gravitational coupling strength), and ω is a parameter ("Brans-Dicke coupling"). In more general scalar-tensor theories [31,32,33], the Brans-Dicke coupling is promoted to a function of φ, ω = ω(φ).

The conformal cousin of the Husain-Martinez-Nuñez solution
A solution of Brans-Dicke gravity was generated, but not interpreted, by Clifton, Mota, and Barrow [242] by conformally transforming the Husain-Martinez-Nuñez solution, g with φ −→φ = 2ω+3 16π ln φ. This is the inverse of the usual transformation from the Jordan frame to the Einstein frame which turns gravity with a scalar field non-minimally coupled to the Ricci scalar into GR with a scalar field with canonical kinetic energy but non-minimally coupled to matter. The Clifton-Mota-Barrow solution is where There are singularities at r = 2C and at t = 0 (here it must be 2C < r < +∞ and t > 0). The scale factor of the spatially flat FLRW background is The solution was interpreted in Ref. [243]. We rewrite the 2-parameter line element as where and is the areal radius. It is useful to study the area of the 2-spheres of symmetry: we have The critical value r 0 lies in the physical spacetime region r 0 > 2C if Θ < −1. R has the limit R(t, r) = r a(t) For Θ < −1, the areal radius R(r) has a minimum at r 0 , the area 4πR 2 of the 2-spheres of symmetry is minimum there, and there is a wormhole throat joining two spacetime regions. Since for α = ± √ 3/2, the condition Θ < −1 requires α = + √ 3/2 (this is a necessary but not sufficient condition for the throat to exist). The sufficient condition Θ < −1 constrains the Brans-Dicke parameter as [243] ω < 1 2 For −3/2 < ω < ω 0 the solution can be interpreted as a cosmological Brans-Dicke wormhole. The region 2C < r < r 0 is not a FLRW region and the scalar field is finite and non-zero at r 0 : The wormhole throat is exactly comoving with the cosmic substratum, which is relevant for the problem of cosmic expansion versus local systems [35] and disappears in the limit C → 0. Let us study the existence and location of the apparent horizons of this spacetime. The relation between differentials , (4.230) turns the line element into where (4.232) Starightforward manipulations yield where H ≡ȧ/a. The inverse metric in coordinates (t, R, θ, ϕ) is The apparent horizons are located by the roots of g RR = 0, or There are solutions which describe apparent horizons with the "S-curve" phenomenology of the Husain-Martinez-Nuñez solution of GR. Eq. (4.235) is satisfied also if the right hand side is time-independent, H = γ/t = 0, γ = 0, β = √ 3, ω = 0, which produces a static Brans-Dicke solution describing an inhomogeneity in a Minkowski background.
Other cases include: a) ω ≥ ω 0 and b) α = − √ 3/2. In both cases there are no wormhole throats and no apparent horizons and the Clifton-Mota-Barrow spacetime contains a naked singularity.

The Brans-Dicke solutions of Clifton, Mota, and Barrow
The next class of solutions of the Brans-Dicke field equations that we consider is that of Clifton, Mota, and Barrow given by the time-dependent and spherically symmetric line element [242] ds 2 = −e ν(̺) dt 2 + a 2 (t)e µ(̺) (d̺ 2 + ̺ 2 dΩ 2 ) , (4.242) where The matter source is a perfect fluid with energy density and pressure ρ (m) and P (m) and equation of state P (m) = (γ − 1) ρ (m) with γ =const. [242]. m is a mass parameter and α, φ 0 , a 0 , ρ and t 0 are related). ̺ is the isotropic radius related to the areal radius Schwarzschildr bỹ and α is real for ω < −2 and for ω > −3/2 but, for brevity, we require that ω 0 > −3/2 and β ≥ 0. The line element (4.242) reduces to the spatially flat FLRW one when m → 0. For γ = 2, setting ω = (γ − 2) −1 produces β = 0 and the spacetime becomes static (but the scalar still depends on time). The values γ = 2 and γ = 4/3 of the γ-parameter yield β = 1/2 and a(t) ∼ √ t independent of the parameter ω. To interpret the solution physically, the apparent horizons were studied in Ref. [244]. Again, the analysis proceeds by rewriting the line element using the areal radius and solving the equation g RR = 0 numerically. It was found that, according to the parameter values, several kinds of behaviours are possible. The "S-curve" familiar from the Husain-Martinez-Nuñez solution of GR is reproduced in a certain region of the parameter space, but different behaviours appear for other combinations of the parameters [244]. In certain regions of the parameter space, the Clifton-Mota-Barrow spacetime contains a naked singularity created with the universe and the metric and scalar field cannot be obtained as regular developments of Cauchy data. Later, this singularity is covered by black hole apparent horizons. In other regions of the parameter space, pairs of black hole and cosmological apparent horizons appear and bifurcate, or merge and disappear, as in the Husain-Martinez-Nuñez solution of GR and in the Clifton solution of f (R c c ) gravity [248,249]. Due to the larger parameter space involved with respect to what seen thus far, the Clifton-Mota-Barrow class of spacetimes exhibits the most varied and rich phenomenology of apparent horizons seen here (including some new one reported in [244]).

Clifton's solution of f (R c c ) gravity
Metric f (R c c ) gravity is described by the action where f (R c c ) is a non-linear function of its argument and S (matter) is the matter part of the action. As usual, R c c denotes the Ricci scalar of the metric g ab . The Jebsen-Birkhoff theorem of GR fails in these theories, as well as in scalar-tensor gravity, adding to the variety of spherical solutions [250]. Of particular interest are black holes in these higher order gravity theories.
A solution of vacuum f (R c c ) = (R c c ) 1+δ gravity was found by Clifton [248] and is given by the line element in terms of the isotropic radius, where Since Solar System tests require δ = (−1.1 ± 1.2)·10 −5 [248,251,252] and it must be f ′′ (R c c ) ≥ 0 for local stability [253,254,255,256], we assume 0 < δ < 10 −5 . Once δ is fixed, two classes of solutions exist, corresponding to the sign of C 2 qr. The line element (4.253) becomes FLRW if C 2 → 0. If δ → 0 (in which case the theory reduces to GR), the line element (4.253) reduces to the Schwarzschild one provided that C 2 qr > 0. In principle both positive and negative values of r are possible according to the sign of C 2 , but we impose that r > 0, C 2 > 0 and we take the positive root in the expression q = ± √ 1 − 2δ + 4δ 2 . Then q ≃ 1 − δ as δ → 0. The Clifton solution is conformal to the Fonarev spacetime [193], which is conformally static [257], hence it is also conformally static. Clifton's solution is dynamical and represents a central inhomogeneity in a spatially flat FLRW universe in vacuum (R c c ) 1+δ gravity. In the fourth order field equations of metric f (R c c ) gravity geometric terms play the role of an effective matter which invalidates the Jebsen-Birkhoff theorem of GR. Metric f (R c c ) gravity has an equivalent representation as an ω = 0 Brans-Dicke theory with a special potential for the scalar field degree of freedom f ′ (R c c ) [25,26]. The apparent horizons of the Clifton solution [248] were studied in Ref. [249]. First using the new radial coordinater dr and then using the areal radius the metric (4.253) is rewritten as Then the identities and the metric turns into the form Let us pass now to the new timet defined by to eliminate the time-radius cross-term, where F (t, R) is an integrating factor [249]. The line element is rewritten as By choosing the function β as the cross-term is eliminated and the metric becomes where H ≡ȧ/a and Using the second of these equations, the metric (4.268) becomes The equation g RR = 0 locating the apparent horizons is A q A 2 vanishes at r = C 2 , which describes the Schwarzschild horizon in the GR limit δ → 0. This is a spacetime singularity where the Ricci scalar R c c = 6(Ḣ+2H 2 ) diverges.
In the second case H 2 R 2 = C 2 A q+1 2 we have choosing the positive sign for an expanding universe. When δ → 0, this equation reduces to Two limits are now interesting: the first one is C 2 → 0, in which the central object disappears leaving FLRW space, r =r and R reduce to the comoving and the areal radius, respectively, and eq. (4.271) reduces to R c = 1/H, the radius of the cosmological horizon. The second limit of interest is δ → 0: the theory now degenerates into GR and eq. (4.271) reduces to A 2 = 0 or r = C 2 with H ≡ 0.
Eqs. (4.256) and (4.260) allow one to express the left hand side of eq. (4.271) as where x ≡ C 2 /r. The right hand side of (4.271) is and eq. (4.271) is simply ≃ 0.366. The radii R of the apparent horizons and the time t can be expressed in the parametric form , (4.276) using x as a parameter. This parametric representation of the horizons radii produces the same "S-curve" phenomenology of the Husain-Martinez-Nuñez solution [249].

Other solutions
Few other solutions which, judging from their apparent horizons, can be interpreted as black holes embedded in cosmological backgrounds are known in Brans-Dicke [258,259] and in other theories of gravity (see, e.g., those of Ref. [260] in Einstein-Gauss-Bonnet gravity, of Ref. [261] in higher order gravity, and of Ref. [262] in Lovelock gravity, and brane-world model solutions are known). Einstein-Gauss-Bonnet and Lovelock gravity, in particular, have been studied but they are appropriate in dimension D > 4 and here we have restricted ourselves to D = 4 -the zoo of black objects (Myers-Perry black holes, black strings, black rings, black Saturns, etc.) is much larger in higher dimension [263]. Moreover, most of the analytical solutions known are static or stationary. We did not mention here the variety of stringy and supergravity black holes (although most of them are stationary), which also deserve attention.
A problem is that, when trying to identify the physical mass of a non-asymptotically flat black hole in GR, we have made use of the Misner-Sharp-Hernandez mass and of the Kodama vector, defined in GR and for spherical symmetry. Thus far, beyond GR, the Misner-Sharp-Hernandez mass has been extended only to Einstein-Gauss-Bonnet gravity [264] (see [265] for an attempt to define it in cosmological spaces in f (R c c ) gravity) and is not available in alternative theories of gravity.

Conclusions
In this short review we have considered black holes with evolving horizons and, necessarily, we have given up the concept of event horizon and we have adopted the apparent and trapping horizon as its replacement. At present, this concept seems the best replacement and is widely used in practical (numerical) investigations, but it suffers from the drawbacks of being defined in a foliation-dependent way and of possibly becoming a timelike surface. The validity of the thermodynamics of apparent/trapping horizons needs to be studied better, as it seems that a quantum generalized second law does not hold for apparent and trapping horizons, while it holds for causal horizons, at least for 1 + 1 dilaton gravity (to which GR reduces in spherical symmetry) and for conformal vacua and coherent states [29].
We have focused on spacetimes representing, at least for part of their temporal extent, black holes embedded in cosmological backgrounds, of which a few solutions are known. We have considered only spherically symmetric spacetimes in D = 4 dimensions, with the goal to provide relatively simple explicit analytical examples that could be used as toy models for various investigations of the thermodynamics of apparent/trapping horizons, spherical accretion [60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75] (possibly of interest for primordial black holes), quantization of black hole areas [176,177,178], the issue of cosmological expansion versus local dynamics [35], and other theoretical topics. In addition, these solutions are of interest in themselves and, after all, they are not so simple. Even in the context of GR, where the McVittie spacetime has been known since 1933, this type of solution is poorly understood.
When asymptotic flatness is given up and non-stationary matter is allowed outside the black hole, there is no Jebsen-Birkhoff theorem and there is no "general" solution of the Einstein equations, as is instead the case for vacuum asymptotically flat solutions of GR (for this case, it is well known from a variety of studies that the Kerr-Newman black hole is the most general stationary, axisymmetric, asymptotically flat, vacuum solution [2,3,1]). The analytical solutions that can be given explicitly are probably pathological in this sense, rather than being generic (in some mathematical sense to be specified). The rather bizarre phenomenology of the apparent horizons in the solutions examined begs the question of whether they are even physically meaningful. They host naked singularities for part of the time and it is known that naked singularities form during gravitational collapse, but they do not seem to be generic and "typical" choices of initial data result in black holes rather than naked singularities [266,267]. It could well be that the solutions examined are exceptional rather than typical in the landscape of possible solutions. The "comoving mass" solution is a late-time attractor in the generalized McVittie class, but it is not a generic cosmological black hole solution with spherical symmetry [35,179]. No definitive statement can be made at present. Moreover, spacetimes hosting naked singularities cannot be obtained as the development of regular Cauchy data. Nevertheless, these spacetimes with time-evolving horizons tell us something about the theory of GR which may be important when trying to move beyond the paradigm of stationary and asymptotically flat vacuum black holes which has characterized black hole research thus far.
Things become even more uncertain when moving beyond the context of GR. To compound our ignorance, sometimes solutions of the field equations are discovered but no attempt is made to interpret them. Until not long ago, one had to be almost apologetic about working outside of GR, but things have changed and nowadays relevance for alternative theories of gravity is often used as a justification, or seen as an added value, for theoretical work. It is recognized that Einstein gravity will fail at some point and research in alternative theories is timely and important. We have only a handful of exact cosmological black holes in alternative gravity and we can only take a glimpse into this unexplored area. Some of the phenomenology of apparent horizons found in GR repeats itself in scalar-tensor gravity, but preliminary research has unveiled a wider range of behaviours.
At the moment, one could think of classifying cosmological black hole solutions of the field equations of a theory of gravity in two ways: 1. based on the type of matter filling the background FLRW universe (e.g. dust, general perfect fluid, imperfect fluid, or scalar field); 2. based on the technique used to generate the solution from a known "seed" metric (when applicable), for example conformally transforming the Schwarzschild or Kerr solution in some coordinate system [95], or applying a Kerr-Schild transformation [217,216,100]; 3. based on the phenomenology of the apparent horizons.
GR solutions with a perfect fluid include the (generalized) McVittie class (and its special case, the Schwarzschild-(anti-)de Sitter black hole). Solutions with a scalar field as a source include the Husain-Martinez-Nuñez, and the (generalized) Fonarev solutions. The phenomenology of apparent horizons distinguishes at least between the McVittie/Lemaître-Tolman-Bondi type with two appearing/disappearing (one black hole and one cosmological) apparent horizon, the Husain-Martinez-Nuñez/Clifton type with its "S-curve" phenomenology [183,249], and Lemaître-Tolman-Bondi spaces with a "double S-curve" [83], but other behaviours appear in the Clifton-Mota-Barrow solutions of Brans-Dicke theory due to the wider range of parameters.
In general, in both GR and alternative gravities, it is rare to find an explicit analytical expression of the proper (areal) radius of an apparent horizon of a dynamical cosmological black hole solution which could be used, for example, to investigate its thermodynamics. To the best of our knowledge such an expression is available only for the "comoving mass" subclass of the generalized McVittie class of solutions of GR. For other analytical solutions of the field equations the apparent horizon can only be obtained numerically or is given by an implicit analytical expression which is not very useful in practice.
The black hole solutions considered here represent eternal black holes which have not been created in a collapse process, but have existed forever or are created together with the universe in the Big Bang.
We can identify several other open problems: one is the possible generalization of the Misner-Sharp-Hernandez mass to non-spherical systems in GR, and its extension to alternative theories of gravity. Since the Misner-Sharp-Hernandez mass seems to be intimately connected with the Kodama vector, the problem can be seen as generalizing this vector to non-spherical systems and beyond GR. Another important problem concerns the thermodynamics of apparent/trapping horizons. There is a fairly large literature on this subject (both for black hole and cosmological horizons), which we cannot review here, and the "tunneling" method has been applied to the computation of the Hawking temperature of apparent horizons, including time-dependent ones (see the review [148]). However, contrary to the case of stationary black holes in which several independent calculational methods are available and produce the same result, for time-dependent apparent horizons only the Parikh-Wilczek "tunneling" method seem to deliver results and it is necessary to confirm these results with independent methods of calculation.
To conclude, cosmological black hole spacetimes as examples of time-evolving black holes disclose some of the complications of gravity and exhibit puzzling phenomenology. They can be useful for various areas of research in gravitational physics and quantum field theory in curved space but, unfortunately, we know too little about them. It is auspicable that in the near future more research is devoted to obtaining new solutions of this kind and, above all, understanding their physics.