1 Introduction

Physics near the event horizon of a black hole has been one of the most fascinating topics of research ever since the discovery of General Relativity (GR) more than a century ago. Indeed, a singularity indicates the limits of applicability of a theory, and the event horizon, which cloaks the central singularity of a black hole might hold the key to understanding the most fundamental aspects of gravity. These include the elusive quantum aspects of gravity, which, many believe, should smoothen this singularity, which is often an inevitable end-state of gravitational collapse in GR.

Studies related to singularities and event horizons assume extreme importance and relevance given the fact that it is now commonly believed that the centers of most galaxies are inhabited by supermassive black holes. While such studies were purely of theoretical interest till a few years back, the advent of the Event Horizon Telescope (EHT) has proved to be a game changer. The phenomenal advances in observational studies near the event horizon of black holes has ushered in a new era where one can meaningfully compare theoretical results with EHT data. In fact, immediately after the EHT collaboration announced its first results on the radio source M\(87^*\) [1,2,3], a flurry of activities have started, one of the most important being the constraining of various solutions of GR and many other gravity theories that are either black holes or can mimic one. The standard way to do this is to compare the theoretically obtained shadow with the one reported by the EHT for M\(87^*\). Indeed, this has been shown to put stringent constraints on the parameters of the underlying theory [4,5,6,7,8].

In astrophysical scenarios, the spacetime geometry around a supermassive compact object is typically modelled by a static Schwarzschild or a stationary Kerr black hole solution of GR. While Birkhoff’s theorem guarantees that the former is the unique solution under the assumption of spherical symmetry in vacuum, the latter is perhaps more interesting as a rotating solution. Such a solution is axially symmetric, and thus more general than the idealized Schwarzschild solution. While the Kerr black hole has been extensively studied in the light of EHT data, a couple of its mimickers were recently studied in the same context in [4], where it was pointed out that such alternatives are indeed viable, within some particular range of parameters. Indeed, the nature of spacetime itself in the presence of strong gravity is still far from being understood in sufficient details. In such a situation, these studies assume importance due to the fact that these can rule out significant regions of parameter space so that one can focus on more meaningful theoretical models. In this context, one has to remember that as of now, experiments are indicative of the existence of dark, compact objects at galactic centres, which are supermassive in nature. While black holes offer the simplest and most elegant explanation for the nature of such objects, more exotic possibilities cannot be ruled out either, given the fact that there might be fundamentally new degrees of freedom in the strong gravity regime. Such exotic objects, have rapidly gained popularity of late, as it is now commonly believed that possible alternatives to black holes and its associated theorems might in fact be necessary. For a recent review on the relevance of these, see [9]. To wit: currently, the role of quantum physics in such regimes are not fully understood, and whether quantum effects can regularize the singularity and remove the horizon is a matter of debate [10]. It is therefore of interest to further the understanding of black hole mimickers and this is what we intend to do in this paper.

In astrophysical contexts, axially symmetric solutions of Einstein’s equations are interesting, as these are more general than idealized spherically symmetric cases. In the simplest scenario, one can envisage a static solution, which can have higher multipole terms in its potential, apart from the usual monopole term. Properties of such spacetimes are well studied in the literature, and exact solutions are known [11,12,13,14,15] (for more details, we refer the reader to the monographs [16] and [17]). Here, we will focus on the Zipoy–Voorhees spacetime [12, 13, 16,17,18,19], whose metric is popularly known as the \(\gamma \)-metric. These are singular and vacuum solutions of Einstein’s equations without an event horizon and are characterized by two parameters, which we call m and \(\gamma \). While m relates to the mass, \(\gamma \) has two special values, namely \(\gamma =0\) representing flat space and \(\gamma = 1\) being the Schwarzschild solution. For all other values of \(\gamma \), the space-time is axially symmetric, with \(\gamma \) being a measure of the deviation from spherical symmetry or from the Schwarzschild solution (for some recent works on various aspects of the \(\gamma \) metric, see [20, 21].)

Our purpose here is to check how far the \(\gamma \)-metric can act as a black hole mimicker. To this end, we first study the shadow cast by the \(\gamma \)-metric, and compare it with the EHT result. Our main observation is that while for \(\gamma < 1/2\), there will in principle be no shadow (this contradicts the recent results of [22, 23]), while the EHT data indicates that all values of \(\gamma \ge 1/2\) are allowed. This is verified by studying the gamma metric in the limit \(\gamma \rightarrow \infty \), where it coincides with the Chazy–Curzon solution in spherical coordinates [17]. We, therefore, come up with an interesting conclusion : the gamma metric, with \(\gamma \ge 1/2\) is an unconstrained axially symmetric solution of vacuum GR, consistent with current EHT data. In this sense, the shadow of the \(\gamma \)-metric adds to the list of black hole mimickers.

In the later part of this paper, we focus on geometrically thin accretion disk images for the \(\gamma \)-metric. We consider the Novikov-Thorne model, and study massive particles on the “equatorial” plane. The images of the disks are then numerically obtained using ray-tracing technique and compared with those of Schwarzschild black hole. Our conclusion here is that, in the presence of light rings, the accretion disk images of the \(\gamma \)-metric can very closely mimic those of the Schwarzschild black hole, but in the absence of such rings it is drastically different.

This paper is organized as follows. In the next Sect. 2, we review some properties of the \(\gamma \)-metric, and also the metric in the limit \(\gamma \rightarrow \infty \), in subsection 2.1. Next, in Sect. 3, we study the shadow cast by the \(\gamma \)-metric. This is followed by the Sect. 4, where we constrain the \(\gamma \)-metric using EHT results, substantiating our discussion above. Section 5 is devoted to the analysis of thin accretion disks in the \(\gamma \)-metric background, followed by the concluding Sect. 6. Throughout this paper, we will work in units where the gravitational constant G and the speed of light c is set to unity.

