Quasinormal modes of black holes in f(Q) gravity

In this work, we have studied the quasinormal modes of a black hole in a model of the type f(Q)=∑nanQ-Q0n\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f(Q)=\underset{n}{\sum }a_{n}\left( Q-Q_{0}\right) ^{n} $$\end{document} in f(Q) gravity by using a recently introduced method known as Bernstein spectral method and confirmed the validity of the method with the help of well known Padé averaged higher order WKB approximation method. Here we have considered scalar perturbation and electromagnetic perturbation in the black hole spacetime and obtained the corresponding quasinormal modes. We see that for a non-vanishing nonmetricity scalar Q0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q_0$$\end{document}, quasinormal frequencies in scalar perturbation are greater than those in electromagnetic perturbation scenarios. On the other hand, the damping rate of gravitational waves is higher for electromagnetic perturbation. To confirm the quasinormal mode behaviour, we have also investigated the time domain profiles for both types of perturbations.


I. INTRODUCTION
In recent years, there has been an increase in interest in modified gravity theories (MGT) in trying to answer various unresolved cosmological puzzles, such as the accelerated expansion of the universe and the formation of dark matter.Some of these theories add higher powers of the scalar curvature R, the Riemann and Ricci tensors, or their derivatives to General Relativity (GR).Some of these efforts include the f (R) and f (G) theories, where R and G are the Ricci scalar and the Gauss-Bonnet topological invariant, respectively [1,2].Certain viable f (R) gravity models have been proposed [3], which support the unification of earlytime inflation and late-time acceleration.Viable f (R) gravity models can also be used to solve dark matter issues [4][5][6].In this paper, we examine the quasinormal modes of black holes in the recently suggested f (Q) theory of gravity, where Q is the non-metricity scalar [7,8].The non-metricity Q of the metric mathematically describes the variation in the length of a vector in a parallel transport process, and it is the key geometric variable describing the properties of the gravitational effect.
Generally, the gravitational effects in the space-time manifold can be described using three geometrical objects: curvature R, torsion T , and non-metricity Q.In GR, gravitational effects are assigned to space-time curvature.Another two scenarios, torsion, and non-metricity offer the equivalent representation of GR, and the associated gravity is known as the teleparallel and symmetric teleparallel equivalent of GR.The f (R) theory is a curvature-based extension of GR with zero torsion and non-metricity.Likewise, the f (T ) theory with zero non-metricity and curvature is an extension of torsion-based gravity (the teleparallel equivalent of GR) [9][10][11][12][13][14][15].Finally, the f (Q) gravity theory generalizes GR's symmetric teleparallel (ST) equivalent with zero torsion and curvature.The expanding history of the Universe in f (Q) gravity attracts instant attention as one of the key drivers for this extension [16][17][18][19].Many observational data sets have been used to test the f (Q) gravity in recent years, as seen in the study of Lazkoz et al. [20].The authors used data from the expansion rate, Type Ia Supernovae (SNe Ia), Quasars, Gamma-Ray Bursts, Baryon Acoustic Oscillations (BAO), and the Cosmic Microwave Background (CMB) to constrain f (Q) gravity.Mandal et al. [21] also shown the validity of f (Q) cosmological models in terms of energy conditions.The authors of such a study established the so-called embedding approach, which allows non-trivial contributions from the non-metricity function to be included in the energy conditions.Furthermore, we can witness an increasing interest in f (Q) gravity in the study of astrophysical objects.Wang et al. [22] also studied the static and spherically symmetric solutions for f (Q) gravity in an anisotropic fluid.In [23], black holes in f (Q) gravity have been studied.The authors of [24] investigated the usage of spherically symmetric configurations in f (Q) gravity.Hassan et al. [25] used linear equation of state (EoS) and anisotropic relations to examine wormhole geometries in f (Q) gravity.They discovered exact solutions to the linear model and confirmed the presence of a modest quantity of exotic matter necessary for a traversable wormhole through VIQ.Mustafa et al. have also derived wormhole solutions from the Karmarkar condition and demonstrated the possibility of creating traversable wormholes while maintaining the energy conditions in [26].Furthermore, Banerjee et al. [27] examined the behaviour of energy conditions under constant redshift function by assuming some particular shape functions and confirmed that wormhole solutions could not exist for f (Q) = Q + αQ 2 model.Ref. [28] recently examined a class of static spherically symmetric solutions in f (Q) gravity.Traversable wormholes with charge and non-commutative geometry are studied in Ref. [29].
Black holes are the cleanest objects in the Universe and are highly associated with the generation of gravitational waves.Quasinormal modes are a fascinating and important aspect of black hole physics.They are oscillations of a black hole that are damped over time and characterized by complex frequencies.The term "quasinormal" refers to the fact that these modes are not exactly normal modes, which would oscillate indefinitely.Instead, they die out due to the presence of dissipative effects, such as the emission of gravitational waves.Quasinormal modes are basically complex values that correspond to the emission of gravitational waves from compact and massive objects in the Universe [30][31][32].The real component of the quasinormal modes signifies the emission frequency, while the imaginary component pertains to its decay.Quasinormal modes are significant because they encode information about the black hole's properties, such as its mass, angular momentum, and the properties of the surrounding spacetime.Additionally, the study of quasinormal modes provides insight into the nature of black holes and the strong gravity regime, which is difficult to probe using other means.These modes are important for understanding the structure and evolution of black holes, as well as their role in astrophysical phenomena such as gravitational wave signals.
The properties of gravitational waves and quasinormal modes of black holes have been extensively explored in various modified gravity theories in recent years .In a recent study, quasinormal modes, black hole shadow, greybody bounds etc., are studied extensively in dyonic modified Maxwell black holes [64].In another study, the quasinormal modes have been investigated for the Kerr-like black bounce spacetime under the scalar field perturbations [65].Another research analyses the quasinormal modes of hairy black holes produced by gravitational decoupling for massless scalar fields, electromagnetic fields, and gravitational perturbations [66].In this work, the equation for the effective potential of these three perturbations is determined within the spacetime of the hairy black holes.The time evolution for the three perturbations is also studied, and the quasinormal mode frequencies are calculated using the Prony method based on the time-domain profiles.In a recent publication, the authors demonstrate how quasinormal modes are produced by the perturbations of massive scalar fields in a curved background through the use of artificial neural networks [67].They design a specific algorithm for the feed-forward neural network method to calculate the quasinormal modes that meet specific boundary conditions.To verify the accuracy of the method, they examine two black hole spacetimes whose quasinormal modes are well-established: the 4D pure de Sitter (dS) and fivedimensional Schwarzschild anti-de Sitter (AdS) black holes.Apart from black holes, quasinormal modes are also studied for the wormholes in different frameworks [68,69].In this study, we shall investigate the quasinormal modes of black holes in f (Q) gravity.One may note that the studies of black hole solutions in this theory of gravity are not very old.Until now, very few papers have dealt with the black holes in the framework of f (Q) gravity [22,70].Quasinormal modes, as well as thermodynamics of black holes in f (Q) gravity, are still unexplored.Therefore, being motivated by the previous studies, in this work, we shall investigate the quasinormal modes and time domain profiles of scalar and electromagnetic perturbations of a static black hole in the framework of f (Q) gravity theory.There are two main prime objectives of this study.The first one is investigating the behaviour of quasinormal modes in a recently obtained black hole solution in f (Q) gravity.This study will enable us to see whether the presence of nonmetricity can be probed via quasinormal modes.Moreover, the impacts of the black hole parameters on the quasinormal spectrum will also be investigated thoroughly.The second one is verifying a newly introduced method to calculate the quasinormal modes, which is known as the Bernstein spectral method.To verify the method in this study, we shall use a well-known method, the 6th-order Padé averaged WKB approximation method.
The paper is organized as follows: In Sec.II, we briefly review the f (Q) gravity.Sec.III is devoted to studying the field equations in f (Q) gravity and possible vacuum black hole solutions.In Sec.IV, we briefly investigate the properties of the black hole solution with the help of associated scalars.Sec.V deals with the scalar and electromagnetic perturbation and associated quasinormal modes.In Sec.VI, the time evolution profiles of the perturbations are investigated.Finally, we conclude the paper with a brief discussion of the results and future prospects in Sec.VII.
Throughout the whole paper, we have used G = c = 1.
Weyl geometry is a significant extension of Riemannian geometry, which serves as the mathematical foundation for GR.According to Weyl geometry, during a parallel transport over a closed path, an arbitrary vector will not only change the direction but also the length.As a result, in Weyl's theory, the covariant derivative of the metric tensor is non-zero, and this feature can be represented mathematically in terms of a new geometric quantity named non-metricity Q.Thus, the non-metricity tensor Q γµν can be defined as the covariant derivative of the metric tensor g µν with regard to the general affine connection Γ σ µγ , and it can be expressed as [7,8,71], In this situation, the general affine connection is represented by a Weyl connection and is divided into two independent components as shown below, where the first term is the usual Levi-Civita connection of the metric g µν , as defined by the standard formulation, The second term, which represents the disformation tensor due to the non-metricity of space-time is written as, Furthermore, as a function of the disformation tensor, the contraction of the non-metricity tensor gives the non-metricity scalar, The Weyl geometry can be extended by accounting for space-time torsion, yielding the Weyl-Cartan spaces with torsion.The general affine connection in the Weyl-Cartan geometry can be divided into three independent components as, where the third term is contortion, described in terms of the torsion tensor Further, the relation between curvatures tensor R ρ σµν and Rρ σµν corresponding to the connection Γ γ µν and Γ γ µν is, and the scalar curvature relation, where ∇ is the covariant derivative operator associated with the Levi-Civita connection Γ γ µν .In Symmetric Teleparallel Gravity (STG) case, the curvature-free, and torsion-free requirements restrict the affine general connection.The curvature-free condition necessitates that the Riemann tensor R ρ σµν Γ be zero.So because the Riemann tensor disappears, the parallel transport represented by the covariant derivative ∇ and its corresponding affine connection Γ γ µν is path independent.In addition to the requirement of zero curvature, this theory imposes a torsionless restriction on the connection, i.e.T γ µν = 0, such that gravitation is completely ascribed to non-metricity in STG.The general affine connection is symmetric in its lower indices because the torsion tensor disappears.
As explained previously, in order to obtain the STG, two constraints must be added to the generic affine connection: R ρ σµν Γ = 0 and T γ µν = 0.These constraints permit the choice of a coordinate system {y µ }, in which the affine connection Γ γ µν (y µ ) disappears i.e.L γ µν = −Γ γ µν , resulting in the so-called coincident gauge [7].Thus, in any other coordinate system {x µ }, the affine connection takes the form, Because there exists a coordinate system {y µ } in which the affine connection disappears, we may always imagine that we are working in this special coordinate system with metric as the only fundamental variable.In this case, the covariant derivative ∇ γ reduces to the partial derivative ∂ γ .As a result, in the coincident gauge coordinate, we obtain while in another coordinate system, Finally, the non-metricity scalar can be used to express the action for the f (Q) gravity as [7], Here, L m denotes the Lagrangian density of matter and τ denotes τ = det(g µν ) = √ −g.Similarly to f (R) gravity, f (Q) will be responsible for the deviation from GR, where, for example, if the function f (Q) is considered to be −Q, we recover the so-called Symmetric Teleparallel Equivalent to GR (STEGR).Since the symmetry of the metric tensor g µν , we can just derive two independent traces from the non-metricity tensor Q γµν , Also, it will be convenient to introduce the non-metricity conjugate defined as, The scalar of non-metricity is calculated as follows: Further, the stress-energy momentum tensor for cosmic matter content is determined by The gravitational field equations derived by varying action (14) with regard to the metric g µν are shown below, For the purpose of simplicity, we designate f Q = df dQ .Again, in the absence of hypermomentum [71], we get the connection field equations by varying the gravitational action (14) with regard to the connection,

