1 Introduction

A possibility of whether graviton can be massive has been a question for a long time in physics community since the birth of the theory of propagating massless graviton known as the Einstein’s general relativity. There were many attempts to generalize the Einstein’s theory to describe propagating massive graviton. In 1939, Fierz and Pauli first successfully constructed a linear theory of massive gravity [1] by adding a suitable mass term parameterized by graviton mass to the linearized version of general relativity. Though adding the mass term seems to break the general covariance, it can be recovered via the use of the so-called Stückelberg trick (see [2] for a review on how to apply the trick to the massive gravity model). However, thanks to the nonvanishing mass, the linear theories of massive gravity suffer from the so-called van Dam–Veltman–Zakharov discontinuity [3, 4] which corresponds to the fact that predictions from the massive gravity in the massless limit do not agree with the predictions from the general relativity. Even though a hope for fixing such discontinuity problems lies in a nonlinear regime of massive gravity [5], generic nonlinear generalization from the Fierz–Pauli linear massive gravity usually possesses the so-called Boulware-Deser ghost mode in the theory [6]. Fortunately, a special class of nonlinear theory of massive gravity has been successfully found by de Rham, Gabadadze, and Tolley, dubbed as dRGT massive gravity theory [7, 8]. The dRGT massive gravity has proved itself useful in cosmology since its cosmological solution incorporates the graviton mass as an effective cosmological constant [9, 10], hence the accelerating expansion of the universe can be explained in this context.

Four dimensional static spherically symmetric black holes in dRGT massive gravity with and without U(1) charge are constructed explicitly in [11]. The effective cosmological constant emerges naturally as a result of graviton mass term. Therefore, these black holes can be considered as generalized Schwarzschild/Reissner–Nordström black holes with positive and negative cosmological constant. The temperature, entropy and thermodynamic stability of these dRGT black holes are also explored [11]. Stability of the Schwarzschild-de Sitter black hole in dRGT massive gravity is studied in Ref. [12]. Stability of black holes in bi-gravity extension of the dRGT model is studied for certain specific conditions in Refs. [13, 14]. Beyond spherical symmetry, a cylindrically symmetric solution with an event horizon is possible even in standard general relativity. The cylindrical black hole (or black string) is proposed by Lemos with the presence of negative cosmological constant [15]. In higher dimensions, black branes and strings can exist and their stability can be explored using the QNMs [16]. As a generalization of Lemos’s black string, static solution with cylindrical symmetry is found in dRGT massive gravity [17]. Unlike the black string in general relativity, there exists dRGT black string solution even with a positive cosmological constant effectively contributed by the graviton mass.

On the other hand, linear stability of black hole closely relates to a study of field perturbation exterior to the black hole’s horizon. The evolution of field perturbation is dictated by gravitational interaction between black hole and the test-field itself. The existence of event horizon requires that the field mode near the black hole must satisfy with ingoing wave condition. Thus, the perturbation modes are not normal modes but rather quasinormal modes (QNMs). The quasinormal spectrum will be complex frequencies that can be uniquely determined by the black hole’s mass, charge and angular momentum. For instance, the QNMs of neutral and charged black holes in dRGT massive gravity are computed [18]. In our previous work, we investigate a rich structure of dRGT massive gravity by calculating the effect of massive charged scalar perturbation on charged dRGT black hole [19]. While stable modes are discovered in some parameter space, we also find that some dRGT black holes suffer from superradiant instability [19]. Interestingly, the QNMs can be categorized into three families, (a) the near event horizon mode, (b) the near cosmic horizon mode, (c) all-region mode [19]. Similar behaviours are found for the QNMs of black string in the dRGT background, the near-horizon modes of the event horizon and cosmic horizon for positive cosmological constant case are purely diffusive (or negative imaginary) but we will present more detailed analysis elsewhere. In this article, we will focus on the all-region QNMs of the black string. We refer interested readers to Refs. [20,21,22] for very nice reviews on black hole QNMs.

Black hole, black brane and black string in asymptotically anti de-Sitter (AdS) space with negative cosmological constant are also interesting from the holographic duality viewpoint. They are dual to the thermal field theory on the boundary whose temperature is identified with the Hawking temperature of the spacetime. Existence of horizon implies finite temperature of the dual field theory and simultaneously implies that any fluctuations in such background would inevitably produce disssipations into the horizon. Such dissipations are encoded in the QNMs of the fluctuations of the spacetime. Due to the holographic duality, these QNMs are exactly the poles of the retarded thermal Green function of the dual field theory. They represent the relaxation time of the perturbed thermal field system in the dual picture to reach the thermal equilibrium again. Massive graviton in the bulk induces the breaking of diffeomorphism-invariance on the boundary and originates the momentum dissipation in the dual hydrodynamics [23, 24].

For cylindrically symmetric set-up, toroidal, cylindrical and planar black holes in general relativity are stable against small perturbation in all channels, scalar, vector and tensor [25]. However, the stability of black string in the dRGT massive gravity has never been addressed. Therefore in this paper, we attempt to investigate this problem systematically. In Sect. 2, the basic set-up for constructing neutral black string in dRGT massive gravity is introduced. Then, we explore the effects of linear term \((\gamma )\) on the horizon structure of black strings in both positive and negative cosmological constant scenario. In Sect. 3, we examine massive scalar perturbation on the dRGT black string background. The Klein–Gordon equation on curve spacetime and boundary conditions are discussed. It is possible to obtain semi-analytical formula for the quasinormal frequency in a small-universe near-extremal limit, this is explored in Sect. 4. The QNMs of dRGT black string with positive cosmological constant are calculated and discussed in Sect. 5. High-momentum behaviour of QNMs with positive cosmological constant is numerically studied in Sect. 5.1 and analytically approximated using the first-order WKB method in Sect. 5.2. In Sect. 6, we present the QNMs of dRGT black string with negative cosmological constant. Finally we summarize our results in Sect. 7.

2 Formalism

We shall consider neutral black string in dRGT massive gravity couples with massive scalar field. This is described by the action [8] (with \(c=8\pi G=1\)),

$$\begin{aligned} S= & {} \frac{1}{2}\int d^4 x \sqrt{-g}\left[ R + m^2_g \mathcal {U}(g,\phi ^a)\right. \nonumber \\&\left. - g^{\mu \nu }\nabla _{\mu } \Phi \nabla _{\nu } \Phi - m_s^2\Phi ^2 \right] , \end{aligned}$$
(1)

where \(m_g\) and \(m_s\) can be considered as the graviton mass and the scalar field mass. The ghost-free massive graviton self-interacting potential is defined by,

$$\begin{aligned} \mathcal {U}(g,\phi ^a) = \mathcal {U}_2 + \alpha _3 \mathcal {U}_3 + \alpha _4 \mathcal {U}_4 \end{aligned}$$
(2)

where,

$$\begin{aligned} \alpha _3&= \frac{\alpha -1}{3}, \end{aligned}$$
(3)
$$\begin{aligned} \alpha _4&= \frac{\beta }{4} + \frac{1-\alpha }{12}. \end{aligned}$$
(4)
$$\begin{aligned} \mathcal {U}_2&= [\mathcal {K}]^2 - [\mathcal {K}^2], \end{aligned}$$
(5)
$$\begin{aligned} \mathcal {U}_3&= [\mathcal {K}]^3 - 3[\mathcal {K}][\mathcal {K}^2] + 2[\mathcal {K}^3], \end{aligned}$$
(6)
$$\begin{aligned} \mathcal {U}_4&= [\mathcal {K}]^4 - 6[\mathcal {K}]^2[\mathcal {K}^2] + 8[\mathcal {K}][\mathcal {K}^3] + 3[\mathcal {K}^2]^2 - 6[\mathcal {K}^4]. \end{aligned}$$
(7)