2 The \(\gamma \)-metric and its properties

We consider a particular solution of the Weyl class which consists of a family of static, axially symmetric vacuum solutions of Einstein’s equations. The solution is called as the Zipoy–Voorhees spacetime [12, 13, 16,17,18,19], whose metric is popularly known as the \(\gamma \)-metric. The Weyl class of metrics have a generic line element of the form [16,17,18]

$$\begin{aligned} ds^2 = -e^{2U}dt^2 +e^{-2U}\left( e^{2k}(d\rho ^2+dz^2) + \rho ^2d\phi ^2\right) , \end{aligned}$$
(1)

in cylindrical coordinates \((t,\rho ,\phi ,z)\). The \(\gamma \)-metric, for a particular solution of \(U(\rho ,z)\), \(k(\rho ,z)\) and after making a transformation from \((\rho ,z)\) to \((r,\theta )\) coordinates, is written as

$$\begin{aligned} ds^2&= g_{tt}dt^2+g_{rr}dr^2+g_{\theta \theta }d\theta ^2+g_{\phi \phi }d\phi ^2 \nonumber \\&= -A(r)dt^2 + \frac{1}{A(r)}\nonumber \\&\quad \times \left[ B(r,\theta )dr^2 + C(r,\theta ) d\theta ^2 +(r^2-2mr)\sin ^2\theta d\phi ^2 \right] , \end{aligned}$$
(2)

where the functions ABC are given as

$$\begin{aligned} A(r)&= \left( 1-\frac{2m}{r} \right) ^\gamma ,\nonumber \\ B(r,\theta )&=\left( \frac{r^2-2mr}{r^2-2mr+m^2\sin ^2\theta } \right) ^{\gamma ^2-1},\nonumber \\ C(r,\theta )&= \frac{(r^2-2mr)^{\gamma ^2}}{(r^2-2mr+m^2\sin ^2\theta )^{\gamma ^2-1}}. \end{aligned}$$
(3)

The transformation between the Weyl coordinates \((\rho ,z)\) of Eq. (1) and the Erez–Rosen coordinates \((r,\theta )\) of Eq. (2) is [16, 18]

$$\begin{aligned} \rho ^2=\left( r^2-2mr\right) \sin ^2\theta ,\quad z=\left( r-m\right) \cos \theta . \end{aligned}$$
(4)

The metric is characterized by two parameters m and \(\gamma \), where m is related to the mass, and \( \gamma \) measures deformation of the spacetime from spherical symmetry. The Arnowitt–Deser–Misner (ADM) mass of the spacetime is \( M=m\gamma \), and the corresponding quadrupole moment is \(Q=\frac{\gamma M^3}{3}(1-\gamma ^2)\) [18]. The monopole M and the quadrupole Q are the only independent components of multipole moments, as all higher order components can be expressed in terms of M and Q. When \( \gamma =0 \), the \(\gamma \)-metric represents the flat Minkowski spacetime, and for \( \gamma =1 \), it reduces to the spherically symmetric Schwarzschild solution. For all other values of \(\gamma \), it deforms the spacetime from spherically symmetric to axially symmetric, with \(\gamma <1 (\gamma >1)\) representing a prolate (oblate) spheroid.

The \(\gamma \)-metric is an interesting model to understand the directional behaviour of naked singularities [24, 25]. Since this is a vacuum solution, the Ricci scalar for the metric vanishes. Whereas the expression of the Kretschmann scalar, \({\mathcal {K}}=R_{\alpha \beta \delta \lambda }R^{\alpha \beta \delta \lambda }\) (\(R_{\alpha \beta \delta \lambda }\) is the Riemann curvature tensor), reads [25,26,27]

$$\begin{aligned}&{\mathcal {K}}\nonumber \\&=\frac{16m^2\gamma ^2}{r^{2(\gamma ^2+\gamma +1)}(r-2m)^{2(\gamma ^2-\gamma +1)} (r^2-2mr+m^2\sin ^2\theta )^{3-2\gamma ^2}}\nonumber \\&\quad F(r,\theta ), \end{aligned}$$
(5)

where

$$\begin{aligned} F(r,\theta )= & {} m^2\sin ^2\theta \left[ 3m\gamma \left( \gamma ^2 +1\right) (m-r)\right. \nonumber \\&\left. +\gamma ^2\left( 4m^2-6mr+3r^2\right) +m^2 \left( \gamma ^4+1\right) \right] \nonumber \\&+ 3r(\gamma m+m-r)^2(r-2m). \end{aligned}$$
(6)

It can be seen that the Kretschmann scalar diverges at \(r=0\) for all \(\gamma >0\). So, there is a curvature singularity at \(r=0\) for all \(\gamma >0\). Moreover, the nature of the surface \(r=2m=\frac{2M}{\gamma }\) is quite interesting. As for the Schwarzschild case, i.e., for \(\gamma =1\), this surface marks the location of the event horizon, which also represents the infinitely red-shifted surface for observers at spatial infinity. From Eq. (5), the expressions of the Kretschmann scalar along the polar axis (\( \theta =0,\pi \)) and \( \theta =\pi /2 \) become [25]

$$\begin{aligned} {\mathcal {K}}|_{\theta =0}&={\mathcal {K}}|_{\theta =\pi } =\frac{48m^2\gamma ^2(\gamma m+m-r)^2}{r^{4+2\gamma }(r-2m)^{4-2\gamma }}, \end{aligned}$$
(7)
$$\begin{aligned} {\mathcal {K}}|_{\theta =\frac{\pi }{2}}&=\frac{16m^2\gamma ^2}{r^{2\gamma ^2 +2\gamma +2}(r-2m)^{2\gamma ^2-2\gamma +2}(r-m)^{6-4\gamma ^2}}~H(r), \end{aligned}$$
(8)

where