III. FIELD EQUATIONS AND A SPECIAL CASE OF VACUUM BLACK HOLE SOLUTION
In this study, we shall consider the following static and spherically symmetric ansatz to obtain a static black hole solution.
where the term dΩ 2 = dθ 2 + sin 2 θdϕ 2 .For this ansatz, the nonmetricity scalar Q can be written as [22] Q Following Ref. [22], we start with a constant nonmetricity scalar Q = Q 0 , for which (22) can be written as Now, the components of the field equation can be written as Now, we shall consider the vacuum case with ρ = p r = p t = 0, in the above field equations which reduces them to From the Eq. ( 28), it is possible to have These conditions suggest the following form of f (Q) as given by, where a n represents arbitrary model parameters.Hence, we have seen that for f (Q) gravity to have nontrivial spacetime solutions, the model should satisfy the conditions (30).Now, we fix e a(r) = 1 − 2M r , which is the normal Schwarzschild black hole metric function with black hole mass M .With this assumption, Eq. ( 23) results where one must have Q 0 < 0, and the black hole metric is found to be [22]

IV. PROPERTIES OF THE BLACK HOLE SOLUTION
In this section, we perform a brief study on the black hole's properties by studying the associated scalars.The simplest scalar associated with the black hole spacetime is the Ricci scalar which is given by the following expression: This shows that the Ricci scalar of the black hole solution depends on the constant nonmetricity scalar Q 0 , and for a fixed mass and at a fixed distance, it varies linearly with the nonmetricity scalar.Another scalar associated with the black hole spacetime is the Ricci squared scalar.The expression for the Ricci squared scalar can be given by, Unlike the Ricci scalar, Ricci squared varies nonlinearly with the nonmetricity scalar.Finally, we calculate the Kretchmann scalar, which is given by One may note that this black hole has a physical singularity at r = 0, and the nonmetricity scalar significantly modifies the scalars associated with the black hole spacetime.For a better visualisation of how the nonmetricity scalar impacts these black hole space-time scalars, we have plotted them with respect to r for different values of Q 0 .We plotted the Ricci scalar on the first panel of Fig. 1.One can see that for smaller values of the nonmetricity scalar parameter Q 0 , the Ricci scalar becomes negative outside the event horizon of the black hole.The scalar vanishes at a large distance from the black hole spacetime, showing Minkowski spacetime.On the second panel, we have considered Ricci squared scalar R µν R µν with respect to r.One can see that the Ricci squared is positive, and when Q 0 decreases, it increases slowly.Finally, we plot the Kretschmann scalar on the third panel of Fig. 1.This scalar is also positive, and with an increase in the value of the nonmetricity scalar Q 0 , it decreases slowly.The analysis of the scalars shows that the nonmetricity scalar modifies this black hole spacetime from the standard Schwarzschild black hole spacetime.The black hole solution is unique in nature, and it has a physical singularity at r = 0, which can't be avoided.

