1 Introduction

In 2016, the announcement from the Laser Interferometer Gravitational-Wave Observatory (LIGO) about the detection of the gravitational signals during the collision of two BHs sent ripples across the research community. The earthshaking news not only confirmed one of Einstein’s predictions in his theory of General Relativity (GR) also opened a new window to explore the dynamics of BH Physics and other space-time-warping phenomena in the universe [1]. After that, in 2019, the Event Horizon Telescope (EHT) discover the first high-resolution image of supermassive BHs in the center of M87 supergiant elliptical galaxy, providing the strong proof of the existence of BHs [2,3,4,5,6,7]. Theoretically, there is a dark central region known as BH surrounded by a ring-shaped lump of radiation, which is a so-called photon ring. A BH is commonly surrounded by a bright illuminated accretion matter. Several achievable observational features of BH were studied for a long time, and concluded that the BH shadow profile provide us significant information about the various accretion matters around the BH geometry [8].

It is well-known that when light emits from the accretion disk and passes by the BH, its trajectory is deflected due to intense gravitational lensing by the BHs. This phenomenon significantly reduces the observed intensity within the image plane, so-called critical curve (a special curve with both heuristics shadow and ring). The dark interior region of the critical curve is called the BH shadow. This discovery has significant ramifications. It not only enables us to verify the background geometry using powerful gravitational deflection but also serves as a test of the Kerr hypothesis on the formation of BHs compared to its numerous competitors [9]. The gravitational lensing allows us to analyze other stellar objects to evaluate the phenomena of shadow caster [10].

The effect of gravitational lensing in stellar objects surrounded by a luminous accretion flow matter may yield the complex sequence of numerous light-ray patterns and plays a significant role in the formation of BH profile. Generally, the trend of the critical curve is interpreted as: when the beam of light experienced by the distant observer traced backward, would have asymptotically led to a confined photon ring. The bound orbits for Schwarzschild BH occur at \(r\equiv 3M\) (in which M represents the mass of BH) and the apparent radius of the critical curve, i.e., impact parameter (IP) \(b\equiv 3\sqrt{3}M\approx 5.2M\) [11]. Further, the number of orbits (in the IP region) depend only on background geometry, and the optical image of BH shadow is not only its function but also depends upon the physical features of the jets and matter dynamics of the illuminating accretion flow around the compact objects.

By a comprehensive review, it is found that the contribution of lensing ring (light rays that cross the equatorial plane two times) to the total luminosity is appreciable, while the photon ring makes a small contribution [11]. Moreover, when the accretion matter is spherically symmetric, the observed angular size of the central dark region was found to be determined not only by the gravitational red-shift but also by the emission flow [12]. Zhu et al. [13] studied the BH shadow with the help of two different charged and slowly spinning BHs in the framework of Einstein-Æther gravity, and concluded that the influence of Æther parameter can change the shape of the shadow. Recently, the study of rotating BHs shadows is analyzed under isotropic fluid and cosmological constant in dark-matter era [14]. Cunha et al. [15] explored the BH shadow features and lensing effect, which is surrounded by thin accretion disk and found that any observer near the equatorial plane could receive different smears of the sky. The study of BH shadow is still debatable, whether the shadow profile is simply determined by spherically symmetric space-time metric or influenced by the flow matters around the BHs. In [16], authors found that the Schwarzschild BH shadow features under simple model spherical accretion and concluded that the edge location of the shadow does not depend upon the inner radius at which the accreting gas stops emission of radiating.

The observed appearance of BH shadow and its consequences are extensively examined in GR as well as in other extended gravitational theories [17,18,19,20,21]. Due to growing interest in exploring the mysterious aspects of accreting matters around the BHs it has still been a fascinating topic in the research community. Guo et al. [22] discussed the static and infalling accretion matter of a charged BH under a perfect fluid radiation field within Rastall gravity. They concluded that the BH shadow image with infalling accretion flow is dark as compared to the static one. Guo et al. [23] also discussed the observable characteristics of the charged BH under various accretion flow matters within Rastall gravity. They considered three positions of radiating accreting disk relative to the BH and compared the theoretically obtained two-dimensional profile of BH with EHT capture. Many other observational features of a BH and some characteristics are widely investigated in literature (for detail see references [24,25,26,27]).

The NC space-time provides a more realistic way to analyze the BH or wormhole solutions, as it provides an alternative way to quantum theory of gravity [28, 29]. Particularly, several methods are considered to implement the NC geometry in extended gravitational theories by modifying the matter source, i.e., Dirac delta function replaced by Gaussian or Lorentzian configuration in space-time structure [30,31,32,33]. However, the observed specific intensity and optical appearance of BHs by different accreting matters have rarely investigated in NC space-time. Recently, some authors discussed the BH shadows and different accreting emissions profiles in the framework of NC geometry. They concluded that the emission of profiles is gradually varies with the variations of NC operator [34]. In this scenario, Zeng et al. [35] studied the observational features of charged four-dimensional Gauss-Bonnet BH with clouds of strings and NC space-time. They analyzed the BH shadow and emission of specific intensities for different values of model parameters and found the darker region during the infalling movement as compared to static one due to Doppler effect.

