No constraints for f ( T ) gravity from gravitational waves induced from primordial black hole ﬂuctuations

. Primordial black hole (PBH) ﬂuctuations can induce a stochastic gravitational wave background at second order, and since this procedure is sensitive to the underlying gravitational theory it can be used as a novel tool to test general relativity and extract constraints on possible modiﬁed gravity deviations. We apply this formalism in the framework of f ( T ) gravity, considering three viable mono-parametric models. In particular, we investigate the induced modiﬁcations at the level of the gravitational-wave source, which is encoded in terms of the power spectrum of the PBH gravitational potential, as well as at the level of their propagation, described in terms of the Green function which quantiﬁes the propagator of the tensor perturbations. We ﬁnd that, within the observationally allowed range of the f ( T ) model-parameters, the obtained deviations from general relativity, both at the levels of source and propagation, are practically negligible. Hence, we conclude that realistic and viable f ( T ) theories can safely pass the primordial black hole constraints, which may oﬀer an additional argument in their favor.


Introduction
Modified gravity is one of the two main avenues that one can follow in order to describe the early and late phases of Universe's acceleration [1,2], and compared to the other alternative, namely the introduction of the inflaton/dark-energy concept [3][4][5], has the additional advantage of being closer to the quantum description of gravity [6].Although the simplest way to build novel classes of gravitational theories is to start from the standard, curvature, formulation of gravity and extend it in various ways [7][8][9][10][11], one can equally well follow the alternative, torsional, formulation and extend it suitably.In particular, since the basic torsional theory, namely the teleparallel equivalent of general relativity (TEGR) [12,13] uses the torsion scalar T as the Lagrangian, one can construct torsional modifications extending it such as in f (T ) gravity , in f (T, T G ) gravity [41,42], in f (T, B) gravity [43], in scalar-torsion theories [44,45], etc.
On the other hand, primordial black holes (PBHs), firstly introduced in the early '70s [46][47][48], have rekindled the interest of the scientific community given the fact that they can solve a number of fundamental issues of modern cosmology.In particular, they may indeed constitute a viable candidate for dark matter [49,50], and explain the largescale structure formation process through the Poisson fluctuations they can seed [51,52].At the same time, depending on their mass they can give access to a wide variety of physical phenomena from the early universe up to late times [53].
Hence, given the motivation behind modified gravity theories, one can use the aforementioned GW portal associated to PBHs in order to constrain them.Specifically, there have been some first attempts in [74,75], where the authors study the primordial scalar induced GWs within the context of Hordenski gravity and non-canonical Higgs inflation, as well as in [76] where the scalar induced GWs from PBH Poisson fluctuations were studied within f (R) gravity theories and in particular within Starobinsky inflation.
In this manuscript, we focus on f (T ) modified gravity, and we study the scalar induced GWs from Poisson fluctuations of ultralight PBHs (i.e. with m PBH < 10 9 g), which evaporate before BBN and transiently dominate the energy content of the universe before their evaporation [59,[77][78][79].Thus, the main goal of the work is to examine whether such an analysis will impose constraints on the various specific f (T ) models, similarly to other observational investigations [80][81][82][83][84][85][86][87][88].
The plan of the work is as follows: In Sec. 2 we review the calculation of the PBH gravitational potential of Poisson distributed PBHs within general relativity, and in Sec. 3 we perform the extended analysis, extracting the PBH gravitational potential in the framework of f (T ) gravity.In Sec. 4 we focus on three mono-parametric f (T ) gravity models and we present the formalism for the computation of the relevant scalar induced gravitational wave signal.Then, in Sec. 5 we investigate the modifications of the GW signal within f (T ) gravity by accounting for GW source and GW propagation effects, showing that f (T ) theories of gravity safely pass the constraints imposed by the gravitational-wave signal induced from PBH Poisson fluctuations.Finally, Sec. 6 is devoted to conclusions.

The primordial black hole gravitational potential in general relativity
In this section, we briefly review the calculation of the PBH gravitational potential in the case of general relativity, following [71,76].We first present the background and perturbation equations, and then we extract the power spectrum of the PBH gravitational potential.

Background evolution and scalar perturbations
We consider a flat Friedmann-Lemaître-Robertson-Walker (FLRW) background geometry with metric where a(t) is the scale factor.Additionally, we consider that the Universe is filled with hydrodynamic fluid matter described with energy-momentum tensor T m µν = diag(−ρ, p, p, p), where ρ and p are the total matter (i.e.including radiation, baryonic and dark matter) energy density and pressure.Thus, the Friedmann equations are where G is the Newton constant (throughout this paper we work in units where c = 1), Λ is the cosmological constant, H = ȧ/a is the Hubble parameter and dots denote derivatives with respect to the cosmic time t.Note that we have introduced the total background energy density and pressure, ρtot and ptot , which include the total matter sector alongside with the cosmological constant.
In order to proceed to the investigation of the scalar perturbations it proves convenient to introduce the conformal time η, defined as dt ≡ adη, and hence the conformal Hubble parameter reads as H ≡ a /a = aH, where primes denote derivatives with respect to η. Restricting ourselves to scalar perturbations in the Newtonian gauge, we can write the perturbed metric as where Ψ and Φ are the two Bardeen potentials [89].Additionally, we include perturbations around the background stress-energy tensor of the total matter content of the Universe (matter and radiation) which we express as follows: where δ ≡ δρ/ρ is the relative energy density perturbation, δu i ≡ υ i /a is the velocity perturbation and Π i j is the (dimensionless) anisotropic stress.In this context, one obtains the following perturbed equations for Φ and Ψ [90]: (2.9) In the time period that we focus on, Π i j is negligible and therefore Φ ≈ Ψ, which we consider to be the case from now on.Elaborating on the above equations, and using also the (total) conservation equation ρ = −3H(ρ + p), one obtains with w ≡ p/ρ being the equation-of-state parameter and c 2 s ≡ p /ρ the sound speed square of the total matter content sector.Hence, during the period of PBH domination, Φ is the potential arising from the PBH distribution.

The PBH gravitational potential power spectrum
In order to proceed we assume conventionally that PBHs are formed in the radiationdominated era.Considering that PBHs are randomly distributed, their energy density is inhomogeneous while the total background energy density is homogeneous, and thus their energy density perturbations can be viewed as isocurvature Poisson fluctuations with the associated Poissonian power spectrum for the PBH density contrast being read as [71] where we have assumed monochromatic PBH mass function [91].Additionally, k UV ≡ a/r is the Ultraviolet (UV) cut-off scale, which is related to the mean PBH separation scale, since at smaller scales the PBH fluid description is not valid.Introducing then the density parameter for the PBH fluid as Ω PBH ≡ ρ PBH ρtot we can find that during the radiation-dominated era, Ω PBH ∝ a, and therefore if the initial PBH abundance is sufficiently large then PBHs can dominate.Hence, the isocurvature PBH perturbations in the radiation-dominated era will be converted to adiabatic curvature perturbations in the subsequent PBH dominated era [92,93].
In order to relate Φ and δ PBH , we introduce the uniform-energy density curvature perturbation of each fluid, ζ i [94].At super-horizon scales, ζ r and ζ PBH are separately conserved, hence at the PBH formation time we can neglect the adiabatic contribution associated to the radiation fluid on the scales we are interested in and obtain Thus, during the PBH-matter dominated era, where w = 0 and Φ is constant in time [95] using the fact that on super-horizon scales ζ −R [94], where R is the comoving curvature perturbation, as well as the relation between R and Φ in GR [95], we finally obtain that [71,76] On the other hand, at sub-horizon scales one can determine the evolution of δ PBH by solving the evolution equation for the matter density perturbations with k d ≡ H d being the comoving scale exiting the Hubble radius at PBH domination time. 3 The primordial black hole gravitational potential in f (T ) gravity In this section we perform the calculation of the PBH gravitational potential in the framework of f (T ) gravity.We first review the relevant background and perturbation equations and then we proceed to calculate the associated power spectrum.

Background evolution and scalar perturbations
In the torsional formulation of gravity one uses the tetrad or vierbein fields e A (x µ ) as the dynamical variables instead of the metric tensor.They form an orthonormal basis for the tangent space at each point x µ of the manifold, that is e A • e B = η AB , where Greek and Latin indices run in coordinate and tangent space respectively and η AB = diag(−1, 1, 1, 1) is the Minkowski metric for the (flat) tangent space.One can express them in the coordinate basis e A = e µ A ∂ µ and thus construct the metric as In teleparallel gravity one describes gravity using the torsion of spacetime instead of its curvature.In this spirit, instead of the familiar Christoffel connection, which is the unique connection whose torsion vanishes, one can introduce the Weitzenböck connection , which is a connection whose curvature vanishes [12].The torsion tensor is given by and its contraction provides the torsion scalar as Using T as a Lagrangian gives rise to the teleparallel equivalent of general relativity (TEGR), since variation in terms of the tetrads leads to the same field equations with general relativity [12].
One then can generalize TEGR by extending T to an arbitrary function of T as the Lagrangian, resulting to f (T ) gravity, whose action is [14] with |e| = det(e A µ ) = √ −g , and where we have included the total matter Lagrangian L m for completeness.Variation of the action (3.4) with respect to the tetrad e A µ yields the following field equations: with f T ≡ ∂f /∂T , f T T ≡ ∂ 2 f /∂T 2 and where em T ρ ν denotes the total matter energymomentum tensor.Note that for convenience we have introduced the super-potential tensor S µν ρ Applying f (T ) gravity in a cosmological framework we impose the FLRW metric (2.1), which in turn arises from the tetrad e A µ = diag(1, a, a, a).Inserting this ansatz into (3.5)one obtains the familiar Friedmann equations where we have defined the effective energy density and pressure due to the f (T ) modification as [14] ρ with equation-of-state parameter given by Thus, ρtot ≡ ρ+ρ f(T) and ptot ≡ p+p f(T) .Lastly, note that in FLRW geometry according to (3.3) the torsion scalar becomes simply T = 6H 2 .
Proceeding to the perturbation level and focusing on scalar perturbations, one can write the perturbed tetrad fields as follows: where the symbols e A µ and ēA µ are used for the perturbed and the unperturbed tetrad fields correspondingly.Then, one can impose the following ansatz for the scalar contributions of the perturbed tetrad fields: and In the above expressions, the scalar perturbations Φ and Ψ are introduced, which are functions of space x and time t.With such a choice, one can match the tetrad perturbations with a perturbed metric in the Newtonian gauge [17,97], namely Following then [14], and writing f (T ) = T + F (T ), one can expand the gravitational equations of motion (3.5) to linear order and obtain the (00), (0i), (ij) and (ii) component perturbation equations, which read as follows: ) and where the total matter content of the Universe is expressed as in (2.5).
In the context of this work, all the models of f (T ) gravity that we are going to investigate (see Sec. 5) are the viable ones, characterised by F T 1, F T T 1 and F T T T 1 [82].Furthermore, as we already mentioned in the previous section, the anisotropic stress can be neglected, and therefore by virtue of (3.17) one obtains that Φ ∼ Ψ.Therefore, under these conditions it is straightforward to show that from the form of the aforementioned equations one can use for the evolution of Φ the GR equation (2.10), namely: On the other hand, we stress out that the effect of f (T ) gravity will be taken into account at the level of the evolution of the matter density contrast, as we will show in the next subsection, as well as at the level of the tensor perturbations (see Sec. 4).
Regarding the evolution of the matter density contrast, it is given via the growth equation, which is derived from (3.15) and (3.16) at subhorizon scales, and assuming matter domination it is given by [82,98] In this expression the quantity is the effective Newton's constant and primes denote derivatives with respect to the conformal time.Thus, comparing (2.14) and (3.20), one can see that the effect of f (T ) gravity at the level of the matter perturbations at sub-horizon scales is essentially captured by the modification of the Newton's constant, which is related to the modification of the gravitational field.As we can see, in the limit F (T ) → const.= −2Λ all the above equations recover the ones of ΛCDM cosmology.