V. PERTURBATIONS AND QUASINORMAL MODES
We have obtained the black hole solutions in f (Q) gravity.Now, in this section, we shall deal with two types of perturbations in the black hole spacetime viz., massless scalar perturbation and electromagnetic perturbation.Here we shall assume that the test field i.e. scalar field or electromagnetic field has negligible impact on the black hole spacetime.To obtain the quasinormal modes, we derive Schrödinger-like wave equations for each case by considering the corresponding conservation relations on the concerned spacetime.The equations will be of the Klein-Gordon type for scalar fields and the Maxwell equations for electromagnetic fields.To calculate the quasinormal modes, we use two different methods viz., Bernstein spectral method and Padé averaged 6th order WKB approximation method.
Taking into account the axial perturbation only, one may write the perturbed metric in the following way [41] here the parameters p 1 , p 2 and p 3 define the perturbation introduced to the black hole spacetime.The metric functions g tt and g rr are the zeroth order terms, and they depend on r only.

A. Scalar Perturbation
We start by taking into account a massless scalar field in the vicinity of the previously established black hole.Since we assumed that the scalar field's effect on spacetime is negligible, the perturbed metric Eq. ( 37) may be reduced to the following form: Now, one can write the Klein-Gordon equation in curved spacetime for this scenario as follows: This equation describes the quasinormal modes associated with the scalar perturbation.It is possible to decompose the scalar field as where we have used spherical harmonics and l and m are the associated indices.The function ψ l (t, r) is the radial time-dependent wave function.One can use this equation and Eq. ( 39) to have here r * is expressed as which is known as the tortoise coordinate.V s (r) stands for the effective potential having the following explicit form: Here, l is referred to as the multipole moment of the black hole's quasinormal modes.