In [36], Kiselev derived a new solution of Einstein’s field equations surrounded by various types of matter fields depending upon the BH state parameter \(\omega \). Kiselev derived the Schwarzchild and Reissner-Nordström-de Sitter like BHs solutions in the range of \(-1<\omega <-\frac{1}{3}\), which is lie in quintessence phase and \(\omega \) represent the ratio of pressure and energy density. Recently, several dynamical aspects of Kiselev BH have been discussed, i.e., accretion rate and critical points for the physical conditions [37], strong gravitational lensing [38], thermodynamical properties and phase transition [39]. Motivated from literature, our aim in this work is to interpret the observational features of charged Kiselev BH under NC geometry. Generally, we depict the NC charged Kiselev BH shadow image, which is surrounded by different accreting flow matters and observe the corresponding specific luminosity in the framework of different emission functions for some specific values of BH state parameters.

The plan of this paper is as follows: In Sect. 2, we provide the basic formulation of charged Kiselev BH and Lorentzian distribution and then observe the behavior of light trajectories near the BH for different values of involved model parameters. We observe the BH shadows image and photon rings with thin accretion disks flow models in Sect. 3. In Sect. 4, we analyze the shadows of BH with a static spherical accretion model, and the infalling spherical accretion model will be discussed in Sect. 5. Finally, we conclude the discussion in the last section.

2 Light deflection by charged Kiselev BH in NC geometry

The static, spherically symmetric metric of charged Kiselev BH is given as [36]

$$\begin{aligned} ds^{2}=-D(r)dt^{2}+\frac{dr^{2}}{D(r)}+d\Omega ^{2}, \end{aligned}$$
(1)

where

$$\begin{aligned} D(r)=1-\frac{2M}{r}-\frac{\lambda }{r^{3\omega +1}}+\frac{Q^{2}}{r^{2}}, \end{aligned}$$
(2)

\(d\Omega ^{2}=r^2(d\theta ^2+\sin ^2\theta d\varphi ^2)\), \(\lambda \) stands for Kiselev BH paremeter, Q denotes the electric charge and \(\omega \) is the equation of state parameter characterizing the matter radiation around the BH. For different classical/topological BH solutions, \(\omega \) can have many possible values (for details one can see reference [36]). We consider \(\omega =1/3\), which develops charged Kiselev BH surrounded by radiation field [40]. In this case, the function D(r) in Eq. (2) turns out to be

$$\begin{aligned} D(r)=1-\frac{2M}{r}-\frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}. \end{aligned}$$
(3)

The NC geometry eliminates the point-like structure in favor of smeared mass configuration in space-time objects. Under Lorentzian distribution, the modified mass-density takes the following form [41]

$$\begin{aligned} \rho _{\phi }=\frac{\sqrt{\phi }M}{\pi ^{\frac{3}{2}}(r^{2}+\pi \phi )^{2}}, \end{aligned}$$
(4)

in which \(\phi \) is the NC operator. The smeared mass configuration can be derived as [34]

$$\begin{aligned} \mathbb {M}_{\phi }= & {} 4\pi \int ^{r}_{0} r^{2}\rho _{\phi }(r)dr = \frac{2M}{\pi }\bigg [\tan ^{{-}1}\big (\frac{r}{\sqrt{\pi \phi }}\big ){-}\frac{r\sqrt{\pi \phi }}{r^{2}{+}\pi \phi }\bigg ],\nonumber \\= & {} M{-}\frac{4\sqrt{\phi }M}{\sqrt{\pi }r}{+}\mathcal {O}(\phi ^{\frac{3}{2}}). \end{aligned}$$
(5)

In this way, the charged Kiselev BH defined in Eq. (3) can be rewritten in the terms of smeared mass configuration as follows

$$\begin{aligned} D(r)=1-\frac{2\mathbb {M}_{\phi }}{r}-\frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}, \end{aligned}$$
(6)

which leads to modify the BH formulation as

$$\begin{aligned} D(r)=1-\frac{2M}{r}+\frac{8\sqrt{\phi }M}{\sqrt{\pi }r^{2}}-\frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}. \end{aligned}$$
(7)

The metric (7) reduces to Schwarzschild BH in the absence of NC operator, charge and Kiselev parameter. Note that in this analysis, we only focus the shadow/optical observation features with the solution of NC distribution of a static spherically symmetric charged Kiselev BH. The observational characteristics of Schwarzschild BH had already been discussed in [23], where authors considered the same accreting flow matter around the BH. The location of horizons can be obtained by putting \(D(r)=0\), as given below

$$\begin{aligned} r_{h}&=M+\sqrt{M^{2}-Q^{2}+\lambda -\frac{8M\sqrt{\phi }}{\sqrt{\pi }}},\nonumber \\ r_{c}&=M-\sqrt{M^{2}-Q^{2}+\lambda -\frac{8M\sqrt{\phi }}{\sqrt{\pi }}}, \end{aligned}$$
(8)

here \(r_{h}\) and \(r_{c}\) represent the event (or outer) horizon and the cosmological (or Cauchy) horizon, respectively. In order to investigate the motion of photons around the BH, we need to find how light trajectories travel in the space-time. So, the motion of photons satisfy the Euler–Lagrange equation:

$$\begin{aligned} \frac{d}{d\kappa }\bigg (\frac{\partial \mathcal {L}}{\partial \dot{x}^{\xi }}\bigg )=\frac{\partial \mathcal {L}}{\partial x^{\xi }}, \end{aligned}$$
(9)