The Power Spectrum of the PBH Gravitational Potential in f (T ) gravity
We shall now repeat the procedure of subsection 2.2 but in the context of f (T ) gravity.
Once again, we will use the uniform-energy density curvature perturbation of each fluid, ζ i , in order to relate Φ and δ PBH .For ζ we will use the usual definition [94]: The relevant fluids for our analysis are the radiation and the PBH-matter ones.Their corresponding energy-momentum tensors are both conserved, i.e. the (background) continuity equation ρ i = −3H(ρ i + pi ) holds for each fluid, and thus ζ i is expressed as: By substituting the equation of state w i for each fluid we obtain: Moreover, we introduce the isocurvature perturbation defined as: On superhorizon scales, ζ r and ζ PBH are conserved separately [94], like the isocurvature perturbation S. Thus, in the PBH-dominated era, ζ ζ PBH = ζ r + S/3 S/3.Since S is conserved, it can be calculated at formation time t f .Consequently, neglecting the adiabatic contribution associated to the radiation fluid at the PBH formation time, since it does not play any role at the scales considered here, from Eq. (3.26) we obtain that S = δ PBH (t f ).Hence, as in the case of GR we find that Further, as we show explicitly in Appendix A, in the context of f (T ) gravity the property ζ −R is valid too at super-horizon scales as it does in GR (see e.g.[94]), by requiring that F T 1, with R being the comoving curvature perturbation defined in the usual way During a matter-dominated era, such as the one driven by PBHs, Φ can be neglected since it is proportional to the decaying mode, thus we obtain R = −ζ = (5/3)Φ.Therefore, combining with (3.27), we deduce that Let us now focus on sub-Hubble scales.We can determine the evolution of δ PBH by solving the evolution equation of the matter density perturbations in f (T ) gravity (3.20).At the background level, the Friedmann equation (3.6) can be expressed as where ρf(T) is given by (3.8).Since at the epochs we focus on, namely before Big Bang Nucleosynthesis (BBN), we expect that deviations from ΛCDM are negligible, we can neglect the effective fluid contribution, writing the Friedmann equation as In this equation, s ≡ a/a d and a d denotes the time at the transition from the radiation to the PBH domination era, while we have assumed that Ω r,f = 1 since PBHs are considered to be formed in the radiation era [71].Note that the scale factor is normalised at one at formation time, i.e. a f = 1.At the perturbation level, we can treat the gas of PBHs as a matter fluid, and by using s as the time variable, the growth equation (3.20) Hence, making an interpolation between Eq. (3.29) and Eq.(3.32), as in the case of GR, and using the expression for the PBH matter power spectrum in Eq. (2.11), we straightforwardly extract the following PBH gravitational potential power spectrum: In the above expression, ξ(a) is defined as where δ PBH (a) is the solution of Eq. (3.31).As verified numerically, ξ(a) has a mild dependence on the comoving scale k, and thus for practical reasons we will consider ξ(a) as k independent.
4 Scalar induced gravitational waves in f (T ) gravity Since we have calculated the power spectrum of the gravitational potential of the initially Poisson-distributed PBHs, we can now proceed to the extraction of the stochastic gravitational wave background induced from the PBH Poisson fluctuations.