B. Electromagnetic Perturbation
Now we move to the electromagnetic perturbation, where one needs to utilise the standard tetrad formalism [41,72] in which a basis e µ a is defined related with the black hole metric g µν .This basis satisfies, e One can express tensor fields in terms of this basis as shown below: One can write the conservation equation as which can be further rewritten in terms of the spherical polar coordinates as In the above expressions, a vertical rule and a comma denote the intrinsic and directional derivative with respect to the tetrad indices, respectively.Using Eq.s ( 45) and ( 46), and the time derivative of Eq. ( 48) one gets, where F = F (t)(ϕ) sin θ.Using the Fourier decomposition (∂ t → − iω) and field decomposition F(r, θ) = F(r)Y ,θ / sin θ, where Y (θ) is the Gegenbauer function one can write Eq. ( 49) in the following form: Using the definition ψ e ≡ r |g tt | F and Eq. ( 42), it is possible to express Eq. ( 50) in the Schrödinger like form as where the potential has the following explicit form: One may note that due to the behaviour of electromagnetic perturbation as shown above, the potential has a simple form and a comparison with Eq. (43) shows that an additional term is absent in the case of the electromagnetic perturbation.This comparison allows us to combine both of them in a compact form or representation as shown below: where the new parameter s stands for spin.In the case of scalar perturbation, we have spin s = 0 which reduces Eq. ( 53) to Eq. ( 43).In the case of electromagnetic perturbation, we have s = 1, which results in the vanishing of the second term of the right-hand side of Eq. ( 53) and reduces it to Eq. ( 52).