\(\alpha \) and \(\beta \) are free parameters. \(K^{\mu }_{\nu } = \delta ^{\mu }_{\nu } - \sqrt{g^{\mu \sigma }f_{ab}\partial _{\sigma }\phi ^a\partial _{\nu }\phi ^b}\). \([\mathcal {K}]=\mathcal {K}^{\mu }_{\mu }\) and \([\mathcal {K}^n]=(\mathcal {K}^n)^{\mu }_{\mu }\). We will work in unitary gauge for which the four Stückelberg fields take the form \(\phi ^a = x^{\mu }\delta ^a_{\mu }\). The fiducial metric \(f_{ab}\) is chosen to be,

$$\begin{aligned} f_{ab} = diag(0,0,\alpha _g^2 h_0^2,h_0^2), \end{aligned}$$
(8)

for the coordinates \((t,r,z,\varphi )\) where \(\alpha _g\) and \(h_0\) are arbitrary constants.

2.1 Field equations

1. Einstein’s equations

$$\begin{aligned} R_{\mu \nu } - \frac{1}{2} R g_{\mu \nu } = -m^2_g X_{\mu \nu } + T_{\mu \nu }^\Phi \end{aligned}$$
(9)

where the full expression for \(X_{\mu \nu }\) can be found in [11, 17]. The energy-momentum tensor is given by,

$$\begin{aligned} T_{\mu \nu }^{\Phi }= & {} \nabla _{\mu } \Phi \nabla _{\nu } \Phi - \frac{1}{2}g_{\mu \nu }\left( g^{\sigma \rho }\nabla _{\sigma } \Phi \nabla _{\rho } \Phi + m_s^2\Phi ^2 \right) \end{aligned}$$
(10)

2. Scalar-field equation

$$\begin{aligned} \nabla _a \nabla ^a \Phi - m_s^2 \Phi = 0 . \end{aligned}$$
(11)

2.2 Static solutions

The Einstein field equation (9) admits cylindrically symmetric solution or “dRGT black string” in the absence of scalar field, i.e., \(T^{\Phi }_{\mu \nu }=0\). The dRGT black string solution is defined as [17],

$$\begin{aligned} ds^{2}&= -f(r)dt^2 + f^{-1}dr^2 + r^2\alpha _g^2dz^2 + r^2 d\varphi ^2, \end{aligned}$$
(12)

where,

$$\begin{aligned}&f(r) = \alpha _m^2 r^2 - \frac{4M}{\alpha _g r} + \gamma r + \epsilon , \end{aligned}$$
(13)
$$\begin{aligned}&\alpha _m^2 = m^2_g(1+\alpha +\beta )\equiv -\frac{\Lambda }{3}, \end{aligned}$$
(14)
$$\begin{aligned}&\gamma = -\frac{\alpha _m^2h_0(1+2\alpha +3\beta )}{1+\alpha +\beta }, \end{aligned}$$
(15)
$$\begin{aligned}&\epsilon = \frac{\alpha _m^2h_0^2(\alpha +3\beta )}{1+\alpha +\beta }. \end{aligned}$$
(16)

M is mass per unit length in z direction of black string. One can clearly see that the graviton mass \(m_g\) naturally generates the effect of cosmological constant. With these definitions, \(\Lambda <0\) case associates with \(\alpha _m^2>0\) while \(\Lambda >0\) is obtained by letting \(\alpha _m^2<0\). This is the unique character of black string in dRGT massive gravity. Since there is no de-Sitter (dS) analogue of black string in standard four dimensional general relativity [15]. The presence of linear term \(\gamma \) and constant term \(\epsilon \) makes black string with positive \(\Lambda \) possible. As will be shown later, these two branches of black string solutions have different horizon and asymptotic structures. Therefore, to determine their QNMs, we need to consider each branch of solution separately. Note that, in the limit, \(\gamma =0,\epsilon =0,\alpha _g=\alpha _m\), this metric (12) reduces to four dimensional black string solution in general relativity found by Lemos [15].

We shall now investigate the root structure of dRGT black string with negative \(\alpha _m^2\). Generally speaking, the metric function (13) has three zeros. With a proper parameter choice, it is possible that all three roots will be real number. More specifically, for \(\alpha _m^2<0\) there are two positive real roots and one negative real root. These two positive roots will be treated as black string’s event horizon \(r_h\) and cosmological horizon \(r_c\), where \(r_h<r_c\). As an example, the roots structure of metric (13) is displayed in Fig. 1. In this plot, the black string mass M, cosmological constant \(\alpha _m^2\), \(\alpha _g\) and \(\epsilon \) are fixed to be \(M=1,\alpha _m^2=-\,0.01,\alpha _g=1,\epsilon =0\) respectively. Each curve represents different values of \(\gamma \). These black strings shown in this figure have two positive roots. The inner root is black string’s event horizon and the outer root is cosmic horizon. As \(\gamma \) increases, the cosmic horizon increases, but the event horizon slightly decreases. In our previous work [19], we showed that the charged dRGT black hole has different properties ranging from naked singularity to a black hole and extremal black hole. However, this is not the case for neutral black string with positive \(\Lambda \). Without the charge term, the metric (13) cannot develop a naked singularity or extremal scenario (for charged black string [17], naked singularity can exist just like in the black hole case). Note that throughout this work, \(\epsilon \) will be fixed to zero [19]. However, this setting prevents us from having black string solution with negative and vanishing \(\gamma \).

An example of black string solution with \(\alpha _m^2>0~(\Lambda <0)\) can be found in Fig. 2. In this plot, we set the black string mass, cosmological constant and other constant to be \(M=1,\alpha _m^2=4,\alpha _g=1,\epsilon =0\) and consider the effect of \(\gamma \) on the spacetime. As can be seen from the plot, the horizon structure dramatically changes from the positive \(\Lambda \) case. There is only one positive real root, i.e., black string’s event horizon, for vanishing \(\epsilon \) [26]. Unlike the positive \(\Lambda \) case, the black string exists with all possible value of \(\gamma \) (negative, zero and positive). We can see that the black string becomes large when \(\gamma \) is more negative. It is clear from the metric (13) that, the black string will always have at least one horizon. This is because the positive \(\alpha _m^2\) term will be dominant at large r. Without charges, one should not expect naked singularity scenario in this case as well.

Fig. 1
figure 1

The behaviour of metric function f(r) plotted against radial coordinate r for varying \(\gamma \) with \(\alpha _m^2=-\,0.01,M=1,\alpha _g=1,\epsilon =0\)

Fig. 2
figure 2

The behaviour of metric function f(r) plotted against radial coordinate r with differing \(\gamma \). For demonstration, we fix \(\alpha _m^2=4,M=1,\alpha _g=1,\epsilon =0\)

3 Scalar perturbations on dRGT black string

Now we shall consider a massive scalar field as a test field on dRGT black string background (12). This can be described by the Klein–Gordon equation (11). We make the following ansatz for the scalar field \(\Phi \),

$$\begin{aligned} \Phi&= \frac{\phi (r)}{r}e^{-i\omega t+ikz+i\lambda \varphi }, \end{aligned}$$
(17)

where \(\omega ,k,\lambda \) are the frequency, the wave number and angular quantum number of the scalar perturbation. Unlike spherically symmetric set-up where we normally obtain (decoupled) radial and angular part of the Klein–Gordon equation, in this case, however, we solely get a second order ordinary differential equation,