Tensor Perturbations
The perturbed metric in the Newtonian gauge, assuming as mentioned above zero anisotropic stress, can be recast as where we have multiplied by a factor 1/2 the second-order tensor perturbation as it is standard in the literature.Then, by Fourier transforming the tensor perturbations and taking into account the two polarization modes of the GWs in f (T ) gravity [99], namely the × and the + as in GR case, the equation of motion for the tensor modes h k reads as with s = (+), (×).In this equation the modified dispersion due to the f (T ) effects is quantified by the term [100] while the source function S s k is given by ) where the polarization tensors e s ij (k) are defined as [2] e As we have mentioned above, since we focus on second-order effects, in this work we assume that the background evolution is close to that of the ΛCDM scenario -note that we consider PBH domination eras before BBN time.Considering also the fact that in the time period that we investigate the Universe is matter (i.e.PBH) dominated, we have c 2 tot ≈ w tot w PBH = 0. Hence, for the time evolution of the potential Φ given by Eq. (3.19), we obtain The solution of the above equation is a superposition of a constant and a decaying mode.In the late-time limit, one can neglect the decaying mode, and write the solution for the Fourier transform of Φ as Φ k (η) = T Φ (η)φ k , where φ k is the value of the gravitational potential at some reference time (which here we consider to be the time at which PBHs dominate the energy budget of the Universe, x d ) and T Φ (η) is a transfer function, defined as the ratio of the dominant mode between the times x and x d .Consequently, Eq. (4.4) can be written in a more compact form as where and the contraction e s ij (k)q i q j ≡ e s (k, q) can be expressed in terms of the spherical coordinates (q, θ, ϕ) of the vector q as e s (k, q) = 1 √ 2 q 2 sin 2 θ cos 2ϕ for s = (+) 1 √ 2 q 2 sin 2 θ sin 2ϕ for s = (×) . (4.9) Finally, the solution of Eq. ( 4.3) for the tensor modes h s k can be obtained using the Green's function formalism where one can write for h s k that and where the Green's function G s k (η, η) is the solution of the homogeneous equation with the boundary conditions lim η→η G s k (η, η) = 0 and lim η→η G s, k (η, η) = 1.Having extracted above the tensor perturbations, the next step is to derive the tensor power spectrum, P h (η, k), for the different polarization modes, which is defined as the equal-time correlator of the tensor perturbations through the following relation: where s = (×) or (+).Finally, after a long but straightforward calculation one acquires that P h (η, k) for the (×) and (+) polarization states can be recast as [101-104] 2 I 2 (u, v, x)P Φ (kv)P Φ (ku) .
(4.13)The two auxiliary variables u and v are defined as u ≡ |k − q|/k and v ≡ q/k and the kernel function I(u, v, x) is given by In the above expressions x = kη, and we use the notation since the function F (q, k − q, η) depends only on the modulus of its first two arguments.