$$\begin{aligned} H(r)&= m^4\left( \gamma ^4+3\gamma ^3+4\gamma ^2+3\gamma +1\right) \nonumber \\&\quad -3m^3r\left( \gamma ^3+4\gamma ^2+5\gamma +2\right) \nonumber \\&\quad + 3m^2r^2\left( 2\gamma ^2+6\gamma +5\right) -6mr^3(\gamma +2)+3r^4. \end{aligned}$$
(9)

From Eqs. (5), (7) and (8), we see that, for \(\gamma <2~(\ne 0,1)\), the Kretschmann scalar diverges at \(r=2m=\frac{2M}{\gamma }\) for the whole domain of \(\theta \), i.e., for \(0\le \theta \le \pi \). So there is a curvature singularity at \(r_s=2m=\frac{2M}{\gamma }\) for the mentioned range of \(\gamma \). On the other hand, for \(\gamma >2\) case, the Kretschmann scalar becomes zero at \(r=2m\) at \(\theta =0,\pi \) and diverges at \(\theta =\pi /2\). Therefore, no singularity exists on the polar axis in the range \(\gamma >2\), whereas along the equatorial direction, the curvature singularity does exist at \(r=r_s\) for all values of \(\gamma >0~(\ne 1)\) [22]. Moreover, the surface \(r=r_s\) also represents an infinitely red-shifted surface for any value of \(\gamma >0\) (\(\gamma =1\) being the Schwarzschild one) for an observer at spatial infinity. For example, if \(\nu \) is the frequency of a light ray emitted from a source at rest at a finite radial distance r, and \(\nu _{\infty }\) is the corresponding frequency of the same light ray received by an observer at rest at infinity, then we can write

$$\begin{aligned} \nu _{\infty }=\nu \sqrt{A(r)}=\nu \left( 1-\frac{2m}{r}\right) ^{\gamma /2}. \end{aligned}$$

As \(r \rightarrow r_s=2m\), we have \(\nu _{\infty } \rightarrow 0\), i.e., the light received by the observer is infinitely red-shifted, irrespective of the value of \( \gamma >0 \). Again from Eq. (5), we observe that \({\mathcal {K}}\rightarrow \infty \) when \(r^2-2mr+m^2\sin ^2\theta =0\) for the range \(0<\gamma \le \sqrt{3/2}\). Therefore, in addition to the singularities at \(r=0\) and \(r_s=2m\), there exists another singular surface inside the \(r_s=2m\) singularity, for this range of \(\gamma \). Since we are interested only in the region exterior to the outermost singular surface \(r_s=2m\), the singularities internal to this surface are immaterial to our analysis. For associated literature regarding the properties of this spacetime, for example the global structure, motion of test particles, accretion disk properties etc., we refer the reader to [18, 19, 26,27,28,29].

2.1 \(\gamma \rightarrow \infty \) limiting case

We have also considered the spacetime in the limiting case of \(\gamma \rightarrow \infty \), keeping the ADM mass M fixed and finite. The resulting spacetime corresponds to the Chazy–Curzon solution of GR in spherical coordinates [17, 30,31,32]. We call it as gamma-infinite (GI) spacetime for simplicity. The corresponding line element of the GI spacetime is

$$\begin{aligned} ds^2&= -\exp \left( -\frac{2M}{r}\right) dt^2 +\exp \left( \frac{2M}{r}\right) \nonumber \\&\quad \times \left[ \exp \left( -\frac{M^2\sin ^2\theta }{r^2}\right) \left( dr^2 + r^2 d\theta ^2\right) + r^2\sin ^2\theta d\phi ^2\right] . \end{aligned}$$
(10)

The Kretschmann scalar of this GI spacetime reads

$$\begin{aligned} \tilde{{\mathcal {K}}}&=\frac{16M^2\left[ 3(M-r)^2r^2 +M^2(M^2-3Mr+3r^2)\sin ^2\theta \right] }{r^{10}}\nonumber \\&\quad \times \exp \left[ \frac{2M(M\sin ^2\theta -2r)}{r^2}\right] . \end{aligned}$$
(11)

From the above expression, we see that