C. Behaviour of Potentials
Here we shall briefly study the behaviour of the perturbation potential for the black hole defined above.Since the behaviour of the potential is highly associated with the quasinormal modes, one can have a preliminary idea related to the quasinormal modes from the potential behaviour.
One may note that the potential associated with an electromagnetic perturbation in the normal coordinate system does not depend on the model parameter Q 0 , and it is completely identical to Schwarzschild's case.However, in the case of quasinormal modes, we shall see Q 0 dependency due to the tortoise coordinates, which exhibit Q 0 dependency.Here, we shall only study the potential associated with scalar perturbation.
On the first panel of Fig. 2, we have shown the dependency of the scalar potential with respect to the multipole moment l.On the second panel, we show the impact of black hole mass on the scalar potential.One can see that the potential shows a different behaviour here which is basically due to smaller values of the parameter Q 0 and multipole moment l.We observe that the peak of the potential increases initially with an increase in the value of M .However, at a large distance r, a smaller value of M corresponds to a large value of potential.To have a clear idea of how Q 0 impacts the scalar potential, we have plotted the potential for different values of Q 0 in Fig. 3. On the first panel, we choose smaller values of Q 0 and l = 3.On the second panel, we choose large values of Q 0 with l = 1.One can see that for l = 3, on the first panel, we observe a normal behaviour of potential.Here, with an increase in the value of Q 0 , potential decreases.While on the second panel, we observe peaks on the potential curve for small values of r.With an increase in r, the peaks vanish and the potential decreases to 0 at r = 2 and beyond this point, potential increases again.In this case also, for smaller values of Q 0 , the potential has a larger value.
This analysis shows that the model parameter Q 0 has a significant impact on the potential behaviour.This suggests that the parameter Q 0 may have noticeable impacts on the quasinormal mode spectrum.