The gravitational-wave energy-density spectrum
In this subsection we calculate the energy density associated to the scalar induced GWs, focusing only on subhorizon scales.Consequently, after a lengthy but straightforward calculation the GW energy density can be recast as [105] ρ which is simply the sum of a kinetic term and a gradient term.The overall bar denotes an oscillation averaging on sub-horizon scales, performed to deduce the envelope of the gravitational-wave spectrum.The GW spectral abundance is just the GW energy density per logarithmic comoving scale, i.e.
Considering a matter-dominated era driven by PBHs, where w = 0, the transfer function T Φ is constant in time (see the discussion after Eq. (4.6)) , and we normalise it to one at PBH domination time, namely T Φ (x d ) = 1.This forces the source term S s k to be constant in time and consequently at sub-horizon scales, where k H, from Eq. ( 4.3) we obtain that h s k 4S s k k 2 .Finally, the tensor modes have a mild dependence on time and therefore the kinetic term in the expression for the GW energy density (4.15) can be neglected.Therefore, we straightforwardly obtain that where the brackets stand for an ensemble average.At the end, by combining Eq. (4.17), Eq. (4.16) and Eq.(4.12) and taking into account from Eq. (4.13) that the (×) and (+) polarization modes give an equal contribution, we find that In order to compute the contribution of the induced GWs to the energy budget of the Universe at the present epoch, one should evolve Ω GW (η, k) from a reference conformal time η * up to today.To do so, one has that where we have taken into account that Ω GW ∼ a −4 , and where the index 0 refers to the present time.Then, taking into account that the energy density of radiation can be recast as ρ r = π 2 15 g * ρ T 4 r and that the temperature of the radiation bath, T r , scales as where g * ρ and g * S stand for the energy and entropy relativistic degrees of freedom.