$$\begin{aligned} \frac{f}{r} \phi '' + \frac{f'}{r} \phi ' + \left( \frac{\omega ^2}{r f}-\frac{m_s^2}{r}-\frac{f'}{r^2}-\frac{\lambda ^2}{r^3}-\frac{k^2}{\alpha _g^2r^3}\right) \phi&= 0, \end{aligned}$$
(18)

where \(f'=df/dr\). This Eq. (18) can be recast into the Schrödinger-like equation,

$$\begin{aligned} \frac{d^2\phi }{dr_*^2} + \left[ \omega ^2 - V(r)\right] \phi&= 0, \end{aligned}$$
(19)

where the effective potential is given by,

$$\begin{aligned} V(r)&= f\left( m_s^2 + \frac{\lambda ^2}{r^2} +\frac{k^2}{\alpha _g^2 r^2} + \frac{f'}{r} \right) . \end{aligned}$$
(20)

We have introduced the tortoise coordinate \(r_*\),

$$\begin{aligned} \frac{dr_*}{dr}&= \frac{1}{f}, \end{aligned}$$
(21)

with \(-\infty<r_*<0\), where \(r_*\rightarrow -\infty \) near the event horizon and \(r_*\rightarrow 0\) at infinity if \(\alpha _m^2>0\). For \(\alpha _m^2<0\), \(r_*\rightarrow -\infty ,\infty \) as r approaches black string’s event horizon \(r_h\) and cosmological horizon \(r_c\), respectively. Note that, in the massless limit \(m_s\rightarrow 0\), Eq. (19) resembles those found in [25]. In Fig. 3, we display the effective potentials for positive and negative \(\alpha _m^2\). In both cases, the effective potentials increase with \(\gamma \). For \(\alpha _m^2<0\), the potential vanishes at the event horizon and cosmic horizon. In contrast, the potential increases with \(r_*\) and diverges as \(r_*\rightarrow 0\) when \(\alpha _m^2>0\).

Fig. 3
figure 3

The effective potential plotted against tortoise coordinate \(r_*\) for \(M=1,\alpha _g=1,\epsilon =0\) with varying \(\gamma \). (Left) \(\alpha _m^2=-0.25,\lambda =k=m_s=0\) (right) \(\alpha _m^2=0.003,\lambda =1,k=0,m_s=0.03\)

The main task of this paper is to solve the radial equation (19) with appropriate boundary conditions. If \(\alpha _m^2<0\), effective potential vanishes at both ends (see Fig. 3), i.e., the event and cosmic horizon. At these points, two independent solutions are \(\phi \sim e^{\pm i\omega r_*}\). With the presence of horizons, only ingoing modes are allowed at the event horizon \(\phi _{r_h}\sim e^{-i\omega r_*}\) and outgoing modes at the cosmic horizon \(\phi _{r_c}\sim e^{i\omega r_*}\). For \(\alpha _m^2>0\) case, the condition at the event horizon is still the same. However, since the effective potential diverges as r approaches infinity, we require the wave function to vanish at infinity. To satisfy these quasinormal boundary conditions, the associated quasinormal frequencies \(\omega \) will be discrete complex numbers. With the sign convention used in this work, scalar perturbations will be stable if Im\((\omega )<0\) (decaying) and unstable if Im\((\omega )>0\) (growing).

Since each branches of solution has different boundary condition, we shall apply two distinct numerical methods for calculating the QNMs. For positive \(\Lambda \) branch, we implement the so-called asymptotic iteration method (AIM). This method was firstly used to solve an eigenvalue problem [27]. An improve version of AIM was developed and applied to compute the QNMs of asymptotically flat Schwarzschild and Schwarzshild de-Sitter black holes [28]. Recently, the QNMs of black holes in dRGT massive gravity [18, 19] and Lifshitz gravity [31] are calculated using improved AIM technique. When \(\Lambda \) is negative, the spectral method (see Ref. [32] and references therein) is used to numerically calculate the QNMs. This method is recently used for computing the QNMs of charged dRGT black hole in asymptotically AdS space [19].

4 Near extremal dRGT black string: \(\alpha _m^2<0\)

One can obtain an analytic formula for the quasinormal frequency, if the effective potential in Eq. (19) can be written in some well-known form. This method was suggested by Ferrari and Mashhoon in 1984 [33]. It turns out that our potential can be put into those exact forms by considering the near extremal limit. By near extreme, we mean the spacetime which the event horizon \(r_h\) is very close to the cosmic horizon \(r_c\), i.e., \(\displaystyle {\frac{r_c-r_h}{r_c}\ll 1}\) [34]. Therefore in this section, we shall compute the QNMs of dRGT black strings with \(\alpha _m^2<0\) in the near extremal limit.

In general, the metric function (13) has three distinct roots at \(r_h\), \(r_c\) and \(r_0\). They can be written in the following form

$$\begin{aligned} f(r)&= \frac{\alpha _m^2}{r}(r-r_0)(r-r_h)(r-r_c). \end{aligned}$$
(22)

We define surface gravity associated to the event horizon

$$\begin{aligned} \kappa _h&\equiv \left. \frac{1}{2}\frac{df}{dr} \right| _{r=r_h}, \nonumber \\&= \frac{\alpha _m^2(r_h-r_0)(r_h-r_c)}{2r_h}. \end{aligned}$$
(23)

Three zeroes of f can be expressed in term of \(\alpha _m^2,\gamma ,\epsilon \) and M as follow

$$\begin{aligned}&r_0+ r_h + r_c = -\frac{\gamma }{\alpha _m^2}, \end{aligned}$$
(24)
$$\begin{aligned}&r_h^2 + r_hr_c + r_c^2 + \left( r_h + r_c\right) \frac{\gamma }{\alpha _m^2} = -\frac{\epsilon }{\alpha _m^2}, \end{aligned}$$
(25)
$$\begin{aligned}&r_hr_c\left( r_h+r_c\right) +\frac{\gamma }{\alpha _m^2}r_hr_c = -\frac{4M}{\alpha _g\alpha _m^2}. \end{aligned}$$
(26)

In the near extremal limit, we shall consider the following approximation

$$\begin{aligned}&r_0 \sim -\left[ 2r_h + \frac{\gamma }{\alpha _m^2}\right] , \end{aligned}$$
(27)
$$\begin{aligned}&\alpha _m^2 \sim -\left[ \frac{2\gamma }{3r_h}+\frac{\epsilon }{3r_h^2}\right] , \end{aligned}$$
(28)
$$\begin{aligned}&M \sim \frac{\alpha _gr_h}{6}\left[ \epsilon +\frac{\gamma r_h}{2}\right] , \end{aligned}$$
(29)
$$\begin{aligned}&\kappa _h \sim -\frac{(r_h-r_c)}{2r_h^2}\left[ \gamma r_h + \epsilon \right] . \end{aligned}$$
(30)

Since r varies between \(r_h\) and \(r_c\), thus \(r-r_0\sim r_h-r_0\). Therefore we can write down the metric function in the near extremal limit as

$$\begin{aligned} f&\sim \frac{\alpha _m^2}{r_h}(r_h-r_0)(r-r_h)(r-r_c), \nonumber \\&\sim -\frac{(r-r_h)(r-r_c)}{r_h^2}\left[ \gamma r_h + \epsilon \right] . \end{aligned}$$
(31)

With the definition of the tortoise coordinate (21), this allows us to express the radial coordinate r as a function of \(r_*\). Now it reads,

$$\begin{aligned} r&= \frac{r_c e^{2\kappa _h r_*}+ r_h}{1+e^{2\kappa _h r_*}}. \end{aligned}$$
(32)

Thus the metric function (31) can be expressed in the tortoise coordinate as

$$\begin{aligned} f&\sim \frac{\left( r_h\kappa _h\right) ^2}{\left( \gamma r_h+\epsilon \right) \cosh ^2(\kappa _h r_*)}. \end{aligned}$$
(33)

By the virtue of these approximations, the Klein-Gordon equation (19) in the near extremal limit becomes

$$\begin{aligned} \frac{d^2\phi }{dr_*^2} + \left[ \omega ^2 - \frac{V_0}{\cosh ^2(\kappa _h r_*)} \right] \phi&= 0, \end{aligned}$$
(34)

where

$$\begin{aligned} V_0&= \frac{\kappa _h^2}{\left( \gamma r_h + \epsilon \right) }\left[ m_s^2 r_h^2 + \lambda ^2 + \frac{k^2}{\alpha _g^2}\right] . \end{aligned}$$
(35)

This potential is the well-known Pöschl–Teller potential [37]. The bound states of this type of potential were studied extensively in [37]. The Pöschl–Teller technique was used to investigate the QNMs of near extremal Schwarzschild-dS [34, 38], d dimensional Reissner-Nordström-dS black hole [39] and black hole in massive gravity [40].

When one applies the boundary conditions of QNMs to (34), the associated quasinormal frequency \(\omega _n\) is given by [33, 34]

$$\begin{aligned} \omega _n&= \kappa _h\left[ \sqrt{\frac{V_0}{\kappa _h^2}-\frac{1}{4}} -\left( n+\frac{1}{2}\right) i\right] ,\quad n=0,1,2\ldots \end{aligned}$$
(36)

It is clear that the quasinormal frequency is quantized by an integer n. The most fundamental mode is characterized by \(n=0\), i.e., the most slowly decaying modes. As an example, the first three lowest modes of the quasinormal frequencies are shown in Table 1. In this table, we fix the mass and other parameters to be \(M=1,\alpha _g=1,\epsilon =0,\lambda =0,k=0\). At each mode n, the \(\omega _n\) is compared between the massless and massive scalar case. From the formula, it should not be surprised that the quasinormal frequencies are purely imaginary for the massless scalar case. It is clear from Table 1 that increasing mode number n affects only imaginary part of the frequencies. Comparing to the black hole case studied in Refs. [35, 36], the behaviour of black string QNMs is quite different. As the scalar mass increases, the imaginary part of the \(\omega \) is not decreasing to the quasi-resonant modes as found in the black holes.

We remark that, the formula (36) is reliable only when the black string’s event horizon is very close to the cosmological horizon, \(\displaystyle {\frac{r_c-r_h}{r_c}\ll 1}\). In the appropriate limit [41], Eq. (36) resembles those found in [34, 38].

Table 1 The QNMs of dRGT black strings with \(\alpha _m^2<0\) in the near extremal limit from formula (36). For demonstration purpose, we set \(M=1,\alpha _g=1,\epsilon =0,\lambda =0,k=0\)

5 QNMs of dRGT black strings with positive \(\Lambda \): improved AIM

In this section, we compute the QNMs of black strings in dRGT massive gravity with \(\Lambda >0\) by using improved asymptotic iteration method [28]. To do this, we shall introduce the radial coordinate transformation \(r=1/x\). The radial wave equation (19) becomes,

$$\begin{aligned} \phi '' + \frac{A'}{A}\phi ' + \left[ \frac{\omega ^2}{A^2} - \frac{1}{A}\left( \lambda ^2+\frac{k^2}{\alpha _g^2}+\frac{m_s^2}{x^2} - xf'\right) \right] \phi&= 0, \end{aligned}$$
(37)

where,

$$\begin{aligned} A(x)&= \alpha _m^2 - \frac{4Mx^3}{\alpha _g} + \gamma x + \epsilon x^2. \end{aligned}$$
(38)

Throughout this section, we shall denote derivative with respect to radial coordinate x by \(^\prime \). Next we define the following [18, 28, 42],

$$\begin{aligned} e^{i\omega r_*}&= \left( x-x_1\right) ^{\frac{i\omega }{2\kappa _1}}\left( x-x_2\right) ^{\frac{i\omega }{2\kappa _2}}\left( x-x_3\right) ^{\frac{i\omega }{2\kappa _3}}, \end{aligned}$$
(39)

where \(x_i\) is real root of f(r). The surface gravity associated with each horizon is given by \(\displaystyle {\kappa _i = \frac{1}{2}\frac{df}{dr}\Big |_{r=r_i}}\). In this section, the black string event horizon and cosmological horizon are defined by \(x_1\) and \(x_2\) whereas the negative real root is \(x_3\). To deal with the singularity at the cosmic horizon, we make the following substitution,

$$\begin{aligned} \phi&= e^{i\omega r_*}u(x). \end{aligned}$$
(40)

The wave equation (37) now reads,

$$\begin{aligned}&Au'' + \left( A'-2i\omega \right) u' - \left[ \lambda ^2+\frac{k^2}{\alpha _g^2}+\frac{4Mx}{\alpha _g}+\frac{\gamma }{x}\right. \nonumber \\&\quad \left. +\frac{\left( 2\alpha _m^2+m_s^2\right) }{x^2}\right] u =0. \end{aligned}$$
(41)

At the black string’s horizon, the divergent behaviour is scaled out by using,

$$\begin{aligned} u(x)&= \left( x-x_1\right) ^{-\frac{i\omega }{\kappa _1}}\chi (x). \end{aligned}$$
(42)

The final equation becomes,

$$\begin{aligned}&\chi ''(x) = \lambda _0(x,\omega )\chi '(x) + s_0(x,\omega )\chi (x), \end{aligned}$$
(43)
$$\begin{aligned}&\lambda _0 = \frac{i\alpha _g \omega }{M\left( x-x_1\right) \left( x_1-x_2\right) \left( x_1-x_3\right) } - \frac{\left( A'-2i\omega \right) }{A}, \end{aligned}$$
(44)
$$\begin{aligned}&s_0 = \frac{1}{A}\left[ \lambda ^2+\frac{k^2}{\alpha _g^2}+\frac{4Mx}{\alpha _g}+\frac{2\alpha _m^2+m_s^2+\gamma x}{x^2} \right. \nonumber \\&\quad \left. \phantom {\left[ \lambda ^2+\frac{k^2}{\alpha _g^2}+\frac{4Mx}{\alpha _g}+\frac{2\alpha _m^2+m_s^2+\gamma x}{x^2} \right. }+ \frac{\alpha _g\omega ^2}{M\left( x-x_1\right) \left( x_1-x_2\right) \left( x_1-x_3\right) }\right] \nonumber \\&\quad + \frac{i\alpha _g\omega f'}{2Mf\left( x-x_1\right) \left( x_1-x_2\right) \left( x_1-x_3\right) }\nonumber \\&\quad + \frac{\alpha _g\omega \left( 2iM\left( x-2x_1\right) \left( x_1-x_2\right) \left( x_1-x_3\right) +\alpha _g\omega x\right) }{4M^2x\left( x-x_1\right) ^2\left( x_1-x_2\right) ^2\left( x_1-x_3\right) ^2}. \end{aligned}$$
(45)

Equation (43) together with the coefficients (44), (45) form a core behind an improved AIM algorithm. We choose to omit the details of the numerical method used in this paper. For interested readers, we refer to Refs [18, 28] for more details on AIM.

Fig. 4
figure 4

The comparison between the QNMs in the near extremal limit obtained via formula (36) and numerical method of AIM

5.1 Results

The QNMs of the dRGT black strings with \(\alpha _m^2<0~(\Lambda >0)\) are calculated using Mathematica notebook modified from the one provided in Ref [43]. As a preliminary check, the results from Table 1 (massive scalar case) are reproduced by using an improved AIM. The first three lowest modes are then compared between the analytic expression (36) and AIM. The results are in a good agreement as displayed in Fig. 4. Note that we shall focus only on the modes with lowest imaginary parts for most situations.

Now we turn our attention to a non-extremal case. In Table 2, we compute the QNMs for massive scalar perturbations for various values of \(\gamma \). At each fixed \(\gamma \), the value of angular quantum number \(\lambda \) is varied from \(0-2\). We observe that as \(\gamma \) increases, the real part and imaginary part of the quasinormal frequencies also increase (in magnitude for Im\((\omega )\)). An increasing in Re\((\omega )\) corresponds to a more oscillation of the scalar wave which means it has more energy. From Im(\(\omega \)), it turns out that the mode with high energy decays faster than the low energy mode. This is because as \(\gamma \) gets larger gravity becomes more attractive. Stronger gravity implies that the scalar wave tends to decay faster. We also see that, Re\((\omega )\) increases as angular quantum number increases while Im\((\omega )\) slightly decreases (in magnitude). In the fourth column, the quasinormal frequencies from the third-order WKB approximation (see appendix of Ref. [19] for details, see also a more accurate sixth-order WKB method in Ref. [29] and the semi-analytic method in Ref. [30]) are shown. The results show a good agreement between the two methods (AIM and WKB).

Table 2 The QNMs for massive scalar perturbations of a neutral black string for \(M=1,\alpha _g=1,\epsilon =0,\alpha _m^2=-0.25,m_s=0.5,k=0\)
Table 3 The QNMs for massive scalar perturbations of a neutral black string for \(M=1,\alpha _g=1,\epsilon =0,\gamma =1.9,m_s=0.5,k=0\)
Table 4 The QNMs for scalar perturbations of a neutral black string for \(M=1,\alpha _g=1,\epsilon =0,\alpha _m^2=-0.16,\gamma =1.5,\lambda =1\)

The effect of \(\alpha _m^2\) on the QNMs is investigated in Table 3. Let’s remark that as \(\alpha _m^2\) becomes more negative the difference between black string’s event horizon \(x_1\) and cosmic horizon \(x_2\) is smaller. When \(\alpha _m^2\) is larger (in magnitude), the real part of \(\omega \) is smaller for each fixed \(\lambda \). However, the behaviour of the imaginary part of \(\omega \) is not trivial as \(\alpha _m^2\) increases. In our chosen parameters, Im\((\omega )\) increases as \(\alpha _m^2\) changes from \(-\,0.16\) to \(-\,0.25\) and then decreases as \(\alpha _m^2\) gets larger. For fixed \(\lambda =2\), however, Im\((\omega )\) gets lower with increasing \(\alpha _m^2\). We also find that at each value of \(\alpha _m^2\), Re\((\omega )\) increases with \(\lambda \) whereas Im\((\omega )\) decreases marginally as happens in Table 2 (except when \(\alpha _{m}^{2}=-\,0.16\)). We find a close agreement between the results from AIM and the WKB approximation.

In Table 4, we consider how scalar mass affects the quasinormal frequencies in the presence of the massive gravity effects. The values of scalar mass \(m_s\) are varied from \(0-0.5\). In this table, the angular quantum number \(\lambda \) is set to be unity and we investigate the effect of wave number k on the QNMs instead. One remarkable feature is that there exists the QNMs with zero real part and negative imaginary part, i.e. these modes are purely decaying. For massless scalar perturbations, the purely imaginary frequencies \(\omega \) are found for all values of k considered here. At \(k=2\), a QNM with nonzero real part appears. However, nonzero real parts are found more abundantly when the scalar mass \(m_s=0.25\) and 0.5. The energy naturally goes up along with increasing k while the modes with larger wave number tend to have a shorter lifetime. When comparing the results with WKB approximation, we find a good match between the two methods. However, only the QNMs with nonzero real parts can be obtained from the WKB method. This is natural since WKB requires the positive-energy solution to exist for certain region around the minimum of the effective potential. For example in the massless case, WKB can only find the non-diffusive QNMs for \(k=2\). Interestingly, AIM also yields this value as the first excited mode \((n=1)\), while the most fundamental mode (obtained from AIM) \((n=0)\) is diffusive as shown in the table.

In order to explore the effect of \(m_s\) on the peculiar behaviour of the QNMs found in Table 4, we plot Re\((\omega )\) and Im\((\omega )\) with smaller increment on the scalar mass to elucidate transit of the lowest mode from the purely imaginary to the values with real parts. This is displayed in Fig. 5 and the values are listed in Table 5. The frequencies plotted in this figure are modes with the lowest imaginary part (in magnitude) for both the purely imaginary mode and the non-zero-real-part mode. The transition of the lowest QNM from purely imaginary mode to non-zero-real-part mode occurs between \(m_s=0.3-0.35\). A closer investigation reveals that the transitions appear at \(m_s=0.33,0.30,0.21\) for \(k=0,1,2\) respectively.

Fig. 5
figure 5

The plot of real and imaginary part of quasinormal frequencies for \(M=1,\alpha _g=1,\epsilon =0,\alpha _m^2=-\,0.16,\gamma =1.5,\lambda =1,k=0\) for varying \(m_s\). A subplot is displayed to clarify the effect of \(m_s\) on the non-zero-real-part modes

Table 5 The quasi-normal frequencies of the two lowest modes for scalar perturbations of a neutral black string for \(M=1,\alpha _g=1,\epsilon =0,\alpha _m^2=-0.16,\gamma =1.5,\lambda =1\) with varying \(m_s\)

The transition is part of a bigger picture of the QNM pattern which consists of two distinct branches found by our numerical method. The first one is purely imaginary mode where Im\((\omega )\) increases with the scalar field mass. The second branch exists with non-zero real part. In this case, \(m_s\) affects the quasinormal frequency such that as the scalar mass increases Re\((\omega )\) also increases but Im\((\omega )\) decreases. For low \(m_{s}\), the non-zero-real-part QNMs have larger \(|\mathrm{Im}(\omega )|\) than the purely imaginary mode, and they exist as higher modes. At transition \(m_{s}\), the imaginary parts of these modes become less than the purely imaginary mode and yet both branches always coexist within the QNM chart.

It is known that the perturbation by a massive field often leads to purely real frequency or quasi-resonance modes [35, 36]. In dRGT model with positive \(\Lambda \) (negative \(\alpha _{m}^{2}\)), the effective potential approaches zero in the asymptotic regions, \(V\rightarrow 0\) as \(r_{*}\rightarrow \pm \infty \) as we can see from Fig. 3. This is generic since \(V\sim f(r)\sim 0\) at both \(r_{h}\) and \(r_{c}\). As a result, the effect of the scalar mass vanishes in the asymptotic regions and there is no quasi-resonance modes. Our numerical results confirm non-existence of the quasi-resonance modes for positive \(\Lambda \) case.

Fig. 6
figure 6

The behaviour of QNMs at large \(\lambda ,k\) for \(M=1,\alpha _g=1,\epsilon =0,\alpha _m^2=-0.04,\gamma =0.4,m_s=0.2\). (Left) Plot of Re\((\omega )\) as a function of \(\lambda \), (right) plot of Im\((\omega )\) as a function of \(\lambda \), a small window shows the behaviour of \(k=0\) curve for small \(\lambda \)

One remarkable feature of the QNMs of neutral dRGT black string is the behaviour of the quasinormal frequency at high momentum, \(\lambda ,k \gg 1\). This is illustrated in Fig. 6. We set \(M=1,\alpha _g=1,\epsilon =0,\alpha _m^2=-0.04,\gamma =0.4,m_s=0.2\) and vary \(\lambda ,k\) from 0 to 200. It is apparent that the scalar wave becomes more energetic as \(\lambda \) and/or k get larger. Moreover, the energy of the wave changes rapidly at small wave number k. From the plot of Im\((\omega )\), the higher \(\lambda ,k\) the wave decays slower. Similar to the real part, the effect of \(\lambda \) on the QNMs is significant only at the low k. We expect the similar pattern if \(\lambda \) and k are interchanged in the above figure. This is because from the equation of motion (19), the roles of \(\lambda \) and k to the QNMs are identical up to \(1/\alpha _g^2\) factor.

5.2 WKB approximation for high-momentum modes

From the first-order WKB approximation, it is possible to estimate the quasinormal frequencies in the large-momentum limit. Note that the following calculations were also done in the Schwarzschild case [33] and the Schwarzschild-dS case [38]. The first-order WKB approximation yields the following relation [44, 45],

$$\begin{aligned} \frac{iQ(r_1)}{\sqrt{2Q^{(2*)}(r_1)}}=\left( n+\frac{1}{2}\right) , \end{aligned}$$
(46)

where \(Q^{(2*)}\) is the second derivative of Q with respect to the tortoise coordinate \(r^*\), \(r_1\) minimizes \(Q(r^*)\) (in terms of the tortoise coordinate), and Q is defined as,

$$\begin{aligned} Q\equiv \omega ^2-f\left( m^2_s + \frac{\ell ^2}{r^2}+\frac{f'}{r}\right) , \end{aligned}$$
(47)

where \(\ell ^2\equiv \lambda ^2+\displaystyle {\frac{k^2}{\alpha _g^2}}\). With \(\epsilon =0,\) in the large-momentum limit Q is then given by,

$$\begin{aligned} Q\approx \omega ^2-\frac{f\ell ^2}{r^2} = \omega ^2 - \frac{\ell ^2}{r^2}\left( \alpha ^2_m r^2 - \frac{4M}{\alpha _g r} + \gamma r\right) . \end{aligned}$$
(48)

A local minimum of \(Q(r^*)\) can be found as a solution to the following,

$$\begin{aligned} Q^{(1*)}(r_1)&= \left. \frac{dQ}{dr^*}\right| _{r_1} = \left. f\frac{dQ}{dr}\right| _{r_1}, \\&=-\frac{f(r_{1})\ell ^2}{r^4_1}\left( \frac{12M}{\alpha _g}-\gamma r^2_1\right) , \end{aligned}$$

The local minimum of \(Q(r^*)\) implies either,

$$\begin{aligned} f(r_{1}) = 0 \text { or } \frac{12M}{\alpha _g}-\gamma r^2_1=0. \end{aligned}$$

Since we demand that \(r_1\) minimizes Q, then we should have \(Q^{(2*)}(r_1)>0\). Thus \(Q^{(2*)}\) is now calculated,

$$\begin{aligned} Q^{(2*)}&= \frac{d^2Q}{(dr^*)^2} = f\frac{d}{dr}\left( f\frac{dQ}{dr}\right) , \\&=-f\ell ^2\left[ f'\left( \frac{12M}{\alpha _g r^4}-\frac{\gamma }{r^2}\right) +f\left( -\frac{48M}{\alpha _g r^5}+\frac{2\gamma }{r^3}\right) \right] . \end{aligned}$$

Since the first choice \(f(r_{1})=0\) yields vanishing \(Q^{(2*)}(r_1)\), therefore we must choose the second choice which implies \(r_1=\sqrt{\frac{12M}{\alpha _g \gamma }}\). Note that \(r_1\) is not defined when \(\gamma <0\). Thus we obtain,

$$\begin{aligned}&Q(r_1) = \omega ^2-\ell ^2\left( \alpha ^2_m+\frac{8M}{\alpha _g r^3_1}\right) , \\&Q^{(2*)}(r_1)=\frac{2f^2_1\ell ^2\gamma }{r^3_1}. \end{aligned}$$

By requiring that \(r_1>r_h\), we can solve (46) for \(\omega \) as,

$$\begin{aligned} \omega&=\ell \sqrt{\alpha ^2_m + \frac{8M}{\alpha _g r^3_1}}\sqrt{1-\frac{2i\sqrt{\gamma r_1}\left( n+\frac{1}{2}\right) }{\ell }}, \nonumber \\&= \sqrt{\alpha ^2_m + \frac{8M}{\alpha _g r^3_1}}\left( \ell -i\sqrt{\gamma r_1}\left( n+\frac{1}{2}\right) \right) \nonumber \\&\quad +\mathcal {O}\left( \frac{1}{\ell }\right) ~\text {for large }\ell . \end{aligned}$$
(49)

It is interesting to note that the above approximation (49) has an interesting interpretation in the context of an unstable circular orbit of a null geodesic around a massive object [46]. Since the local minimum of Q, which is \(r_1\), has been proven to also be a radius of a circular null geodesic, then, according to Ref. [46], the real part of the QNMs in (49) can be interpreted as multiples of an angular frequency of the corresponding circular null geodesic while the imaginary part indicates instability timescale of the orbit.

In the large-momentum limit, the imaginary part of the lowest QNMs (\(n=0\)) approaches the asymptotic value,

$$\begin{aligned} \text {Im}~\omega (\ell \rightarrow \infty )=-\frac{1}{2}\left( \alpha ^2_m\sqrt{\frac{12M\gamma }{\alpha _{g}}} + \frac{2}{3}\gamma ^{2}\right) ^{1/2}. \end{aligned}$$
(50)

Numerical results we found in Fig. 6 match the values given by (49) for \(n=0\) with errors less than \(0.1\%\) for large momentum.

6 QNMs of dRGT black strings with negative \(\Lambda \): spectral method

In this section, we will explore the QNMs of scalar perturbation in the black string spacetime in the dRGT model with negative \(\Lambda \). Geometrically, this is not an asymptotically global AdS space since we have cylinder instead of sphere. The boundary topology is cylindrical but otherwise the holographic nature of the space remains the same as the AdS space. We can understand the general aspects of the perturbations by considering (19) in two limits, near-horizon and far-away regions. For the near-horizon region, the equation of motion becomes,

$$\begin{aligned} \frac{d^{2}\phi }{dr_{*}^{2}}=-\omega ^{2}\phi , \end{aligned}$$
(51)

so the solution is simply \(\phi (r)\sim e^{\pm i\omega r_{*}}\). For QNMs calculation, we take the infalling wave \(\phi \sim e^{-i\omega r_{*}}\), the perturbation is leaking into the black string horizon. For the far-away region since \(f(r)\simeq \alpha _{m}^{2}r^{2}=-\Lambda r^{2}/3\), Eq. (19) takes the form,

$$\begin{aligned} \frac{d^{2}\phi }{dr_{*}^{2}}= -\frac{\Lambda }{3}r^{2}\left( m^{2}_{s}-\frac{2\Lambda }{3}\right) \phi , \end{aligned}$$
(52)

and has the power-law solution \(\phi \sim r^{\alpha }\) where

$$\begin{aligned} \alpha = \frac{1}{2}\left( -1\pm \sqrt{9-\frac{12 m_{s}^{2}}{\Lambda }}\right) . \end{aligned}$$
(53)

Note that \(m^{2}_{s}\) is bounded by \(m_{s}^{2}\le 3\Lambda /4\). In order for the field to vanish at infinity for the plus sign solution of (53), \(m_{s}^{2}\ge 2\Lambda /3\) is required otherwise it will become a non-normalisable configuration from the viewpoint of holographic duality. Such non-vanishing mode can generate back-reaction to the background and cause transition to other geometry. However, the minus sign solution always exists as normalisable mode without the lower bound on \(m^{2}_{s}\).

In order to numerically calculate the QNMs, we let \(\phi = e^{-i\omega r_{*}}S\) to linearise the equation of motion (19) with respect to \(\omega \) and change the coordinate by \(u\equiv r_{h}/r\) so that the physical region is \(u\in [0,1]\). The resulting equation of motion is,

$$\begin{aligned}&u^{2}\frac{\partial }{\partial u}\left( f(u) u^{2}\frac{\partial }{\partial u}S(u)\right) +2iw~u^{2}\frac{\partial }{\partial u}S(u) \nonumber \\&\quad -\left[ m^{2}+\left( \lambda ^{2} +\frac{k^{2}}{\alpha _{g}^{2}}\right) u^{2}-f'(u) u^{3}\right] S(u) = 0, \end{aligned}$$
(54)

where \(w\equiv \omega r_{h}, m\equiv m_{s}r_{h}\) are dimensionless parameters.

We shall now calculate the QNMs by using the spectral method (see Ref. [32] and references therein). This is the method of expanding the solution of the quasinormal mode equation with series of orthonormal functions. Each basis function is regular both at the horizon and the infinity so that the solution satisfies the regular boundary conditions, i.e. only normalizable modes are considered. First expand the solution for positive integer N

$$\begin{aligned} S(u) = \sum _{n=0}^{N}b_{n}T_{n}(2u-1), \end{aligned}$$
(55)

where \(T_{n}\) is the Chebyshev polynomials of the first kind. The larger value of N gives the more accurate approximation of S(u). Substituting the expansion into the equation of motion (54), we obtain the linear equation of coefficients \(b_{n}\). By dividing the domain of interest \((2u-1) \in [-1,1]\) into a finite number of grid points and solve the system of linear equations of coefficients \(b_{n}\), the quasinormal frequencies will be determined. In the numerical calculation, we use the Gauss-Lobatto grid points

$$\begin{aligned} u_{k}=\frac{1}{2}\left( 1+\cos \left( \frac{k\pi }{N}\right) \right) , \end{aligned}$$
(56)

where \(k=0,1,\ldots ,N\) and solve the generalized eigenvalue problem to obtain the quasinormal frequencies w for a given N. An excellent example Mathematica code of the spectral method for the calculation of the QNMs is given by Yaffe in Ref. [47].

Table 6 The QNMs \(w=\omega r_{h}~(r_{h}=10.6266)\) of dRGT black string with \(M=1,\alpha _g=1,\epsilon =0,\Lambda =-0.01, \gamma =0, k=0, m=0.2,\lambda =0, 1, 2\)

6.1 Results

For \(\gamma =0\) as shown in Table 6 for each different angular momentum \(\lambda =0, 1, 2\) state at zero k momentum, the lowest lying modes tend to live even longer as the momentum \(\lambda \) increases while the energy also naturally increases. It is apparent from the equation of motion (19) that the effect of \(\lambda \) and k to the QNMs are similar upto the factor of \(\alpha _{g}^{2}\) so we expect the similar dependence on k. Another notable feature is the convergence of the lowest QNMs for high momentum \(\lambda , k \gg 1\) to a simple power-law asymptotic relation

$$\begin{aligned} \text {Im}~w = e^{-\beta }(\text {Re}~w)^{-0.204}, \end{aligned}$$
(57)

where e.g. \(\beta = 0.08, 0.06, 0.0156\) for \(m=0, 0.2, 0.4\) respectively. This result is consistent with what found in Ref. [48] (with \(d=3\) modulo topology of the boundary, the exponent \(0.204\simeq (d-2)/(d+2)=1/5\)) using the WKB method. Figure 7 shows this asymptotic relation for the lowest modes. The high momentum limits in both \(\lambda \) and k thus give the same asymptotic behaviour governed by (57). Naturally for the infinite momentum behaviour, it has a universal power-law for any mass \(m_{s}\). The only effect of \(m_{s}\) shows up at low momentum as different values of \(\beta \). These high-momentum modes have \(\text {Im}(w)\ll \text {Re}(w)\) and thus can be thought of as quasi-particles with long lifetime in the dual gauge theory.

Fig. 7
figure 7

The QNMs \(w=\omega r_{h}~(r_{h}=10.6266)\) of dRGT black string with \(M=1,\alpha _g=1,\epsilon =0,\Lambda =-0.01, \gamma =0, m=0.2\) and \((\lambda , k)=(0,0), (1,0), (2,0), (3,0), (5,0), (7,0), (20,1), (20,100), (20,200)\) (presented by colored dots from the vertical axis out to both sides). In addition to the lowest modes, the next-to-lowest modes appear for sufficiently high momentum \((\lambda , k)\) in a repeated pattern at exactly the same values of \(\text {Re}~w\)

Table 7 The normal-mode frequencies \(w=\omega r_{h}~(r_{h}=3.122)\) of dRGT black string with \(M=1,\alpha _g=1,\epsilon =0,\Lambda =-0.01, \gamma =0.4,m=0\). The momentum parameter \(\ell \) is defined as \(\ell ^{2}\equiv \lambda ^{2}+k^{2}/\alpha _{g}^{2}\)

For positive \(\gamma \), the large-momentum behaviour is different in a remarkable way. Around \(\ell =10\) (defined in Sect. 5.2), some QNMs start to become normal modes with \(\text {Im}(w)\sim 0\), as \(\ell \) grows more QNMs turn to normal modes as shown in Table 7. An interpretation is that the scalar perturbations with high momentum live far away from the horizon in a circular orbit, shielded by the effective potential wall and do not feel the presence of the black string horizon.

Fig. 8
figure 8

The QNMs \(w=\omega r_{h}~(r_{h}\in (0.25,5.8))\) of dRGT black string with \(M=1,\alpha _g=1,\epsilon =0,\Lambda =-0.01, \lambda =k=0, m=0\) and \(\gamma = 0.1,0.2,0.4,0.7,1,2,3,5,10,20,60\) (presented by colored dots from the vertical axis out to both sides) for \(N=300\)

Table 8 The diffusive QNMs for negative \(\gamma , M=1,\Lambda =-\,0.01, \epsilon = 0, \alpha _{g}=1,m=0,\lambda =k=0\) at \(N=300\). T is the corresponding Hawking temperature of the black string

Interestingly, similar to the black hole in dRGT massive gravity studied in Ref. [19] when the horizon size is comparable to the characteristic “AdS” radius \(R_{\Lambda }\equiv \sqrt{1/|\Lambda |}\), the converging QNMs become scarce as shown in Table 6 and Fig. 7. On the other hand for black string with small \(r_{h}< R_{\Lambda }\), the tower of QNMs appears with excellent convergence. Figure 8 shows towers of the QNMs of black string for \(R_{\Lambda }=10\). The real and imaginary parts appear to be almost on the straight line (closer examination shows that it is not exactly a straight line) for each positive value of the massive gravity parameter \(\gamma \). The larger \(\gamma \) the more energetic the QNMs become while the lifetime given by the inverse of \(\text {Im}(w)\) gets larger at first then smaller for very large \(\gamma \). As explained in Ref. [49], the \(\gamma r\) term in the spacetime metric represents the constant radial force \(-\gamma \hat{r}\) originated from the self-interaction of the massive graviton. For positive (negative) \(\gamma \), the massive gravity force is attractive (repulsive). Naturally for positive \(\gamma \), the attractive gravity is stronger so we expect the QNMs to be more energetic, i.e., having larger \(\text {Re}(w)\) values.

In contrast to the \(\gamma \ge 0\) case, the lowest QNMs for negative \(\gamma \) contain two converging diffusive modes (i.e. modes with zero \(\text {Re}(w)\) and negative \(\text {Im}(w)\)) as shown in Table 8. Negative \(\gamma \) implies antigravity force from massive-graviton self-interaction, the larger the magnitude, the stronger the repulsive force. Stronger antigravity leads to larger imaginary parts of the diffusive QNMs and thus shorter lifetime. In constrast to the QNMs of black string in conventional gravity (see e.g., Ref. [25]), the magnitude of the lowest diffusive QNMs of black string with negative \(\gamma \) in massive gravity increases non-linearly with the horizon (\(r_{h}\)) and Hawking temperature (T) even though \(r_{h}\) still increases linearly with T, i.e., \(r_{h}\sim T\sim \gamma \sim m_{g}^{2}\). The approximate power-law relation is given by \(\text {Im}(w)\sim -|\gamma |^{1.89}\) for the parameters in Table 8.

Holographically, quasinormal frequencies are the poles of the retarded thermal Green function in the dual gauge theory living on the boundary [50]. Larger negative imaginary values of QNMs imply that the relaxation time to thermal equilibrium of the perturbed gauge matter in the dual picture is shorter while the temperature is higher due to the antigravity (since \(\gamma \) is negative) effect of massive graviton self-interactions in the bulk. Scalar perturbation in the bulk is dual to the perturbation of certain scalar operator \(\mathcal {O}\) of the gauge theory on the boundary. The generation of such operator could break scale invariance of the dual field theory. The perturbed scalar operator will settle to the thermal equilibrium value in the timescale of the relaxation time given by \(1/\text {Im}(\omega )\) of the QNMs.

For holographic massive gravity models (see e.g., Refs. [23, 24, 51,52,53]), the choice of fiducial metric such as (8) breaks diffeomorphism invariance of the dual gauge theory in the corresponding directions (\(z,\varphi \) in our case). The diffeomorphism-breaking term is effectively the graviton mass term in the gravity theory side. It generates the momentum dissipation in the broken directions in the dual field theory. In coherent regime with \(\gamma =0\), the sound poles of QNMs are found as in Table 6 and Fig. 7. Turning on the massive graviton paramater \(\gamma >0\) does not generate diffusive modes of QNMs, the sound poles simply move away from the imaginary axis as shown in Fig. 8. On the other hand, for \(\gamma <0\), the momentum dissipation appears as there are two purely-imaginary modes found, see Table 8. In contrast to the planar geometry case [24], both diffusive poles move down the imaginary axis as \(\gamma \sim m_{g}^{2}\) becomes more negative.

7 Conclusions

In this paper, we have studied the quasinormal modes of massive scalar perturbations on neutral black string background in dRGT massive gravity. The unique characteristics of dRGT black string are the cosmological constant and the linear term \(\gamma \) in the metric (13) that are generated naturally from the graviton mass via self-interactions. The \(\Lambda \)-like term allows one to consider the black string metric separately in either positive and negative \(\Lambda \) cases. This point makes the black string in dRGT massive gravity different greatly from the cylindrical black object in general relativity [15]. This is because in standard general relativity there is no black string in asymptotically de-Sitter spacetime. We have also examined the horizon structure in each cases. It is found that for positive \(\Lambda \), the black string possesses two horizons which we define as event horizon and cosmological horizon. In the negative \(\Lambda \) case, the black string only has one event horizon. In a generic set-up where the \(\epsilon \) term in the metric (13) is nonzero, the number of horizon structure in negative \(\Lambda \) case will be different from ours [17].

For the QNMs of a black string with positive \(\Lambda \), we start our investigation by considering the spacetime metric (12) in the near extremal limit [34]. In this limit, the cosmological horizon is taken to be very close to the event horizon. It turns out that we can derive the quasinormal frequency analytically by using the Pöschl–Teller technique [37]. The frequencies are labeled by the mode number n. As n increases, the scalar perturbations decay faster. We find the purely imaginary modes when the scalar mass is turned off. In Sect. 5, we perform a fully numerical computation for the QNMs of black string. The numerical technique called improved asymptotic iteration method is implemented [28]. Firstly, the results from numerical AIM and analytical formula are compared in the near extremal limit. We find a perfect agreement between the two methods. The results show that the modes with higher oscillation (higher Re\((\omega )\)) decay faster than the modes with lower oscillation (lower Re\((\omega )\)). In addition, the purely imaginary modes are discovered when the scalar field is massless. However as the scalar mass increases (also the wave number k), the quasinormal frequencies appear to be non-diffusive (nonzero real part). We remark that the third-order WKB approximation gives mostly the same results as those obtained from the AIM. We numerically show that the imaginary parts of \(\omega \) approach asymptotically to a constant at high momentum \(\lambda ,k\gg 1\). Finally, we find no evidence of any instabilities for all dRGT black string with positive \(\Lambda \) since all the perturbation modes investigated in this work decay exponentially.

For negative \(\Lambda \) scenario, the black string that is large comparing to the AdS radius \(R_{\Lambda }\) has few converging modes at low momentum for \(\gamma =0\) (we can choose M to be small to get small black string though). As the momentum of the scalar field gets higher, either with k or \(\lambda \), more converging modes appear at exactly the same \(\text {Re}(w)\) with larger negative \(\text {Im}(w)\). For large momentum limit, we numerically establish asymptotic relation between real and imaginary parts of the QNMs consistent with the WKB method found in Ref. [48]. When the massive graviton effect \(\gamma \sim m_{g}^{2}\) is turned on to positive value, the QNMs becomes sequences of normal modes for high momenta. The higher the momentum, the more number of normal modes appear as shown in Table 7. At zero momentum for a fixed mass, positive \(\gamma \) makes black-string horizon smaller and there is a number of QNMs found. As the effect of massive graviton gets stronger, the energy of the QNMs grows larger, the lifetime becomes longer (i.e., \(\text {Im}(w)\) becomes less negative) at first then starts to get shorter for very large \(\gamma \).

When \(\Lambda<0, \gamma <0\), the massive graviton self-interactions generate antigravity at long distances. With such antigravity, the QNMs remarkably become diffusive with \(\text {Re}(w)=0\). The relaxation time is shorter with stronger antigravity. Holographically, the choice of fiducial metric (8) breaks diffeomorphism invariance in the \(z,\varphi \) directions of the perturbations resulting in momentum dissipation in the dual field theory. We note the importance of sign of \(\gamma \) in the movement of sound poles to collide and form diffusive poles, as \(\gamma \) changes from positive to zero and to negative values. Long-distance gravity changes from attractive to repulsive as the dual hydrodynamics changes from coherent to momentum-dissipative regime.

As a generalization of this work, one can investigate the effect of linear charged scalar perturbation on charged dRGT black string spacetime. With the existence of black string and scalar charges, this could possibly introduce the effect of superradiance [54] on this type of spacetime background. Holographically, unstable charged scalar perturbations in such charged-string (asymptotically AdS) background signals the superconducting phase transition of the dual field system on the cylindrical boundary that resembles e.g., the carbon nanotube.