D. Bernstein Spectral method for Quasinormal modes
The Bernstein spectral method is a numerical technique used to calculate quasinormal modes of black holes and other objects in GR and MGTs [73].The Bernstein method is based on a spectral decomposition of the wave equation into a finite set of basis functions, which allows for highly accurate and efficient computation of quasinormal modes [73].
Such a method works by approximating the wave equation as a linear combination of basis functions, such as Bernstein polynomials or Chebyshev polynomials.The coefficients of this linear combination are then calculated using a collocation method, where the wave equation is evaluated at a set of discrete points in the domain of interest.By doing this, the wave equation is transformed into a matrix eigenvalue problem, which can be solved using standard numerical techniques such as the Arnoldi method or the Lanczos algorithm.The eigenvalues obtained in this way correspond to the complex frequencies of the quasinormal modes [73].
One of the key advantages of the Bernstein method is that it allows for a very high degree of accuracy, since the basis functions are globally defined and can be tailored to the specific problem at hand.Additionally, the method is computationally efficient, as it only requires the solution of a matrix eigenvalue problem, which can be done very quickly using specialized algorithms.The Bernstein method has been applied to a wide range of problems in both astrophysics and mathematical physics, including the calculation of quasinormal modes for black holes of different masses and spins, the computation of wave scattering by black holes and other objects, and the study of gravitational wave generation by binary black hole systems.
In summary, the Bernstein spectral method is a powerful and versatile numerical tool for the calculation of quasinormal modes in GR as well as MGTs and other areas of physics.Its combination of accuracy and efficiency make it a valuable tool for exploring the fundamental properties of black holes and other gravitational systems, and for the study of gravitational waves and their sources.
In this subsection, we shall implement the Bernstein spectral method [73] to obtain the quasinormal modes for the black hole system considered in this study.The method has been thoroughly discussed in Ref. [73].In another recent study, this method has been revisited for asymptotically de Sitter, anti-de Sitter and flat black hole spacetimes.The corresponding quasinormal modes have been compared with those obtained using time domain analysis [74].We define a compact coordinate u as given by, We also define a wavefunction Ψ(u) which is regular in the range 0 ≤ u ≤ 1, as where ω is the quasinormal frequency, and the other unknown parameters can be determined from the characteristic equations [74].Ω c and Ω + satisfy the conditions The quasinormal boundary conditions representing incoming and outgoing waves are Here r ⋆ implies tortoise coordinate.
In our case, we express the wavelike equation in terms of u = 1/r so that u = 0 and u = 1 correspond to the event horizon and cosmological horizon respectively.After that, we factorise the equation to include e iΩc(ω)/u u α(ω) (1 − u) −iΩ+(ω) so that we obtain Eq. ( 55).This makes the new function regular at u = 0 and u = 1.This is of course due to the inclusion of u and (1 − u) in the function.Hence, both u and (1 − u) play a very important role.More specifically, the new function at the boundaries becomes: Ψ(u = 0) = 0 and Ψ(u = 1) = 0, which are both regular expressions providing mathematical feasibility to do the calculations numerically.The unknown terms such as Ω c (ω), α(ω) and Ω + (ω) are obtained from the characteristic equations utilising series expansion near u = 0 and u = 1.
The function ψ(u) is represented as a sum in the above expression are known as the Bernstein polynomials.We use ( 55) in ( 41) and then a Chebyshev collocation grid [74], where p = 0, ..., N , to get a number of linear equations.Now the differential equation is simplified to the numerically solvable eigenvalue problem of a matrix of order 2 with respect to ω.The solution gives the quasinormal frequencies, and by calculating the relevant coefficients C k and explicitly determining the polynomial (58), which roughly approximates the solution to the wave equation, one can get the complete solution of the problem.We compare the eigenfrequencies and related approximating polynomials for various values of N in order to rule out the erroneous eigenvalues that emerge as a result of the polynomial basis' finiteness.For this purpose, we follow Ref. [74] and calculate where α is the angle between two nearby eigenfunctions ψ (1) and ψ (2) .We set a minimum cut-off α to determine the quasinormal modes from the eigenvalues.For better accuracy, one needs to consider a large value of N , which increases the computation time.
E. WKB method with Padé Approximation for Quasinormal modes Apart from the above method, we shall use another well-established method known as the WKB method to estimate the quasinormal modes of the black hole considered in this study.We shall compare the results obtained from both methods to verify our findings.
The first-order WKB method or technique was introduced for the first time by Schutz and Will in Ref. [75].Although this method can give approximate results of quasinormal modes, the error associated with this method is comparatively higher.For this reason, higher-order WKB methods have been implemented in the study of black hole quasinormal modes.The WKB method was upgraded to higher orders in Ref.s [76][77][78].It was stated in Ref. [78] that the WKB technique may be improved by averaging the Padé approximations.Subsequently, it was discovered that this improved the findings of quasinormal modes with more precision [77].The Padé averaged 6th-order WKB approximation approach will be used in this investigation.
Each table shows the quasinormal modes obtained from Padé averaged 6th order WKB approximation method and Bernstein spectral method.The errors associated with the WKB method are represented by the rms error △ rms and ∆ 6 , which is defined as [77] ∆ where the terms ω 7 and ω 5 represent quasinormal modes obtained from 7th and 5th order Padé averaged WKB method.
TABLE I: Quasinormal modes of the black hole with n = 0, M = 1 and Q 0 = −0.1 for the massless scalar perturbation.In Table I, we have listed the quasinormal modes for the massless scalar perturbation with the model parameter Q 0 = −0.1, the mass of the black hole M = 1 and overtone number n = 0.The term ∆ BW represents the percentage deviation of the quasinormal modes obtained via the Bernstein spectral method from the 6th order padé averaged WKB approximation method.Here in the Bernstein spectral method, a large value of N increases the accuracy [79] but it needs a larger computation time.We have used N = 20, 30, 40 and 60 to compare the results for excluding the spurious eigenvalues.We compare the results by extracting the eigenvalues from each set that differ by an amount smaller than the specified cutoff value (around 0.01).For each pair of eigenvalues meeting this criterion, we calculate the squared sin of the angle between their respective eigenfunctions using Eq.61.When the sin value approaches zero or less than a cutoff value (we took it to be 1 × 10 −10 ), it signifies that the eigenfunctions exhibit a negligible difference and approximate a common solution, possibly differing by a constant factor.By employing this technique, we effectively identify the dominant quasinormal modes.We see that excluding spurious eigenvalues and considering a suitable precision of numerical calculations allow us to choose small values of N say around 40 to obtain the desired accuracy instead of considering a large value such as N = 100 which needs a higher computation time.One can see that for l = 1, ∆ BW has a comparatively higher value.Apart from this, the errors associated with the Padé averaged WKB method are also higher.The percentage deviation, as well as the errors, decrease significantly as the multipole moment l increases.The variation of these errors can be associated with the WKB method, which gives less significant results when l − n is smaller.For both methods, it is clearly seen that the quasinormal frequencies and the damping rate increase with the value of the multipole moment l.
In Table II, we have shown the quasinormal modes for electromagnetic perturbation for different l values with M = 1, n = 0, and Q 0 = −0.1.One can see that in this case also, the deviation term ∆ BW has higher values for smaller multipole moment l.As l increases, the deviation, as well as the error terms, decreases noticeably.
From both Tables, it is clear that quasinormal frequencies have lower values for the electromagnetic perturbation than those obtained in the case of massless scalar perturbation.However, the decay rate or damping rate of GW is higher in the case of electromagnetic perturbation.
To see the effect of the model parameter Q 0 on the quasinormal modes, we have plotted the variations of quasinormal modes with respect to Q 0 .On the left panel of Fig. 4, we have plotted the real frequencies of quasinormal modes with respect to the nonmetricity scalar Q 0 for both massless scalar perturbation and electromagnetic perturbation.One can see that the scalar quasinormal mode frequency decreases significantly with an increase in the nonmetricity scalar Q 0 .On the other hand, electromagnetic quasinormal frequencies increase slowly with an increase in the value of Q 0 .Both types of frequencies approach each other when Q 0 is close to zero.One may note that Q 0 has more impact on the scalar quasinormal frequencies than the electromagnetic quasinormal frequencies.On the second panel of Fig. 4, we have plotted the imaginary quasinormal modes with respect to the nonmetricity scalar Q 0 .We observe that the damping rate or decay rate of gravitational waves increases with a decrease in the value of the scalar Q 0 .The decay rate for the scalar perturbation is lower than that of electromagnetic perturbation for smaller values of Q 0 near −1.However, when Q 0 increases to 0, the damping rate decreases and approaches zero.One may note that the variation of the damping rate with respect to Q 0 is nonlinear, unlike the case for real quasinormal modes, where the variation was linear for both types of perturbations.Our analysis shows that the quasinormal mode spectrum carries the signature of nonmetricity, and in the near future, provided sufficient significant observational data is available, it might be possible to distinguish the scalar and electromagnetic quasinormal modes as well as the presence of nonmetricity in theory from a standard Schwarzschild black hole in GR.