Constraints on f (T ) gravity
In this section we use the portal of the scalar induced GWs from PBH Poisson fluctuations presented above, in order to derive constraints on f (T ) gravity.

Mono-parametric f (T ) Models
Since we will perform specific calculations we have to specify the form of the function f (T ).In particular, we consider the following three f (T ) gravity models depending on one free parameter, denoted here as β.

Theoretical parameters
In this subsection we discuss the theoretical parameters that are involved in the analysis.These parameters include the mass of the PBH m PBH , the initial PBH abundance at formation time Ω PBH,f , and of course the parameter β of the mono-parametric f (T ) model at hand.For the PBH mass range we consider that the PBHs are formed after the end of inflationary era and evaporate before the BBN time.In particular, one can derive an upper bound on the PBH mass by accounting for the current Planck upper bound on the tensor-to-scalar ratio for single-field slow-roll models of inflation, which gives ρ 1/4 inf < 10 16 GeV [107].On the other hand one can extract a conservative lower bound on the reheating energy scale, i.e. ρ 1/4 reh > 4MeV, derived by taking into account the thermalization of neutrino background and the hadron scatterings emitted from PBHs as discussed in [108].Hence, one can show that the relevant PBH mass range is given by [71] 10g < m PBH < 10 9 g. (5.7) We proceed to the range of Ω PBH,f .In order to have a transient PBH domination era, we can set it by requiring that the PBH evaporation time t evap is larger than the PBH domination time t d .Consequently, knowing that during a radiation dominated era we have Ω PBH = ρ PBH /ρ d ∝ a −3 /a −4 ∝ a and demanding that t evap > t d , we obtain that Ω PBH,f > 10 −15 g eff 100 10 9 g m PBH . ( Lastly, concerning the parameter β, according to observational constraints mentioned above from [82,98,106], its value should roughly vary within the following range: depending on the f (T ) model at hand.

Gravitational waves from an early primordial black hole dominated era
We have now all the necessary material in order to investigate the relevant GW signal created in the early PBH dominated era, in the context of f (T ) gravity.Let us mention here that one should discriminate between two decisive effects which introduce deviations from standard GR gravity, namely a) the effect of the source of the induced GWs, which is actually encapsulated within the power spectrum of the PBH gravitational potential P Φ , and b) the effect of the GW propagation, which is encoded within the time evolution of the Green function G k (η, η), which can be viewed as the propagator of the tensor perturbations as it can be seen by Eq. (4.10).