where \(\kappa \) is an affine parameter, \(\dot{x}^{\xi }\) is the four velocity of the photon, “.” is the derivative with respect to the affine parameter and \(\mathcal {L}\) is Lagrangian, which can take the following form

$$\begin{aligned} \mathcal {L}{} & {} =\frac{1}{2}g_{\xi \zeta }\dot{x}^{\xi }\dot{x}^{\zeta }= \frac{1}{2}\bigg (-D(r)\dot{t}^{2}+\frac{\dot{r}^{2}}{D(r)}\nonumber \\{} & {} \quad +r^{2}(\dot{\theta }^{2}+\sin ^{2}\theta \dot{\varphi }^{2})\bigg ). \end{aligned}$$
(10)

The photon motions satisfy the null geodesic equation i.e., \(\mathcal {L}=0\). We only assume that the motion of photons take place in the equatorial plane, so we impose the initial conditions, \(\theta =\pi /2\) and \(\dot{\theta }=0\). Since the metric coefficients are independent of both time t and \(\varphi \), so we have two conserved quantities, which are \(E=-D(r)dt/d\kappa \) and \(L=r^{2}d\varphi /d\kappa \) represent the energy and angular momentum of photons, respectively. Using Eqs. (7), (8) and (9), the time, azimuthal angle and the radial components of four velocity satisfy the following equations:

$$\begin{aligned} \frac{dt}{d\kappa }= & {} \frac{1}{b\bigg (1-\frac{2M}{r}+\frac{8\sqrt{\phi }M}{\sqrt{\pi }r^{2}}- \frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}\bigg )}, \end{aligned}$$
(11)
$$\begin{aligned} \frac{d\varphi }{d\kappa }= & {} \pm \frac{1}{r^{2}}, \end{aligned}$$
(12)
$$\begin{aligned} \frac{dr}{d\kappa }= & {} \sqrt{\frac{1}{b^{2}}-\frac{1}{r^{2}}\bigg (1-\frac{2M}{r}+\frac{8\sqrt{\phi }M}{\sqrt{\pi }r^{2}}- \frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}\bigg )}, \end{aligned}$$
(13)

where the symbol “±”, corresponds to clock-wise (\(+\)) and anti-clock-wise (−) direction of motion of photons, b defines the IP, i.e., \(b\equiv \frac{L}{E}\). From Eq. (13), one can obtain the effective potential \(V_{eff}(r)\) as

$$\begin{aligned} \dot{r^{2}}=\frac{1}{b^{2}}-V_{eff}(r), \end{aligned}$$
(14)

where

$$\begin{aligned} V_{eff}(r)=\frac{1}{r^{2}}\bigg (1-\frac{2M}{r}+\frac{8\sqrt{\phi }M}{\sqrt{\pi }r^{2}}- \frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}\bigg ). \end{aligned}$$
(15)
Table 1 The numerical values of \(r_{h}\), \(b_{ph}\), and \(r_{ph}\), for some specific values of \(\phi \) with \(\lambda =0.1\), \(Q=0.1\) and \(M=1\)
Table 2 The numerical values of \(r_{h}\), \(b_{ph}\), and \(r_{ph}\), for some specific values of \(\lambda \) with \(\phi =0.01\), \(Q=0.1\) and \(M=1\)
Table 3 The numerical values of \(r_{h}\), \(b_{ph}\), and \(r_{ph}\), for some specific values of Q with \(\phi =0.01\), \(\lambda =0.1\) and \(M=1\)

At photon sphere, the motion of photons should have a definite value, i.e., \(\dot{r}=0\) and \(\ddot{r}=0\). Therefore, the radius of photon ring satisfies the following conditions

$$\begin{aligned} V_{eff}(r_{ph})=\frac{1}{b_{ph}^{2}}, \quad V'_{eff}(r_{ph})=0, \end{aligned}$$
(16)

where \(r_{ph}\) and \(b_{ph}\) represent the radius of photon sphere and critical IP, respectively. The quantities \(r_{ph}\) and \(b_{ph}\) can be re-organized as

$$\begin{aligned} r^{2}_{ph}=b^{2}_{ph}D(r), \quad 2b^{2}_{ph}D(r)^{2}-r^{3}_{ph}D'(r)=0. \end{aligned}$$
(17)

The numerical result of \(r_{h}\), \(r_{ph}\) and \(b_{ph}\) of the considering BH taking different model parameters are listed in Tables 1, 2 and 3. From Table 1, it is found that the increasing value of NC parameter \(\phi \), leads to decrease the value of \(r_{h}\), \(r_{ph}\) and \(b_{ph}\), which implies that the photon orbit shrinks inward to BH by increasing \(\phi \). In Table 2, it can be seen clearly, the numerical values of \(r_{h}\), \(r_{ph}\) and \(b_{ph}\) increase with the increasing values of \(\lambda \), whereas in Table 3, the numerical values of \(r_{h}\), \(r_{ph}\) and \(b_{ph}\) decrease with the increasing values of charge Q. In particular, we depict the behavior of \(V_{eff}(r)\) for some fixed values of model parameters in Fig. 1 (left panel). The effective potential vanishes at the event horizon, then it increases to a maximum and reaches the peak position of the photon sphere, \(r_{ph}\) subsequently drops down with respect to BH radius r. In region 1, where \(r>r_{ph}\), the light rays crossed the potential barrier and then reflected in the outward direction. The photons will fall into a singularity if they start the motion from \(r<r_{ph}\) (which do not hold in general).

In region 2 at \(b=b_{ph}\), the angular velocity is non-zero, since the light rays approach the photon sphere and revolve around the BH an infinite number of times. Here, the photon will move in circular intrinsically form and interpret the unstable behavior [11]. In region 3, the light trajectories will continue their motion toward the center and do not account for any potential limits at this location and hence, fall back into the BH singularity. The singularity carries no information about a BH, and it should not influence the creation of the shadow image, since the BH shadow is a space-time topological phenomenon. The evolution of the effective potential for some specific values of \(\phi \), \(\lambda \) and Q are depicted in Fig. 1. One can see that the effective potential is sensitive for all the model parameters.

By adopting the ray-tracing method, we can obtain the deflection of light rays near the BH. Using Eq. (13), we have

$$\begin{aligned} \frac{dr}{d\varphi }=\pm r^{2}\sqrt{\frac{1}{b^{2}}-\frac{1}{r^{2}}\bigg (1-\frac{2M}{r}+\frac{8\sqrt{\phi }M}{\sqrt{\pi }r^{2}}- \frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}\bigg )}.\nonumber \\ \end{aligned}$$
(18)

After that, we introduce a parameter \(\overline{u}=1/r\) and the above equation is modified as

$$\begin{aligned} \Psi (\overline{u}){} & {} = \frac{d\overline{u}}{d\varphi }\nonumber \\{} & {} = \sqrt{\frac{1}{b^{2}}{-}\overline{u}^{2}\bigg (1{-}2M\overline{u}{+} 8M\overline{u}^{2}\sqrt{\frac{\phi }{\pi }}{-}\lambda \overline{u}^{2}{+}(Q\overline{u})^{2}\bigg )}.\nonumber \\ \end{aligned}$$
(19)
Fig. 1
figure 1

The effective potential of the considered BH for some fixed values of state parameters (left panel), for different \(\phi \) and \(\lambda \) (middle panel) and for different Q (right panel) with \(M=1\)

Fig. 2
figure 2

The light rays for different specific values of model parameters in polar coordinates \((r,\varphi )\) with \(M=1\). The purple and red curves classified the \(b>b_{ph}\) and \(b<b_{ph}\) regions, respectively. The dashed green line denotes the radius of photon sphere and solid black disk represents the BH face orientation

The optical appearance of geodesics lines based on the roots of \(\Psi (\overline{u})=0\), such as when \(b>b_{ph}\), the light rays will be deflected at some radial position \(\overline{u}_{i}\), i.e., \(\Psi (\overline{u}_{i})=0\), and do not enter the BH. For the case, \(b<b_{ph}\), the light rays would not have \(\overline{u}_{i}\) and photons will continue to move towards the BH, and finally falls into it. The physical interpretation of light rays is depicted in Figs. 2 and 3, where the light rays are depicted for different model parameters. Clearly, one can see the influence of each parameter, as the size of the BH shadow is changing and change the deflection degree of each trajectories paths as well. It is worth mentioning that when \(b=b_{ph}\), the light rays revolve around the BH, which is an unstable circular photon orbit, and this conclusion is justified from the analysis of Fig. 1 (left panel, region 2).

Fig. 3
figure 3

The light rays for different specific values of model parameters in polar coordinates \((r,\varphi )\) with \(M=1\). The purple and red curves classified the \(b>b_{ph}\) and \(b<b_{ph}\) regions, respectively. The dashed green line denotes the radius of photon sphere and solid black disk represents the BH face orientation

Table 4 The numerical values of direct emission, lensing emission and photon ring emission are related to IP under specific values of model parameters. In the similar manner, we can analyze the numerical values of these emissions for other choices of model parameters. However, the value \(M=1\) is fixed for all cases
Fig. 4
figure 4

The light rays for different values of state parameters in Euclidean polar coordinates \((r,\varphi )\) with \(M=1\). Here, the spacings in IP are 1/5, 1/100 and 1/1000 for direct, lensing and photon ring emissions shown by red, blue and green bands, respectively. The dashed red line represents the photon orbit and solid black disk represents the BH face orientation

3 Thin disk accretion flow

Now, we study the effect of accretion disk on the observed BH shadow and photon rings. We deal this problem by considering optically and geometrically thin disk accretion scenario, considering that the disk radiate uniformly in the rest frame of static world-lines, the disk is located on the equatorial plane and the viewer is at the north pole.

3.1 Direct emission, lensing emission and photon ring emission

In [11], authors defined the important features of a BH surrounded by thin disk accretion flow matters and differentiate the lensing ring from photon ring near the BH. One can classify the total number of orbits made by different light trajectories on its path from its source to the observer as changed the azimuthal angle, i.e., \(n(b)\equiv \varphi /2\pi \), which is classified as [11]:

  • Case (I): \(n<3/4\) corresponds to direct emission, here the light rays crossed the equatorial plane just once.

  • Case (II): \(3/4<n<5/4\) corresponds to lensing ring emission, here the light rays crossed the equatorial plane twice.

  • Case (III): \(n>5/4\) corresponds to photon ring emission, here the light rays crossed the equatorial plane atleast three times.

In this manuscript, we show that the accretion flow scenario for three sets of model parameters. The range of b which is defined the numerical regions of emissions corresponding to three data sets are listed in Table 4. Combining with Table 4 and Fig. 4, we noticed that the influence of each parameter definitely change the behavior of photon ring and lensing ring in space-time, and the corresponding luminosity of these regions also varies significantly, however, the size of the BH which is shown as solid black disk is almost the same for these chosen sets of values.

3.2 Transfer functions

In the case of direct emission, the light drops on the front of the thin disk and for the lensed ring, the light rays revolve around the BH and then drop back off the thin disk, which leads to depicting the additional luminosity from the second crossing between accretion flow and light. While for the photon ring emission, the light again reaches the front side of the thin disk, which produces more brightness from the three times crossing between accretion flow and light. Therefore, the observed specific intensity on the observer’s plane is defined by the gravitational red-shift radiate density, where the influence of absorption and reflection of light can be ignored safely. And hence, whenever any light trajectory traced backward from the observer’s screen intersects the accretion disk plane it will gain additional brightness. From Liouvilles theorem, \(I_{\text {em}}(r)/(\nu _{\text {em}}^{f})^{3}\) is conserved along photon ray, where \(I_{\text {em}}(r)\) and \(\nu _{\text {em}}^{f}\), denote the emission-specific intensity and emission frequency, respectively. The observed intensity of photons with single frequency is defined as [42]

$$\begin{aligned} I_{\text {obs}}(r)= & {} D(r)^{\frac{3}{2}}I_{\text {em}}(r)\nonumber \\= & {} \bigg (1-\frac{2M}{r}+\frac{8\sqrt{\phi }M}{\sqrt{\pi }r^{2}}-\frac{\lambda }{r^{2}}+\frac{Q^{2}}{r^{2}}\bigg )^{\frac{3}{2}}I_{\text {em}}(r). \end{aligned}$$
(20)

Therefore, the total specific intensity can be calculated by integrating over different frequencies:

$$\begin{aligned} I_{obs}(r)= & {} \int I_{\text {obs}}(r)d\nu _{\text {obs}}^{f}(r) \nonumber \\= & {} \int D(r)^{2}I_{\text {em}}(r)d\nu _{\text {em}}^{f} =D(r)^{2} I_{\text {em}}(r), \end{aligned}$$
(21)

where \(I_{\text {em}}(r)\equiv \int I_{\text {em}}(r)d\nu _{\text {em}}^{f}\), represents the total emission of thin disk. The total received luminosity from each intersection can be defined as

$$\begin{aligned} I_{\text {obs}}(r)\equiv \sum \limits _{\begin{array}{c} n \end{array}} D(r)^{2}I_{\text {em}}|_{r}\equiv r_{n}(b), \end{aligned}$$
(22)

where \(r_{n}(b)\) represents the radial position of the \(n ^{th}\) intersection with the disk plane outside the horizon.

Fig. 5
figure 5

The first three transfer functions of the BH under some specific values of the state parameters with \(M=1\). All the panels from left to right are plotted for numerical values of parameters as defined in top row of Table 4

Figure 5 depicts the relationship between \(r_{n}(b)\) and IP function for different values of defined parameters. The demagnification factor dr/db, defines the slope of the transfer function. The red, blue and green lines represent \(n=1\), \(n=2\) and \(n=3\) transfer functions, respectively. When \(n=1\), it represents the direct image, so the direct image profile is just gravity red-shift source profile and its slope is nearly approaches to unity. The second (\(n=2\)) transfer function denotes the lensing ring emission. In this scenario, the obtained profile from the back side of the disk will be (de)magnified, because its slope will have relatively small value. As the value of IP increases, its slope also increases sharply to a very high value. Finally, photon ring corresponds to third (\(n=3\)) transfer function. In this situation, we can see an extremely (de)magnified profile from the front side of the disk as its slope approaches to infinite. Hence, its contribution can be ignored safely to the total flux of the BH profile.

Fig. 6
figure 6

The optical appearance of thin disk accretion with different emission flow profiles near the NC Charged Kiselev BH. The profiles are depicted for \(\phi =0.01, \lambda =0.1, Q=0.1\), and \(M=1\). The details can be found in the main text

Fig. 7
figure 7

The optical appearance of thin disk accretion with different emission flow profiles near the NC Charged Kiselev BH. The profiles are depicted for \(\phi =0.01, \lambda =0.2, Q=0.1\), and \(M=1\). The details can be found in the main text

3.3 Observational features of NC charged Kiselev BH

Here, we consider the three types of inner radii, where the accretion flow matter stops radiating, and observed the features of charged Kiselev BH NC geometry with the help of transfer function and the total observed luminosity equation. We consider three toy models listed below and consider innermost stable circular orbit (ISCO), where the accretion radiation stop. In the first model, we consider that \(I_{\text {em}}(r)\) is a decaying function with second power as

$$\begin{aligned} I^{a}_{\text {em}}(r)= {\left\{ \begin{array}{ll} {(\frac{1}{r-(r_{\text {ISCO}}-1)})^{2},\, \quad \,\text {if}\, r>r_{\text {ISCO}}} \\ 0,\quad \quad \quad \quad \quad \qquad ~~~\text {if} \,r\le r_{\text {ISCO}} \end{array}\right. } \end{aligned}$$
(23)

In the second model, the emission has a sharp peak position at \(r_{ph}\), having comparatively the same center and asymptotic trend as find in model 1. This is defined by

$$\begin{aligned} I^{b}_{\text {em}}(r)= {\left\{ \begin{array}{ll} {(\frac{1}{r-(r_{ph}-1)})^{3}, \quad \text {if}\, r>r_{ph}} \\ {0,\quad \quad \quad \quad \quad ~~~\text {if}\, r\le r_{ph}} \end{array}\right. } \end{aligned}$$
(24)

In the third model, the emission begins at the outside of the horizon and drops off more uniformly to zero as compare to previously defined two models, one can defined as

$$\begin{aligned} I^{c}_{\text {em}}(r)= {\left\{ \begin{array}{ll} {\frac{\frac{\pi }{2}-\tan ^{-1}(r-r_{\text {ISCO}}+1)}{\frac{\pi }{2}-\tan ^{-1}(r_{ph})}, \quad \text {if}\, r>r_{h}} \\ {0,\quad \quad \quad \quad \quad \quad \quad ~~~~~~\text {if}\, r\le r_{h}} \end{array}\right. } \end{aligned}$$
(25)

Now we proceed with the study of BH solution for three emission models and we consider two sets of parametric values, as we mentioned in the caption of Figs. 6 and 7. To this end, we depict the total emitted intensity with respect to radius r, the total observational intensity \(I_{obs}(r)\) as a function of IP, and two-dimensional density plots in celestial coordinates for each model, in left, middle and right panels, respectively. In the first row of Fig. 6 (left panel), the ISCO radius is \(r_{\text {ISCO}}\approx 5.92\), and the emission reaches at peak position when \(r\approx 5.92\), which denotes the radius of ISCO as emission stop position. In the middle panel of Fig. 6 (first row), the profile of direct emission obtains the peak position at \(b\approx 6.94\), the observed lensing ring lies within the small range \(5.08\sim 5.42\), and weak contribution in the observed intensity, while the photon ring is a highly narrow spike at \(b\approx 5.52\), and negligible contribution in total observed intensity. The two-dimensional density plots are showing the observational characteristics of the BH, where \(r_{\text {ISCO}}\) corresponds to the black disk boundary, the bright line in the black disk denotes the lensing ring, and the position of the photon ring shows the significantly weaker ring is moving continuously towards the interior as shown in Fig. 6 (first row right panel).

The second row of Fig. 6 (left panel) reflects the emission at the maximum position when \(r_{ph}\approx 2.73\), while in the middle panel, where direct emission is at peak position when \(b\approx 4.84\), and the limited range of lensing ring is \(4.81\sim 5.09\). Here, the photon ring overlapped in the lensed ring, when \(b\approx 4.81\). The optical appearance is provided in the right panel of the second row (Fig. 6), where the profile of direct emission is depicted with the merging of lensing and photon rings. The contribution of the lensing ring is highly appreciated and the contribution of the photon ring is barely visible.

It is observed that emission has a peak value at \(r_{h}\approx 1.79\) as shown in Fig. 6 (third row and left panel), and direct emission has peak value at \(b\approx 5.08\), and the range of lensed and photon ring emissions are limited to \(4.83\sim 5.24\) as shown in the middle panel. The photon ring is located at \(b\approx 5.2\). The two-dimensional optical appearance is shown in Fig. 6 (third row and right panel), where direct emission shows a narrow but a brighter ring, lensed and photon ring emission has the same reflection, so the latter can be ignored safely. In addition, we found that in all cases, the observed luminosity is generated by mainly direct emission, while the lensed ring has less contribution and the photon sphere has a negligible contribution in total brightness.

In Fig.  7, we depict the profiles for all the three emission flow models, they behaved the same as we analyzed in Fig. 6, the only difference is the numerical values and position of the photon and lensing ring, and intensities/descriptions as well.

4 Shadows of the NC charged Kiselev BH with a static spherical accretion

The observed intensity measured usually in \(\text {ergs}^{-1}\text {cm}^{-2}\text {str}^{-1}\text {Hz}^{-1}\) by the observer, expressed in [43, 44]:

$$\begin{aligned} I_{s}^{\text {obs}}(b)=\int g_{s}^{\text {obs}3}j(\nu ^{\text {obs}}_{em})dl_{pr}, \end{aligned}$$
(26)

where \(g_{s}^{\text {obs}}=\nu ^{\text {obs}}_{o}/\nu ^{\text {obs}}_{em}\) represents the red-shift factor, \(\nu ^{\text {obs}}_{o}\) is the observed frequency of the photon and \(\nu ^{\text {obs}}_{em}\) is the emission photon frequency, \(j(\nu ^{\text {obs}}_{em})\) denotes emission per unit volume, which is obtained in the rest frame of emitter and \(dl_{pr}\) is the infinitesimal proper length. From Eq. (1), the red-shift factor has value \(g_{s}^{\text {obs}}=(D(r))^{\frac{1}{2}}\) and radiation of emitted light is monochromatic having fixed frequency \(\nu _{fr}\), i.e.,

$$\begin{aligned} j(\nu ^{\text {obs}}_{em})\propto \frac{\delta (\nu _{em}-\nu _{fr})}{r^{2}}. \end{aligned}$$
(27)

In addition, the emissivity of light has \(1/r^{2}\) radial coordinate as argue in [44] and proper length in static frame is

$$\begin{aligned} dl_{pr}= & {} \sqrt{(D(r)^{-1}dr^{2}+r^{2}d\varphi ^{2})}\nonumber \\= & {} \sqrt{D(r)^{-1}+r^{2}\bigg (\frac{d\varphi }{dr} \bigg )^{2}}dr. \end{aligned}$$
(28)

Using Eqs. (26)-(28), one can derive the specific intensity as

$$\begin{aligned} I_{s}^{\text {obs}}(b)=\int \frac{D(r)^{3/2}}{r^{2}}\sqrt{D(r)^{-1}+r^{2}\bigg (\frac{d\varphi }{dr}\bigg )^{2}}dr. \end{aligned}$$
(29)

From Eq. (29), we investigate the shadow profile and corresponding intensity of charged Kiselev BH in NC space-time. One can observe that the intensity of light depends upon the trajectory of the light ray, which is related to IP. So, we analyze how the intensity varies with respect to IP. From Fig. 8 (top panel), we see that as b increases the intensity of light ramps first, and then moves towards the peak position at \(b_{ph}=4.85M\) (left panel), \(b_{ph}=4.95M\) (middle panel) and \(b_{ph}=4.82M\) and then drops smoothly to the lower values. This result is physically consistent with Figs. 1, 2 and 3. When \(b<b_{ph}\), the intensity of light, which is coming from accreting matter mostly absorbed by BH, and we see that the observed intensity is small. The light rays move in circles around the BH, when \(b=b_{ph}\), an unstable circular orbit of BH is shown in Fig. 1 (left panel) as a specific example of an effective potential.

For \(b>b_{ph}\), only the refracted light mainly participate in the observed intensity. Moreover, as the value of b increases, the refracted light starts to obtain the smaller values and thus, vanishes when \(b\rightarrow \infty \). In addition, we can also see from Tables 1, 2 and 3, each parameter will affect the observed intensity. The two-dimensional optical appearance is depicted in Fig. 8 (bottom panel), where the shadow is circularly symmetric around the BH, and the bright photon ring lies beyond the BH, which is known as the photon sphere. For different parameters, the radius of the photon sphere is listed in Tables 1, 2 and 3.

Fig. 8
figure 8

Profiles of the specific intensity \(I_{s}^{\text {obs}}(b)\), cast by a static spherical accretion mechanism, under some specific values of BH parameters as three instant examples with \(M=1\). The details can be found in the main text

5 Shadows of the NC charged Kiselev BH with an infalling spherical accretion

When we do not consider the angular velocity of BH, the accreting matter will only have radial velocity towards the BH. Here, Eq. (29) is applicable, but now red-shift factor can be written as

$$\begin{aligned} g=\frac{\mathcal {K}_{\eta }u_{0}^{\eta }}{\mathcal {K}_{\epsilon }u_{e}^{\epsilon }}, \end{aligned}$$
(30)

in which \(\mathcal {K}^{\tau }=\dot{x}_{\tau }\), \(u_{0}^{\tau }\) and \(u_{e}^{\tau }\) denote the four-velocity components of photon, distant observer and infalling accretion flow, respectively. The components of four-velocity of accretion flow are defined as

$$\begin{aligned} u_{e}^{t}=\frac{1}{D(r)}, \quad u_{e}^{r}=-(1-D(r))^{1/2}, \quad u_{e}^{\theta }\equiv u_{e}^{\varphi }\equiv 0. \end{aligned}$$
(31)
Fig. 9
figure 9

Profiles of the specific intensity \(I_{i}^{\text {obs}}(b)\), cast by a infalling spherical accretion mechanism, under some specific values of BH parameters as three instant examples with \(M=1\). The details can be found in the main text

With the help of Eqs. (11)–(13), one can obtain the \(\mathcal {K}_{t}\equiv 1/b\) (is a constant) and \(\mathcal {K}_{r}\) can be obtained from \(\mathcal {K}_{\alpha }\mathcal {K}^{\alpha }=0\), i.e.,

$$\begin{aligned} \frac{\mathcal {K}_{r}}{\mathcal {K}_{t}}=\pm \frac{1}{D(r)}\sqrt{1-\frac{b^{2}D(r)}{r^{2}}}, \end{aligned}$$
(32)

in which the sign “±” denotes for the motion of photons travels towards and move away from the BH, respectively. Using Eq. (32), the red-shift factor defined in Eq. (30) can be simplified as

$$\begin{aligned} g=\frac{1}{u_{e}^{t}+\mathcal {K}_{r}/\mathcal {K}_{e}u_{e}^{r}}. \end{aligned}$$
(33)

Now, the proper distance can be evaluated as

$$\begin{aligned} dl_{pr}=\mathcal {K}_{\epsilon }u_{e}^{\epsilon }d\kappa = \frac{\mathcal {K}_{t}}{g|\mathcal {K}_{r}|}dr, \end{aligned}$$
(34)

where \(\kappa \) denotes the affine parameter along the path of photon \(\epsilon \). Now, we consider the emission of observed intensity is monochromatic, so Eq. (27), is still applicable. So, the infalling accretion flow mechanism can be obtained as

$$\begin{aligned} I^{\text {obs}}_{i}(b)\propto \int g^{3}\mathcal {K}_{t}(r^{2}|\mathcal {K}_{r}|)^{-1}dr. \end{aligned}$$
(35)

With help of Eq. (35), we can analyze the shadow of BH numerically for different involved model parameters for infalling spherical accretion case. From Fig. 9 (top panel), it is found that for different involved parameters, the observed intensity is still close to the critical case \(b\sim b_{ph}\), and reaches the peak position first and then starts to drop with the increasing value of b at \(b= b_{ph}\) and then almost vanishes when \(b\rightarrow \infty \), behave similarly as observed in the static accretion case. The two-dimensional intensity of the image is shown in (see Fig. 9 (bottom panel)). We can observe that the radius of the shadow profile and the position of the photon sphere are the same as compared to the static accretion case. However, a major new feature is found that in the central region, the infalling accretion is darker as compared to static accretion, which is happening due to the influence of the Doppler effect and it is found to be darker nearer the event horizon.

6 Conclusions and discussion

A team of astronomers working on the EHT has captured the first image of the “shadow” and glowing surroundings the supermassive BH at the heart Milky Way galaxy. After that, the scientific team of EHT found an accretion disk structure around the BHs, formed by diffuse material in angular motion. The influence of gravitational and frictional forces changes the thermodynamical properties of the BH material, causing the emission or absorption of radiation from the BHs. In this paper, we study the shadow and observational features of charged Kiselev BH under NC geometry.

For this purpose, we examined some physical properties of the considered BH, such as light deflection, various emissions profiles, shadows and photon rings around the BH, and observed the intensities for static and infalling matter accretion profiles. Using the technique of null geodesic, the behavior of effective potential is shown in Fig. 1. In Fig. 1 (left panel), we depict a specific case of effective potential and mentioned three different regions and locations of the photon sphere as well as the event horizon \(r_{h}\). In the middle and right panels of Fig. 1, the trajectories of the effective potential reflect the influence of each parameter significantly, i.e., as we increase the value of NC operator \(\phi \), the effective potential reaches the peak value and then decreases as moves towards the boundary. The increasing values of \(\phi \) lead to a decrease in the radius of the event horizon, critical IP and radius of the photon sphere, which confirms that the location of the photon sphere was found to be near the BH.

For \(\lambda =0.2\), the effective potential also reaches the peak Value and attains lower values as we increase the values of \(\lambda \). With the increasing values of \(\lambda \), the event horizon’s radius, critical IP and the radius of the photon sphere also increase, so the photon ring is expanded outward the BH (as shown in the middle panel of Fig. 1). The charge also affects the trend of the effective potential. We analyzed the large influence of magnetic charge leads to a maximum ramp of effective potential with a minimum radius of BH. The increasing charge supports decreasing radius of the event horizon, shadow and critical IP values, and hence, the photon ring is also found near the BH with increasing values of Q. Moreover, we have also shown in Figs. 2 and 3, the ingoing and outgoing light trajectories and found that the curves are sensitive for all involved model parameters as they change the behavior of light rays and size of the BH solid disk.

With the help of the ray-tracing method, we study the light trajectories around the BH, according to the total number of orbits such as \(n(b)\equiv \varphi /2\pi \), as shown in Fig. 4. We classified these light trajectories according to the crossing with the equatorial plane, defined as direct emission, lensed ring, and photon ring. Via three transfer functions, as shown in Fig. 5, it is concluded that the photon ring is more highly (de)magnetized than the lensing ring and then compared with three different emission rings. One can see that the light trajectories of the photon ring are highly curved and crossed the disk at least three times, but its role in the total luminosity can be ignored due to its narrow area. The lensing ring occupied more area as compared to the photon ring, but its (de)magnetization is not too high in a specific range of IP, and hence makes a major contribution to the total flux and has a lesser contribution than direct emission. It is worth mentioning here that our result is consistent with the Schwarzschild solution [45], where direct emission always dominates the total observed intensity for all three model functions.

Next, we consider three toy models and observed the optically and geometrically thin disks accretion scenario for the considered BH in NC space-time geometry. From Figs. 6 and 7, one can observe that the emitted intensity from these toy models attain peak positions at ISCO for a distant observer. It is observed from the first model, there are two isolated spikes represent the photon and lensing ring emissions, which are created due to gravitational lensing. Therefore, the main contribution in the total flux is obtained from the direct emission, creating a wide rim, while the lensing ring appeared in the interior of the wide rim and photon ring lies in the innermost region, which is difficult to detect to naked eye.

In the second model, these three emissions are overlapped in a certain range of IP. We observed that in direct emission there are two peaks corresponding to the initial trajectory, which falls off at the position, where lensing and photon ring are merged, and IP develops the large ramp of spike and after that falls off sharply. The total optical evolution is that the lensing and photon rings are merged with the direct emission and the contribution of the lensing ring gain significant appreciation while the photon ring is difficult to detect.

The physical interpretation of the third model is that the directly observed region expanded all the way and moved down toward the outer horizon. In this case, the direct emission again gets contributions at the peak of IP from the ramp spike in the light ray orbit first and in the lensing ring very shortly after, before uniformly dropping off to 0. The two-dimensional intensity appearance is narrow and the luminosity expanded ring, which is made of direct lensing and photon ring emissions can be ignored safely.

Moreover, we also obtain the specific intensity of BH, which is observed by a static observer, in which accretion was considered to be either static or infalling as shown in Figs. 8 and 9, respectively. For both cases, we observed that each parameter has a significant role in the observed intensities. For the infalling accretion, the shadow of the BH turned out to be darker as compared to the static case, which is due to the influence of the Doppler effect.

Based on our analysis, we conclude that the optical appearance of charged Kiselev BH solution with NC distribution depends on the accretion flow matter and space-time geometry. We hope that this theoretical analysis will play a significant role to understand the BH shadow and luminosity in a more comprehensive way.