VI. EVOLUTION OF SCALAR AND ELECTROMAGNETIC PERTURBATIONS ON THE BLACK HOLE GEOMETRY
In the previous section, we have numerically calculated the quasinormal modes and studied their behaviour with respect to the nonmetricity scalar Q 0 .In this section, we shall deal with the time domain profiles of the scalar perturbation and electromagnetic perturbation.To obtain the time evolution profiles, we shall implement the time domain integration formalism [80].For this purpose, we define ψ(r * , t) = ψ(i∆r * , j∆t) = ψ i,j and V (r(r * )) = V (r * , t) = V i,j .Now, one can express Eq.(39) as ψ i+1,j − 2ψ i,j + ψ i−1,j ∆r 2 * − ψ i,j+1 − 2ψ i,j + ψ i,j−1 ∆t 2 − V i ψ i,j = 0.
Using the above iteration scheme and choosing a fixed value of ∆t ∆r * , one can easily obtain the profile of ψ with respect to time t.However, one should keep ∆t ∆r * < 1 so that the Von Neumann stability condition is satisfied during the numerical procedure.On the left panel of Fig. 5, we have shown the time domain profile for the scalar perturbation and on the right panel, we have shown the time domain profile for the electromagnetic perturbation with overtone number n = 0, nonmetricity scalar constant Q 0 = −0.1 and multipole number l = 1, 2, 3 and 4. We can see that the time domain profiles are significantly different for different values of multipole number l.For both cases, with an increase in the value of l, we observe an increase in the frequency.However, the decay rate of the time domain profile seems to increase for the scalar perturbation with an increase in the value of l while for electromagnetic perturbation the variation is very small.Another interesting observation is that for l = 1, both the scalar and electromagnetic time domain profiles have significantly different damping rate.The damping or decay rate seems to higher for the electromagnetic perturbation case than that of the scalar perturbation.
In Fig. 6, we have plotted the time domain profiles for scalar perturbation (on the left panel) and electromagnetic perturbation (on the right panel) with different values of the parameter Q 0 with overtone n = 0 and multipole number l = 2.One may note that the decay rate for the time profile of electromagnetic perturbation increases more rapidly with a decrease in the value of the parameter Q 0 .The results obtained from the time domain profiles are consistent with the previous results of quasinormal modes.

FIG. 2 :
FIG. 2: Variation of scalar potential V s (r) w.r.t.r.On the first panel, we have used Q 0 = −0.1 and M = 1.On the second panel, we have used Q 0 = −1.5 and l = 1.

FIG. 3 :
FIG.3: Variation of scalar potential V s (r) w.r.t.r.On the first panel, we have used l = 3 and M = 1.On the second panel, we have used l = 1 and M = 1.

FIG. 4 :
FIG.4: Variation of quasinormal modes with respect to model parameter Q 0 with M = 1, n = 0 and l = 4 for massless scalar perturbation and electromagnetic perturbation.