The effect of the gravitational-wave source
After solving numerically the Meszaros equation (3.31) for the PBH energy density perturbations, we proceed to the calculation of the PBH gravitational potential power spectrum P Φ , as dictated by Eq. (3.33), by accounting for the three mono-parametric f (T ) models presented in Sec.5.1.The power spectrum P Φ is actually the source of the induced GWs as it can be seen from Eq. (4.13).
In Figs.1-3 we depict the power spectrum of the PBH gravitational potential at PBH domination time, for the three f (T ) models, for various choices of the initial PBH abundance Ω PBH,f , of the PBH mass m PBH , and of the model-parameter β.As it can be clearly seen from Fig. 1, the amplitude of the PBH gravitational potential power spectrum increases with the initial PBH abundance, whereas from Fig. 2 we can notice that the position of the peak of P Φ depends on the PBH mass.Finally, from Fig. 3, where P Φ is plotted for different values of the modified gravity parameter β, one can clearly infer that the deviation from GR is practically indistinguishable.As it was verified numerically for all the monoparametric f (T ) models considered here, the time evolution of the GW spectrum one should take into account the fact that during a sudden transition from the PBH dominated era to the subsequent radiation era, in the case of a monochromatic PBH mass function as the one considered here, the GW is enhanced due a very rapid increase of the time derivative of the gravitational potential Φ which is present in the source term (see Eq. (4.7)) as noted in [109,110].Then, during the radiation-dominated era, the source term is decaying on subhorizon scales and therefore Ω GW stops growing after the moment when the source term has sufficiently decayed.After this point, the scalar induced GWs are propagating as free waves, with their present energy density spectrum given by Eq. (4.20) and η * being a reference time during the radiation-dominated era when GWs start to propagate as free waves.
Having these in mind let us focus on Eq. (4.11), namely (5.11) One can identify the dominant terms of the above equation by taking the ratios between the GR terms and the new f (T ) terms multiplied by the γ T function.In this procedure, we should take into consideration the fact that the γ T function, as it was checked numerically for all the three mono-parametric models studied here, is a negative decreasing function of time, which implies that its absolute value increases with time.Therefore, in order to find the maximum deviation from GR we compute the ratios between the GR and f (T ) terms at a time during radiation domination when the γ T function acquires its maximum value.Being quite conservative we set this time to be the standard matter-radiation equality time at redshift z eq = 3387.
At this point, we should stress that in the comparison of the different terms in Eq. (5.11) one should compute the derivative terms G k (η, η) and G k (η, η).To achieve this we take into account the fact that the solutions of Eq. (5.11) are expected to be trigonometric functions (sines and cosines), in particular Bessel functions in the case of GR.One then expects that G k (η, η), G k (η, η) and G k (η, η) should differ by a phase difference, and this was indeed verified numerically.Consequently, one expects in gen- 1), when comparing the different terms of Eq. (5.11).
Thus, taking the above discussion into account, let us start the identification of the dominant terms by comparing the first two terms in Eq. (5.11), namely the second derivative term G s, k (η, η) and the friction term 2Hγ T G s, k (η, η), and in particular let us consider their ratio G s, k (η, η)/[2Hγ T G s, k (η, η)].Eventually, we find that for the powerlaw f (T ) model and for m PBH = 10 5 g, Ω PBH,f = 10 −3 and β = 0.1, we acquire 1. (5.12) Similar results are obtained for the square-root exponential and the exponential f (T ) models, and by varying the parameters m PBH , Ω PBH,f and β too.
With the same reasoning we can additionally examine the ratio between the k 2 and 2H 2 γ T terms inside the parenthesis of Eq. (5.11).Choosing again η = η eq in order to find the maximum deviation from GR, we straightforwardly find that where k = k evap is the comoving scale exiting the Hubble radius at PBH evaporation time and as a consequence it is the largest scale considered here.In summary, we can safely argue that the f (T ) modifications at the level of the propagation equation (5.11) can be neglected and consequently one obtains that G f (T ) k (η, η) G GR k (η, η). (5.14) Hence, we conclude that the effect of f (T ) gravity on the propagation of gravitational waves, is very small, and practically indistinguishable from GR for realistic f (T ) modelparameter values.We mention that this behavior is different than the case of f (R) gravity, in which the corresponding effect is small but still distinguishable from GR [76], which reveals the different nature and effects of the two gravitational modifications.