$$\begin{aligned} \left. \tilde{{\mathcal {K}}}\right| _{r\rightarrow 0}=\left\{ \begin{array}{ll} 0, &{} \text {if}\ ~~\theta =0,\pi \\ \infty , &{} \text {if}\ ~~\theta \ne 0,\pi \end{array} \right. \end{aligned}$$
(12)

Therefore, there is no curvature singularity at \( r=0 \) for \(\theta =0,\pi \), whereas for all other values of \(\theta \), curvature singularity exists at \( r=0 \). We shall also discuss the lensing properties of this spacetime in sequel.

3 Optical properties and shadows of the \(\gamma \)-metric

To analyze the gravitational lensing and shadows cast by the \(\gamma \)-metric, we need to study the motion of photons in this background. Due to time translational and azimuthal symmetries of this spacetime given by the metric in Eq. (2), we have two constants of motion along the null geodesics, namely, the energy E and the angular momentum L of photons about the axis of symmetry. Therefore, the t and \(\phi \) components of the geodesic equations for photons are given by

$$\begin{aligned} \dot{t}=-\frac{E}{g_{tt}}, \quad \dot{\phi }=\frac{L}{g_{\phi \phi }}. \end{aligned}$$
(13)

where an ‘overdot’ represents differentiation with respect to the affine parameter (unless otherwise specified). From the normalization of four velocities of photons (\(u^{\mu }u_{\mu }=0\)) along null geodesics, we obtain

$$\begin{aligned}&(-g_{tt})\left( g_{rr}\dot{r}^2+g_{\theta \theta }\dot{\theta }^2\right) +\left( \frac{-g_{tt}}{g_{\phi \phi }}\right) L^2=E^2~,\nonumber \\&\quad \text {i.e.,}~ (-g_{tt})\left( g_{rr}\dot{r}^2 +g_{\theta \theta }\dot{\theta }^2\right) +V_{\mathrm{eff}}=E^2 \end{aligned}$$
(14)

where \(V_{\mathrm{eff}}\) is the effective potential for photons and is given by

$$\begin{aligned} V_{\mathrm{eff}}=\frac{L^2}{r^2\sin ^2\theta }\left( 1-\frac{2m}{r}\right) ^{2\gamma -1} =\frac{L^2}{r^2\sin ^2\theta }\left( 1-\frac{2M}{\gamma r}\right) ^{2\gamma -1}. \end{aligned}$$
(15)

On the equatorial plane, \(\theta =\pi /2\) and \(\dot{\theta }=0\). Therefore, on this plane, we have

$$\begin{aligned} (-g_{tt})~g_{rr}\dot{r}^2+V_{\mathrm{eff}}^{\pi /2}=E^2 \end{aligned}$$
(16)

where the effective potential on the equatorial plane takes the form

$$\begin{aligned} V_{\mathrm{eff}}^{\pi /2}=\frac{L^2}{r^2}\left( 1-\frac{2m}{r}\right) ^{2\gamma -1} =\frac{L^2}{r^2}\left( 1-\frac{2M}{\gamma r}\right) ^{2\gamma -1}. \end{aligned}$$
(17)

If a photon arrives at the turning point (\(r_{\mathrm{tp}}\)) of its trajectory, we have \(\dot{r}|_{r=r_{\mathrm{tp}}}=0\), which from Eqs. (16) and (17) gives,

$$\begin{aligned} V_{\mathrm{eff}}^{\pi /2}~|_{r=r_{\mathrm{tp}}}&=\frac{L^2}{r_{\mathrm{tp}}^2}\left( 1-\frac{2M}{\gamma r_{\mathrm{tp}}}\right) ^{2\gamma -1}=E^2 \nonumber \\&\implies ~~ b(r_{\mathrm{tp}})=\frac{L}{E} =r_{\mathrm{tp}}\left( 1-\frac{2M}{\gamma r_{\mathrm{tp}}}\right) ^{(1-2\gamma )/2} \end{aligned}$$
(18)

where \( b(r_{\mathrm{tp}})\) is the impact parameter of a light ray. The maximum of the effective potential marks the position of the photon capture radius which is known as the light ring in an axially symmetric case, and the photon sphere in a spherically symmetric case. On the equatorial plane, the photon capture radius becomes

$$\begin{aligned} r_{\mathrm{ps}}=(2\gamma +1)m=2M+\frac{M}{\gamma }. \end{aligned}$$
(19)

Figure 1 shows the variation of the photon capture radius \(r_{\mathrm{ps}}\) on the equatorial plane and the singularity radius \(r_s\) as functions of \(\gamma \). Note that the photon capture radius exists only for \(\gamma \ge 1/2\). However, for \(\gamma <1/2\), it does not exist as \(r_{\mathrm{ps}}<r_s\) in this case. This is also clear from the effective potential \(V_{\mathrm{eff}}^{\pi /2}\). For \(\gamma >1/2\), \(V_{\mathrm{eff}}^{\pi /2}\) vanishes both at the singularity \(r_s\) and at the spatial infinity with a maximum in between marking the position of \(r_{\mathrm{ps}}\). However, for \(\gamma <1/2\), \(V_{\mathrm{eff}}^{\pi /2}\) diverges at the singularity \(r_s\) (see Fig. 2). Therefore, in this latter case, photons on the equatorial plane will always have turning points outside the singularity. This is also true for off-equatorial photon geodesics. The reason is that, if the geodesics governed by Eq. (16) always have turning points, then so do the ones governed by Eq. (14), as, in the limit \(r\rightarrow r_s\), the effective potential in Eq. (14) still diverges and the coefficient of \(\dot{\theta }^2\) vanishes. Therefore, in the \(\gamma <1/2\) case, as a photon with any non-zero impact parameter always has a turning point, there will be no capturing of photons at all, and hence no shadow will be produced in this case.

Fig. 1
figure 1

Variation of the photon capture radius \((r_{\mathrm{ps}})\) and the singularity \((r_s)\) on the equatorial plane as functions of \(\gamma \). The solid red line corresponds to \(r_{\mathrm{ps}}\) and the solid blue line represents \(r_s\). We have chosen \(M=1\) to obtain the plots. For \(\gamma <0.5\), we see that \(r_{\mathrm{ps}}<r_s\), i.e., the singular surface is outside the photon capture circle

Fig. 2
figure 2

A typical plot of \(V_{\mathrm{eff}}^{\pi /2}\) (in units of \(\bar{L}^2\)) as a function of r (in units of M) for \(\gamma =0.45\). A ray of light (indicated by the black line) comes from a distant source, hits the turning point at \(r=r_{\mathrm{tp}}\), turns back and escapes to infinity

To check the above results, we use our numerical ray-tracing techniques (with certain modifications) discussed in our previous work [33] and produce the images. We shoot photons with different impact parameters from a distant observer towards the lensing objects and integrate the geodesic equations backward in time. If a geodesic has a turning point and escapes to infinity for a given impact parameter, we assign bright point to this. On the other hand, if it is captured by the singularity, we assign dark point to the corresponding impact parameter. However, since at the singularity \(r_s=\frac{2M}{\gamma }\) some of the metric functions become infinite, we cannot exactly touch the surface of singularity due to numerical limitation. We shall have to take a region of tolerance \((\delta r)\) around the singularity. Therefore, we take the inner grid to be at \(r=r_s+\delta r\), where \(\delta r\) is very small. For \(\gamma \ge 1/2\), any photon hitting this surface will be captured by the singularity. Also, while performing ray-tracing, we consider piecewise step size in the affine parameter \(\lambda \). When the radial coordinate reaches below a predefined value \(r_s+0.3\), i.e., when \(r\le (r_s+0.3)\) along a geodesic, we decrease the step size to \(\Delta \lambda =10^{-5}\).

Figure 3 shows the ray-traced shadows of the \(\gamma \)-metric for different \(\gamma \). Note that, for \(\gamma \ge 1/2\), the singularity always casts shadows whose shapes and sizes will depend on the value of \(\gamma \). The shadow becomes more and more prolate as we decrease \(\gamma \). This is similar to the results obtained in [22, 23]. However, for \(\gamma =0.4\), we note that the dark region shrinks as we decrease \(\delta r\) and will vanish in principle in the limit \(\delta r\rightarrow 0\). Therefore, as explained before via a theoretical argument, \(\gamma <1/2\) case does not cast any shadow. This conclusion is in contradiction to the studies of shadows in [22, 23].

The shrinking of the dark region with decreasing \(\delta r\) can be understood in more detail from the analysis of geodesics on the equatorial plane. In Fig. 4, we have shown dependence of the impact parameter on turning point \(r_{tp}\) for equatorial geodesics. Note that the turning point lies very close to the singularity even when the impact parameter is \(\sim \) 3 or 4. If we take the turning point to be \(r_{tp}=r_s+0.001\), then corresponding impact parameters are 2.92 2.13 and 1.15 for \(\gamma =0.45\), 0.4 and 0.3, respectively. Therefore, if we take \(\delta r=0.001\), then this means that we are excluding all those photons having turning points \(r_{tp}\le (r_s+0.001)\), i.e., from the impact parameters space, we are excluding impact parameters \(b\le 2.92\), 2.13 and 1.15 for \(\gamma =0.45\), 0.4 and 0.3, respectively. Photons having impact parameters within the excluded region form dark spots, as they are excluded from the analysis. As a result, we are having the dark region for \(\gamma <1/2\). However, for a given \(\gamma <1/2\), decreasing \(\delta r\) means decreasing the excluded region from the impact parameter space. The excluded region and hence the dark region vanishes in the limit \(\delta r\rightarrow 0\). Therefore, we do not have any shadow in principle in the \(\gamma <1/2\) case.

Fig. 3
figure 3

Ray-traced shadows of the \(\gamma \)-metric. The black region will disappear for \(\gamma <1/2\) in the limit \(\delta r\rightarrow 0\), implying that there is no shadow in this case. The observer’s inclination angle is taken to be \(\theta _o=\pi /2\). The shapes of the dark regions for \(\gamma =0.4\) are numerical artefacts, see text

Fig. 4
figure 4

Impact parameter as a function of the turning point \(r_{tp}\) on the equatorial plane. The vertical dashed line corresponds to turning point \(r_{tp}=r_s+0.001\). The corresponding impact parameters for this turning point are 2.92 2.13 and 1.15 for \(\gamma =0.45\), 0.4 and 0.3, respectively. Here, we have considered \(M=1\)

Fig. 5
figure 5

Dependence of the diameter and deformation of the shadow on \(\gamma \) and \(\theta _o\). Here, \(M=(6.5\pm 0.7)\times 10^9 M_\odot \) and \(D=(16.8\pm 0.8)\) Mpc. The red dashed line shows the inclination angle \(\theta _o=17^\circ \). Note that, for this inclination angle, the diameter and the deformation of the shadow are consistent with the M87\(^*\) observation

We reiterate that the somewhat “strange” nature of the dark regions in Fig. 3e–h are numerical artefacts. Namely, when we take \(\delta r=0.05\) to produce the image in Fig. 3e, we are taking only those rays into consideration whose turning points \(r_{tp}\) lie on and outside the radius \(r_s+0.05\), and all the rays having turning points inside the radius \(r_s+0.05\) are excluded. Photons from these excluded rays produce the dark regions. When we decrease \(\delta r\) and take \(\delta r=10^{-5}\), we exclude fewer rays (only those having turning points inside \(r_s+10^{-5}\)), and more rays take part in image formation (Fig. 3h), thus reducing the size of dark region. Then, in the limit \(\delta r \rightarrow 0\) (when all the rays, in principle, take part in image formation), the shadow disappears. Finally, a word about numerically tracking small values of \(\delta r \sim 10^{-9}\) that we are taking, is in order. As in a numerical procedure, we cannot exactly touch the surface of the singularity (\(\delta r=0\)) as some of the metric functions blow up at the singularity \(r=r_s\), we must take a very small but finite region of tolerance \(\delta r\) around the singularity. Now, we need to track \(\delta r=10^{-9}\) between two successive numerical steps only for photons which reach close to the inner grid \(r_s+\delta r\). This happens only for photons which have turning point around \(r_s+\delta r\), with \(\dot{r}=0\) at the turning point. Therefore, as the photons reach close to their turning points around \(r_s+\delta r\), the magnitude of \(\dot{r}\) becomes small (close to zero). Additionally, with \(\Delta \lambda \) also small, the difference \(|\Delta r| \simeq |\dot{r}|~\Delta \lambda \) in r between two successive numerical steps becomes small. Note that this \(\Delta r\) must be smaller than \(\delta r\) for photons which reach close to \(r_s+\delta r\) in order to track small \(\delta r=10^{-9}\). That this is the case can be readily checked.

4 Constraining the \(\gamma \)-metric using the M87\(^*\) results

We now use the results from M87\(^*\) observation [1] and put possible constraint on the \(\gamma \)-metric. For this purpose, we use the average size of the shadow and its deformation from circularity. To this end, we first denote the horizontal and the vertical axes in the shadow plane by \(\alpha \) and \(\beta \) respectively and define an angle \(\phi \) between the \(\alpha \)-axis and the vector connecting the geometric centre (0, 0) of the shadow with a point \((\alpha ,\beta )\) on the boundary of a shadow. Therefore, the average radius \(R_{av}\) of the shadow is given by [7]

$$\begin{aligned} R_{av}^2=\frac{1}{2\pi }\int _{0}^{2\pi }l^2(\phi )\; d\phi , \end{aligned}$$
(20)

where \(l(\phi )=\sqrt{\alpha (\phi )^2+\beta (\phi )^2}\) and \(\phi =tan^{-1}(\beta (\phi )/\alpha (\phi ))\). Following [1], we define the deviation \(\Delta C\) from circularity as

$$\begin{aligned} \Delta C=\frac{1}{R_{av}}\sqrt{\frac{1}{2\pi }\int _{0}^{2\pi }(l(\phi )-R_{av})^2\; d\phi }. \end{aligned}$$
(21)

Note that \(\Delta C\) is the fractional RMS distance from the average radius of the shadow. According to EHT collaboration, the angular size of the observed shadow is \(\Delta \theta _{sh}=42\pm 3\) \(\mu \)as, and its deviation from circularity (\(\Delta C\)) is less than \(10\%\) [1]. Also, following [1], we take the distance to M87\(^*\) to be \(D=(16.8\pm 0.8)\) Mpc and the mass of the object to be \(M=(6.5\pm 0.7)\times 10^9 M_\odot \). These numbers imply that the diameter of the shadow in dimensionless units should be

$$\begin{aligned} \frac{d_{sh}}{M}=\frac{D\Delta \theta _{sh}}{M}=11.0\pm 1.5, \end{aligned}$$
(22)

where the errors have been added in quadrature. The above quantity must be equal to \(\frac{2R_{av}}{M}\), i.e., \(\frac{d_{sh}}{M}=\frac{2R_{av}}{M}\). In Fig. 5, we have shown the average diameter and the deviation from circularity of the shadow for different \(\gamma \) and the observers inclination angle \(\theta _o\). Here, we have taken \(0.5\le \gamma \le 1\). Note that the size of the shadow is consistent with the M87\(^*\) observations for all \(\gamma \) and \(\theta _o\). However, the deviation from circularity is more than \(10\%\), i.e., \(\Delta C>0.1\) over a small parameter region where \(\gamma \) is close to 0.5 and inclination is high simultaneously. If we restrict the inclination angle to be \(\theta _o=17^\circ \), which the jet axis makes to the line of sight [1], then both the shadow size and the deviation from circularity are consistent with the M87\(^*\) observations for all \(\gamma \) values considered in Fig. 5. For \(\gamma >1\), we have found that the deviation from circularity slowly increases with increasing \(\gamma \) for a given inclination angle. Therefore, for a given \(\theta _o\), the maximum deviation from circularity occur for the GI spacetime. However, with increasing \(\gamma \) (\(>1\)), the shadow size slowly increases from the Schwarzschild value \(6\sqrt{3}M\) for a given higher inclination angle and decreases for a given lower inclination angle. Therefore, for \(\gamma >1\), the maximum and minimum size respectively occur for the GI spacetime with \(\theta _o=90^\circ \) and \(\theta _o=0^\circ \). These maximum and minimum values are respectively given by \(d_{sh}/M\simeq 10.47\) and \(d_{sh}/M\simeq 10.14\) (see Fig. 6). For any other values of \(\gamma \) (\(>1\)) and \(\theta _o\), \(10.14\lesssim d_{sh}/M\lesssim 10.47\). For the GI spacetime with \(\theta _o=17^\circ \), \(d_{sh}/M\simeq 10.18\) and \(\Delta C\simeq 0.002\), which are consistent with the observation. Therefore, we find that, for the inclination angle of \(17^\circ \), the shadow of the \(\gamma \)-metric is always consistent with the M87\(^*\) observations for all \(\gamma \ge 1/2\).

Fig. 6
figure 6

The diameter of the shadow for the GI spacetime

The importance of the above result is that it provides a horizonless static, axially symmetric vacuum solution that can mimic a black hole for the above mentioned range of \(\gamma \). As we have mentioned in the introduction, such exotic black hole mimickers are of current relevance in the light of the EHT data, as the nature of quantum effects in the strong gravity regime is still elusive. In that context, what we have shown here is that the \(\gamma \) metric being consistent with current experimental data makes it an attractive solution for further study. Furthermore, any value of \(\gamma \ge 1/2\) being consistent with the data implies that the parameter is essentially unconstrained. More precise measurements might possibly lead to stronger constraints on \(\gamma \).

5 Accretion disks and their images

We now consider the properties and image of a geometrically thin accretion disk in the \(\gamma \)-metric given in Eq. (2). The disk consists of massive particles moving in stable circular timelike geodesicsFootnote 1 on the equatorial (\(\theta =\pi /2\)) plane and is described by the Novikov–Thorne model [34, 35]. Since the spacetime has time translational and azimuthal symmetries, we have two constants of motion along a timelike geodesic, namely, the specific energy \(\tilde{E}\) (energy per unit mass) and the specific angular momentum \(\tilde{L}\) about the axis of symmetry respectively. Therefore, the geodesic equations corresponding to t and \( \phi \) can be written as,

$$\begin{aligned} \dot{t}=-\frac{\tilde{E}}{g_{tt}}, \quad \dot{\phi }=\frac{\tilde{L}}{g_{\phi \phi }}. \end{aligned}$$
(23)

From the normalization of four velocity \(u^{\mu }\) (i.e. \(u^{\mu }u_{\mu }=-1\)) for massive particles, the radial geodesic equation on the equatorial plane can be written as

$$\begin{aligned} g_{rr}\dot{r}^2=-1-\frac{\tilde{E}^2g_{\phi \phi }+\tilde{L}^2g_{tt}}{g_{tt} g_{\phi \phi }}=\tilde{V}_{\mathrm{eff}}(r,\theta ,\tilde{E},\tilde{L}), \end{aligned}$$
(24)

where \(\tilde{V}_{eff}\) is the effective potential. A stable circular orbit is given by \(\tilde{V}_{\mathrm{eff}}=0\), \(\tilde{V}'_{\mathrm{eff}}=0\) and \(\tilde{V}_{\mathrm{eff}}''<0\). The first two conditions yield the following expressions for the specific energy and the specific angular momentum of the particles moving in the stable circular orbits:

$$\begin{aligned} \tilde{E}&=-\frac{g_{tt}}{\sqrt{-\left( g_{tt}+g_{\phi \phi }\Omega ^2\right) }}, \quad \tilde{L}=\frac{g_{\phi \phi }\Omega }{\sqrt{-\left( g_{tt}+g_{\phi \phi }\Omega ^2\right) }},\nonumber \\ \Omega&=\frac{d\phi }{dt}=\sqrt{-\frac{g_{tt,r}}{g_{\phi \phi ,r}}}, \end{aligned}$$
(25)

where \(\Omega = d\phi /dt\) is the angular velocity of the particles forming the disk. The flux of the electromagnetic radiation emitted from a radial position r of the disk is given by the standard formula [34, 35]

$$\begin{aligned} {\mathcal {F}}(r)=-\frac{\dot{M}}{4\pi \sqrt{-g}} \frac{\Omega '}{(\tilde{E}-\Omega \tilde{L})^2}\int _{r_{in}}^r (\tilde{E}-\Omega \tilde{L})\tilde{L}' dr, \end{aligned}$$
(26)

where \(\dot{M}=dM/dt\) is the mass accretion rate, \(r_{in}\) is the inner edge of the disk, and \(\sqrt{-g}=\sqrt{-g_{tt}g_{rr}g_{\phi \phi }}\) is the determinant of the metric on the equatorial plane. The marginally stable circular orbit is given by \(V_{\mathrm{eff}}''=0 \). This gives

$$\begin{aligned} r_{\pm }=\frac{1+3\gamma \pm \sqrt{5\gamma ^2-1}}{\gamma }M. \end{aligned}$$
(27)

Figure 7 shows the variations of \(r_{\pm }\), the radius of singularity \(r_s=\frac{2M}{\gamma }\) and the photon capture radius \(r_{ps}\) as functions of \(\gamma \). For \(0<\gamma < 1/\sqrt{5}\), \(r_{\pm }\) do not exist. Therefore, in this case, we have a single continuous disk with its inner edge at \(r_{in}=r_s\) and outer edge at some radius \(r>r_s\). For \(\gamma =1/\sqrt{5}\), \(r_{-}=r_{+}\), i.e., the outer edge of the inner disk coincides with the inner edge of the outer disk, thereby giving a single continuous disk. For \(1/\sqrt{5}< \gamma < 1/2\), \(r_s< r_{-}<r_{+}\), implying that there exist stable circular orbits in between the singularity \(r_s\) and \(r_{-}\), and also at radii greater than \(r_{+}\) with no stable circular orbits in between \(r_{-}\) and \(r_{+}\). Therefore, in this case, we have double disk configuration (two concentric disjoint disks). The inner disk extends from the singularity (i.e., \(r_{in}=r_{s}\)) to \(r_{-}\), and the outer disk extends from \(r_{in}=r_{+}\) to some radius \(r>r_{+}\). For \(1/2\le \gamma \le 1\), \(r_{-}\le r_s<r_{+}\). Therefore, in this case, \(r_{-}\) does not exist, and we have a single accretion disk with its inner edge at \(r_{in}=r_{+}\) and the outer edge at some radius \(r>r_{+}\). Although both the roots \(r_{\pm }\) are real for \(\gamma >1\) case, the angular momentum \(\tilde{L}\) of circular orbits with radii \(r\le r_{-}\) becomes imaginary. Therefore, in \(\gamma >1\) case, we have a single disk with its inner edge at \(r_{in}=r_{+}\). Note that, in the limit \(\gamma \rightarrow \infty \), i.e., for the GI spacetime, \(r_s=0\) and \(r_{\pm }=(3\pm \sqrt{5})M\).

Fig. 7
figure 7

Variations of \(r_{+}\) (blue), \(r_{-}\) (green), \(r_s\) (red) and \(r_{ps}\) (black) as functions of \(\gamma \). We have considered \(M=1\) to obtain the plots

Fig. 8
figure 8

The images of accretion disks in a Schwarzschild black hole and the \(\gamma \)-metric (ag). See the text for the disk configurations for different values of \(\gamma \). The outer edge of the outer disk is at \(r = 20M\), and the observer’s inclination angle is \(\theta _{o}=80^{\circ }\). The observer is placed at the radial coordinate \(r=10^4M\), which corresponds effectively to the asymptotic infinity. In order to get rid of the parameters M and \(\dot{M}\), we have normalized the fluxes by the maximum flux observed for the Schwarzschild black hole. Also, we have plotted square-root of the normalized flux for better illustration. The color bars show the values the flux. All spatial coordinates are in units of M. Here, we have taken \(\delta r=10^{-9}\) for \(\gamma < 1/2\) case. The central black strip in this case is a numerical artefact as was the case in Fig. 3e–h, and will disappear in the limit \(\delta r\rightarrow 0\)

We now use our numerical ray-tracing techniques (with certain modifications) discussed in our previous work [33] (see also [36]) and produce the images of accretion disks. The intensity maps of the images of accretion disks for the different cases discussed above are shown in Fig. 8. Note that, when there are light rings in the \(\gamma \)-metric, the image is very similar to that of a black hole, thereby mimicking a black hole. In the absence of light rings, however, the images for the \(\gamma \)-metric differ strikingly from that of the black hole, as shown by Fig. 8d–f. Note that the central black strips in Fig. 8d–f which look “strange” are only numerical artefacts, as was the case in Fig. 3e–h, and will disappear in the limit \(\delta r\rightarrow 0\), as was discussed towards the end of Sect. 3. Namely in this case also, this is a result of exclusion of rays that have turning points within the region of tolerance, and will disappear in the limit when this region is made to approach zero.

We should also point out that here, we have solved null geodesic equations numerically to generate the accretion disk images using the energy flux formalism, based on Novikov-Thorne model of accretion. Our numerical method is only used for ray-tracing. In a more realistic scenario, one needs to incorporate numerical methods to take into account the general relativistic magneto-hydrodynamic (GRMHD) phenomena around the accretions disk, which is not relevant for our analysis here. The main significance of our analysis lies in the fact that it does give us a good example of a black hole mimicker for a particular range of \(\gamma \).

Let us also mention that in Fig. 8, we have plotted the energy flux distributions coming from the accretion disk. The higher the emitted energy flux, the higher will be the brightness of the image. We can see that the brightness of images for \(\gamma <0.5\) cases is much higher as compared to the \(\gamma >0.5\) cases. From Fig. 7 and the discussion before it, we see that for \(\gamma \le 1/\sqrt{5}\), there exists a continuous disk with its inner edge at the singularity (\(r_{in}=r_s\)) and outer edge at some radius \(r>r_s\). Again, for \(1/\sqrt{5}<\gamma <0.5\), we have a double disk configuration (two concentric disjoint disks) where the inner edge of the inner disk coincides with the singular surface (\(r_{in}=r_s\)). On the contrary, for \(\gamma >0.5\) cases, we have single disk configurations where the singularity is covered by photon capture radius (\(r_{ps}\)) and this \(r_{ps}\) lies inside the inner edge of the accretion disk, i.e., \(r_s<r_{ps}<r_{in}\) (similar to black holes). Since, in this scenario, the inner edge of the disk lies outside the singular surface, the energy flux emitted from the inner edge of the disk has the maximum value which is finite. However, when the inner edge of the disk coincides with the singular surface, the emission from the vicinity of the naked singularity becomes unboundedly large. Similar features of diverging flux from accretion disks around naked singularities have been obtained in the literature, see, e.g., [36, 37]. Since the flux becomes unboundedly large, the brightness also becomes high for the images of \(\gamma =0.45\), \(\gamma =0.4\) or \(\gamma =0.3\) (Fig. 8d–f.)

6 Conclusions

The unprecedented advances in observational studies in the current era of the EHT have resulted in the exciting possibility of understanding the nature of singularities in GR, by comparing theory with experimental data. Motivated by the celebrated work of Hawking [38], the importance of such studies lies in the fact that these can ultimately throw light on the nature of strong gravity, effective at horizon scales [39]. It is by now well established that there might be fundamentally new physics at such scales, near a spacetime singularity, due to hitherto unknown quantum effects.

Currently, the best possible scenario to try and uncover such effects is via the study of black hole images. Such studies, which are abundant in the literature by now, point to the fact that horizonless objects might as well mimic black holes. We recall some of the current literature. Shadows of black holes with spherical accretion were studied in [40]. In [41], an accretion disc model of the Schwarzschild black hole was considered, and its shadow was studied in relation to the EHT data. The astrophysical nature of such shadows was studied in [42], where the authors showed that this was a fairly robust indicator of the associated physics, in the sense that it is not affected greatly by astrophysical parameters. A recent review and overview of aspects of shadows in black hole physics appears in [43]. An important aspect of such analyses is the fact that one can use the EHT data to constrain possible solutions of Einstein equations which include exotic compact objects other than black holes. Indeed, a plethora of activities have been reported in this direction in the recent past, and it is by now understood that several singular solutions, as well as horizonless compact objects are consistent with current data. The totality of such results will be pivotal in understanding the correct nature of the fundamental aspects of strong gravity. Currently, the debate is thus ongoing, and with improved data it is expected that stronger constraints will emerge on theoretical models that explain the data. It is in this context that the results presented here assume importance, as we have quantified a new black hole mimicker.

Here, we have carried out the analysis of shadows and accretion disk images of the Zipoy–Voorhees spacetime, characterized by the \(\gamma \)-metric. This is known to be an attractive singular vacuum solution of the Einstein equations without a horizon and with axial symmetry. We have shown that while for the parameter \(\gamma <1/2\), there will be no shadow, the \(\gamma \ge 1/2\) class is essentially unconstrained, i.e., all such values of \(\gamma \) are consistent with current EHT observations. We have further constructed thin accretion disk images in the \(\gamma \) metric background and shown that while for \(\gamma \ge 1/2\), these mimick the Schwarzschild black hole, they are dramatically different from the Schwarzschild case for \(\gamma <1/2\). Our results indicate that at the moment, the \(\gamma \) metric is a perfectly reasonable candidate for the central singularity at the galactic centre from EHT observations of M\(87^*\).

As we mentioned in the introduction, astrophysical black holes are often approximated by the static Schwarzschild or the stationary Kerr solution. The \(\gamma \)-metric on the other hand describes a static, axially-symmetric vacuum solution, and is attractive in its own rights. Since spherical symmetry in static vacuum solutions is by no means a fundamental criterion in black hole physics, our result that the \(\gamma \)-metric is perfectly admissible, and should be an interesting addition to the current literature. Further, we have shown how the thin accretion disk images here might be very similar to that of the Schwarzschild black hole in cases when there are light rings. These cases thus exemplify situations where the horizonless object might mimic a black hole. In the absence of any light ring, however, these images might be very different from the black hole case.

In continuation of this work, it should be interesting to compare other static axially symmetric solutions of GR with the current data from the EHT. It would also be interesting to apply the analysis in this paper to the most recent observational results from the EHT on Sgr\(\mathrm{A}^*\) (see [44] and its followup papers).