Conclusions
Primordial black holes (PBH) can address a number of issues of modern cosmology since they may account for a part or all of the dark matter contribution, they can seed the large-scale structure formation through Poisson fluctuations, and they may constitute the progenitors of the black-hole merging events recently detected by LIGO-VIRGO.Interestingly, PBHs are tightly associated with gravitational wave (GW) signals, providing us the possibility to gain information of the physics of different cosmic epochs, from the very early universe up to later times, depending on the GW production mechanism.
Since their formation and effects are determined by the underlying gravitational theory, one can use them as a novel tool in order to test general relativity and investigate possible modified gravity deviations.
In this work we focused on the primordial scalar induced gravitational waves, generated at second order in cosmological perturbation theory, from PBH Poisson fluctuations, in the framework of f (T ) modified gravity.In particular, we desired to use it as a novel probe to extract constraints on the involved model-parameters.We considered three viable mono-parametric f (T ) models, and we investigated the induced modifications at the level of the gravitational-wave source, which are encoded in terms of the power spectrum of the PBH gravitational potential P Φ , as well as at the level of their propagation, described in terms of the Green function G k (η, η) which can be considered as the propagator of the tensor perturbations.
Our detailed analysis showed that within the observationally allowed range of the parameters of the f (T ) models at hand, the obtained deviations from GR, both at the level of source, as well as at the level of propagation, are practically indistinguishable.Indicatively, regarding the PBH gravitational potential power spectrum we found that the deviation from GR is of the order of 10 −2 for all the monoparametric f (T ) models considered here.This behavior is different than the case of other modified gravity theories, such as f (R) gravity, in which the corresponding effect is small but still distinguishable from GR [76].Hence, we conclude that realistic and viable f (T ) theories of gravity can safely pass the primordial black hole constraints, which may offer an additional argument in their favor.
Finally, one should stress that one can extend our analysis to other modified teleparallel theories of gravity, such as f (T, B) gravity and scalar-torsion theories, whose polarization numbers have been calculated in [111].Especially in the cases where extra polarization modes do appear, one expects to find significant differences, as it was found in [76].This interesting and necessary investigation will be performed in a separate paper.

Figure 1 :
Figure 1: The power spectrum of the PBH gravitational potential at PBH domination time, the three mono-parametric f (T ) models, by choosing different values of the initial PBH abundance Ω PBH,f .In all graphs we have used β = 0.1 and m PBH,f = 10 5 g.The dashed curves correspond to the GR results.

Figure 3 :
Figure 3: The power spectrum of the PBH gravitational potential at PBH domination time, for the three mono-parametric f (T ) models, by choosing different values of the parameter β within its observationally allowed range.In all graphs we have used m PBH = 10 5 g and Ω PBH,f = 10 −3 .The dashed curves correspond to the GR results.
We proceed by relating our solution for δ PBH from (3.31) with Φ, via the sub-Hubble scale approximation of the time-time field equation in f (T ) gravity for the PBH dominated era (equation (3.15)), which is: