Covariant density and velocity perturbations of the quasi-Newtonian cosmological model in $f(T)$ gravity

We investigate classes of shear-free cosmological dust models with irrotational fluid flows within the framework of $f(T)$ gravity. In particular, we use the $1 + 3$ covariant formalism and present the covariant linearised evolution and constraint equations describing such models. We then derive the integrability conditions describing a consistent evolution of the linearised field equations of these quasi-Newtonian universes in the $f(T)$ gravitational theory. Finally, we derive the evolution equations for the density and velocity perturbations of the quasi-Newtonian universe. We explore the behaviour of the matter density contrast for two models - $f(T)= \mu T_{0}(T/T_{0})^{n}$ and the more generalised case, where $f(T)= T+ \mu T_{0} (T/T_{0})^{n}$, with and without the application of the quasi-static approximation. Our numerical solutions show that these $f(T)$ theories can be suitable alternatives to study the background dynamics, whereas the growth of energy density fluctuations change dramatically from the expected $\Lambda$CDM behaviour even for small deviations away from the general relativistic limits of the underlying $f(T)$ theory. Moreover, applying the so-called quasi-static approximation yields exact-solution results that are orders of magnitude different from the numerically integrated solutions of the full system, suggesting that these approximations are not applicable here.


Introduction
In 1998 [1,2], observational evidence for an accelerating universe was discovered within the framework of a Friedmann-Lemaître-Robertson-Walker (FLRW) cosmology. The only way to explain this phenomenon is to introduce an additional dark component to the total energy density [3]. However, the physical properties of Dark Energy (DE) is still not well understood. Moreover, there are many problems that need to be explained such as the inhomogeneity problem and how the primordial fluctuations seeded the formation of structure on large-scales. There are exist several approaches to the theoretical description of these problems. One of these is a modified gravity theories, which provide the very natural gravitational alternative for dark energy and they are extremely attractive in the applications for late-time acceleration [4]. One of these modified gravity theories is f (T ) gravity [5][6][7][8][9][10]. In the simple case f (T ) = T , where T is the torsion scalar. The f (T ) theory can be directly reduced to the teleparallel equivalent of general relativity (TEGR) [11]. The idea of the teleparallel gravity (TG) was originally proposed in 1928 by Einstein after the formulation of general relativity (GR) to address the unification of gravity with electromagnetism by introducing the notion of tetrad (vierbin) field together with the suggestion of absolute parallelism [12][13][14][15]. In this theory, a set of tetrad fields e a (x) µ are considered to be the dynamical object instead of the metric g ab . While (GR) uses the well-known torsionless Levi-Civita connection, TG uses Weitzenböckconnection [16,17] that has no curvature but instead torsion. From Einstein-Cartan-Kibble-Sciama theory of gravitation [18], torsion is an alternative to curvature to describe the gravitational interaction [19,20]. From this view, the TG description of the gravitational interaction is completely equivalent to that of GR in some respects [21]. In (GR), curvature is used to geometrize the gravitational interaction, in other words, the gravitational force is replaced by geometry and particle trajectories are determined by geodesics, rather than the force equation [22]. In (TG), there is no notion of geodesics and the torsion yields a true gravitational force which quite similar to the Lorentz force of electrodynamics [17,23]. Although at the background and perturbation levels TG is completely equivalent to (GR). f (T ) gravity has various cosmological solutions which are consistent with the observational data [24][25][26], especially at the cosmological background. According to these features, f (T ) gravity is assumed as a viable theory both at cosmological and at astrophysical scales. One of the significant advantages of f (T ) gravity is that its field equations are always of second-order in contrast with f (R) gravity, where the field equations are governed by fourth-order equations. In [27], it was shown that T is not a local Lorentz scalar but the lack of this local Lorentz symmetry appears to be of little importance in (TG) when the Lagrangian is just T . However, it is not the case for the f (T ) generalisation of the (TG). If T is not a local Lorentz invariant then f (T ) cannot be either, and this is considered to be one of the disadvantages of the f (T ) gravity.
In this paper, we consider the covariant form of the field equations of f (T ) gravity to study linear cosmological perturbations [28]. There are two approaches to study the cosmological perturbations namely, the metric based approach [29][30][31][32][33][34][35][36][37][38] and the 1 + 3 covariant approach [39,40]. In the 1 + 3 covariant approach, the perturbations are formulated using variables that are covariantly defined in the real universe, and are exactly gauge-invariant by construction [41]. This approach has been used recently to study the cosmological perturbations for different contexts of modified gravity and GR [34,40,42].
This paper is organised as follows: in Sections 2, 3 and 4 respectively, we review the 1 + 3 covariant approach, kinematics quantities in the presence of torsion and we provide the covariant form of the field equations in f (T ) gravity which are required to study the cosmological perturbations. In Sections 5 and 6 respectively, we study the quasi-Newtonian models in the f (T ) gravity and we derive the integrability conditions that describes a consistent evolution of the linearised field equations of the quasi-Newtonian universes. In Section 7, we define the gradient variables that describe the cosmological perturbations and derive the linear evolution equations for matter and torsion perturbations. In Section 8, we analyse the growth of the matter density contrast by considering the power-law f (T ) theory where f (T ) = µT 0 (T /T 0 ) n and the more generalised case, where f (T ) = T + µT 0 (T /T 0 ) n by solving the whole system of perturbation equations numerically. We introduce the so-called quasi-static approximation to admit the approximated solutions on small scales. Section 9 is devoted for discussions and conclusions.

The 1 + 3 covariant approach in f (T ) gravity
In this approach, space-time is split into space and time, where 1+3 refers to the number of dimensions involved in each slice to investigate the deviation from homogeneity and isotropy of the Universe. The 4-velocity field vector of the observer u a associated with the worldlines is defined as where x a is the worldline in terms of local coordinates and τ is proper time measured along the worldlines. In this approach the metric g ab is decomposed into the projected tensor h ab as follows: We consider the covariant form of f (T ) gravity to clearly show the equivalence between teleparallel gravity and General Relativity. This form of field equation is very advisable to define the covariant variables in a gauge-invariant formalism for the study of the cosmological perturbations. Instead of using the torsionless Levi-Civita connection Γ c ab in general relativity, we use the curvature-less Weitzenböckconnection connectioñ Γ c ab in TG . The torsion tensor has different symmetry properties from the curvature case and it can be expressed as [22,43] Each vector e a is described by its components e µ a = 0, 1, 2, 3, in a coordinate basis. The contorsion tensor is expressed as [22] K c ab =Γ c ab − Γ c ab .
The teleparallel Lagrangian density is described by the torsion scalar as follows [22] where the super-potential term is given as [6,22,43] S ab and the contortion tensor can be rewritten as The modified teleparallel action for f (T ) is given by [22,25] where e is the determinant of the tetrad field e µ a i.e., (e = det|e µ a | = √ −g), and the coupling constant κ = 8πG/c 4 §. Varying the action in Eq. (8) with respect to the vierbein vector field e µ a , we obtain where f and f denote the differentiation with respect to T and Θ a A is the matter energy momentum tensor, and all aindices on the manifold run over 0, 1, 2, 3. From the relation between the Weitzenböckconnection and Levi-Civita connection in Eq. (4), one can write the Riemann tensor associated with the Levi-Civita connection and contorsion tensor as [22] and the Ricci scalar is given as The field equations can be written as where G ab = R ab − 1 2 g ab R is the Einstein tensor and g ab is the metric tensor. From Eq.
(12), we consider a covariant version of the field equations of f (T ) gravity with a clear analogy to Einstein's field equations as [5,22,28] f G ab where T ab denotes the usual energy-momentum tensor of the matter fluid expressed as a . § From here on-wards, the geometric units convention where 8πG = c = 1.
The linearised thermodynamic quantities in the presence of torsion are the energy density ρ T , the pressure p T , the energy flux q T a and the anisotropic pressure π T ab , respectively given by The total effective energy density, isotropic pressure, anisotropic pressure and heat flux of standard matter and torsion combination are defined as From Eq. (13), the Freiedmann equations of the effective fluid are presented as follows: One can directly obtain the corresponding thermodynamics quantities such as the effective energy density and the effective pressure of the fluid therefore, the Friedmann Eqs. (28) can be expressed as For simplicity, we have introduced the following dimensionless variables as presented in [10] are HereΩ m is the fractional energy density of effective matter like fluid, Ω m is the normalized energy density parameter of standard matter fluid and X is the fractional energy density of torsion fluid.
In the FLRW spacetime universe, the Raychaudhuri Eq. (21) that governs the expansion history of the Universe can be written as follows due to the resulting non-trivial field equations [49]:

Covariant equations
Given a choice of 4-velocity field u a , the Ehlers-Ellis approach [40,50] employs only fully covariant quantities and equations with transparent physical and geometric meaning [51]. In such a treatment for any scalar quantity X in the background we havẽ thus, by virtue of the Stewart-Walker Lemma [52], any quantity is considered to be gauge-invariant if it vanishes in the background. Therefore, the FLRW background is characterised by the following dynamics, kinematics and gravito-electromanetics [53]: where E ab and H ab are the gravito-electromagnetic fields responsible for tidal forces and gravitational waves. They are the "gravito-electric" and "gravito-magnetic" components of the Weyl tensor C abcd defined from the Riemann tensor R a bcd as The covariant linearised evolution equations in the general case are given by [51,53,54] These evolution equations propagate consistent initial data on some initial (t = t 0 ) hypersurface S 0 uniquely along the reference time-like consistency [55]. They are constrained by the following linearised equations [51,53,54]: These constraints restrict the initial data to be specified and they must remain satisfied on any hypersurface S t for all comoving time t.

Quasi-Newtonian spacetimes
There being no proper Newtonian limit for GR on cosmological scales, recent works on so-called quasi-Newtonian cosmologies [51,56,57] have shown that gravitational physics can be studied to a good approximation [36]. The importance of investigating the Newtonian limit for general relativity on cosmological contexts is that, there is a viewpoint that cosmology is essentially a Newtonian affair, with the relativistic theory only needed for examination of some observational relations. Most of the astrophysical calculations on the formation of large-scale structure in the universe rely on such a limit [57]. In [57], a covariant approach to cold matter universes in quasi-Newton has been developed and it has been applied and extended in [56] in order to derive and solve the equations governing density and velocity perturbations. This approach revealed the existence of integrability conditions in GR.
If a comoving 4-velocityû a is chosen such that, in the linearised form the dynamical, kinematic and gravito-electromagnetic quantities Eq. (33) undergo transformation.
Here v a is the relative velocity of the comoving frame with respect to the observer in the quasi-Newtonian frame, defined such that it vanishes in the background. In other words, it is a non-relativistic peculiar velocity. Quasi-Newtonian cosmological models are irrotational, shear-free dust spacetimes characterised by [51,53]: Therefore, the evolution equations (36) - (40) for this class of spacetimes can be written asθ Due to the vanishing of the shear in the quasi-Newtonian frame, Eq. (37) is turned into a new constraint and using the identity in Eq. (A.1) for any scalar ϕ, Eq. (38) can be simplified as where ϕ is the covariant relativistic generalisation of the Newtonian potential. The shear-free σ ab = 0 and irrotational condition ω = 0 and the gravito electromagnetic(GEM) constraint Eq. (42) result in the silent constraint H ab = 0 . Thus there is no gravitational radiation, which further justifies the term 'quasi − N ewtonian . And Eq. (44) show that q a is irrotational and thus v a it follows that for a vanishing vorticity, there exists a velocity potential ψ such that v a =∇ a ψ .

Integrability conditions
A constraint equation C A = 0 evolves consistently with the evolution equations in the sense that [51,53,56,58] where F and G depend on the kinematic, dynamical and gravito-electromagnetic quantities but not their derivatives [51,53,59]. It has been shown that the non-linear models are generally inconsistent if the silent constraint H ab = 0 is imposed, but that the linear models are consistent. Thus, so a simple approach to the integrability conditions for quasi-Newtonian cosmologies follows from showing that these models are in fact a sub-class of the linearized silent models. This can happen by using the transformation between the quasi-Newtonian and comoving frames.
The transformed linearised kinematics, dynamics and gravito-electromagnetic quantities from the quasi-Newtonian frame to the comoving frame are given as follows [51,53,59,60]: The first and second integrability conditions From Eq. (50), we need to ensure its consistent propagation at all epochs and in all spatial hypersurfaces. The differentiating this equation with respect to a cosmic time t and and together with Eqs. (26), (41) and (48) the identity of (A.2), one we obtain This equation is the first integrability condition for quasi-Newtonian cosmologies in f (T ) gravity and it is a generalisation of the one obtained in [56]. Eq. (57) reduces to an identity for the generalized van Elst-Ellis condition [51,53,57,59] aṡ From Eq. (47) with the time evolution of the modified van Elst-Ellis condition Eq. (58), we obtain the covariant modified Poisson equation in f (T ) gravity as follows as For the case of f = T , Eq. (59) reduces to the one obtained in [51]. By taking the gradient of Eq. (58), one gets using the identity (A.3), the above Eq. (60) can be written as using Eq. (51) together with the shear-free constraint Eq. (41), one can obtain the evolution equation of the 4-acceleration A a aṡ and it reduces to the one obtained in [51]. To check for the consistency of the constraint Eq. (50) on any spatial hyper-surface of constant time t, we take the divergence of Eq. (50) ad by using the identity (A.4), we get by using the constraint Eq.(43) with the identity (A.10) it follows that: By using Eqs. (26), (27) and (41), one obtains which is the second integrability condition, The left and right hand sides of Eq. (65) refer to GR and non-GR contributions. For GR, the right hand side vanishes when f = T and this matches the result obtained in [54] for the second integrability condition. By taking the gradient of Eq. (58) and using Eq. (41), one can obtain the peculiar velocity: which evolves according tȯ

Cosmological perturbations
In Section 6, we showed how imposing special restrictions to the linearised perturbations of FLRW universes in the quasi-Newtonian setting result in the integrability conditions. These integrability conditions imply velocity and acceleration propagation equations resulting from the generalised van Elst-Ellis condition for the acceleration potential in f (T ) gravity. In this section, we obtain the velocity and density perturbations via these propagation equations, thus generalizing GR results obtained in [51].

Definition of vector gradient variables
The covariant vector gradient variable D m a for the total matter fluid and the volume expansion of the fluid can be defined respectively as follows Those two gradient variables are the key to examine the evolution equation for the matter density perturbations. Now we define an extra key variables resulting from the spatial gradient of gaugeinvariant quantities which are connected with torsion fluid for f (T ) gravity, we define F a for the torsion density fluid and B a for the torsion momentum density respectively as follow.
We also define the comoving acceleration A a and the covariant vector gradient variable V m a for the velocity inhomogeneity of the matter as follows

Linear evolution equations
In the following we derive the linear evolution for each covariant variables. The system of equations governing the evolution of the variables defined in the previous subsection are given as follows

Definition of scalar gradient variables
The vector gradient variables we have defined so far contain both a scalar and a vector parts. However, the formation of structures on large scales follows a spherical clustering mechanism, and this could only be studied through the scalar parts. Therefore, we extract the scalar parts of the perturbation vectorial gradients of these quantities by applying a local decomposition [61] a∇aX = X ab = 1 3 where Σ X ab = X (ab) − 3h ab X describes shear whereas X [ab] describes the vorticity. When extracting the scalar contribution the vorticity term vanishes. Based on the above decomposition, by applying the comoving differential operator a∇ a to Eqs. (68)-(70), our scalar variables can be given to linear order as

first-and second-order evolution equations
Due to the above definitions of the scalar gradient variables, here we present firstand second-order evolution equations to demonstrate the growth of perturbations with cosmological red shift aṡ

Harmonic decomposition
The above evolution Eq. (79) -(97) can be thought of as a coupled system of harmonic oscillator differential equations of the form [36,62] X + Aẋ where A, B and C are independent of X and they represent friction (damping), restoring and source forcing terms respectively. To solve Eq. (88), a separation of variables is applied such that Since the evolution equations obtained so far are complicated to be solved, the harmonic decomposition approach is applied to these equations using the eigenfunctions and the corresponding wave number for these equations, therefore we write where Q k (x) are the eigenfunctions of the covariantly defined spatial Laplace-Beltrami operator [36,62], such that∇ The order of the harmonic (wave number) is given by k = 2πa λ , where λ is the physical wavelength of the mode. The eigenfunctions Q are covariantly constant, ieQ k ( x) = 0 . Therefore, the first-order evolution Eqs. (79) -(97) becomė Then, we will study the growth of the matter density contrast with cosmological redshift.
To do this we applied the transformation technique to make the redshift dependent instead of cosmic time. Therefore, our evolution equations can be written as follow: For more simplicity, we introduce here some quantities such as: By using the Friedmann Eq. (28) and the introduced dimensionless variables from Eq. (31), then we rewrite the more generalised form of the evolution equations as follows: For further analysis, in this part we are going to apply the quasi-static approximation to our evolution equations (114) -(116). In this approximation, we assume very slow temporal fluctuations in the perturbations of both the torsion energy density and its momentum compared with the fluctuations of the matter energy density. Therefore, terms involving time derivatives for torsion fluid are neglected, i.e., F = F ≈ 0.Then , the second-order evolution equations (114) -(116) for quasi-static approximations yield as In the above evolution equations (108) -(116) , we have to note that the k dependence appears in the first-order evolution equations and disappeasr in the second-order equations and this exactly the same as the work presented in [63] for GR. The GR can be recovered for the case of f (T ) = T , and we have In the following section, we explore the solutions of the density and velocity contrast in GR and f (T ) gravity models.

Solutions
In this section we will solve the whole system of perturbations equations we obtained so far Eqs. (108) -(116) to explore the growth of the matter density contrast in the GR context and for f (T ) gravity approach for non-quasi-static approximations and quasi-static approximation from Eqs. (117) -(119). The exact solutions of the matter density contrast can be found in the quasi static approximations and the numerical solution will be presented in the non-quasi-static approximations as well. To find those solutions, we consider the power-law f (T ) model where f (T ) = µT 0 (T /T 0 ) n and the more generalized model where f (T ) = T + µT 0 (T /T 0 ) n . These models produce the accelerating expansion of the Universe without invoking the cosmological constant. The cosmological and spherical solutions in f (T ) gravity lead to various viable models that support cosmological observations along with the solar system tests. In this section we investigate the matter density contrast and the velocity contrasts for non-quasi static and quasi-static approximations. We also show how these f (T ) gravity respond to the linear cosmological perturbations and formation of large-scale. We defined the normalized energy density for matter fluid as presented in [10] where ∆ in is the initial value of ∆ m (z) at z in = 1100, since the variation of CMB temperature detected observationally in the order of 10 −5 [64] at z ≈ 1100. In the same manner, we define normalized velocity contrast as

The growth of matter and velocity-density fluctuations in GR limits
Here, we analyse the growth of matter energy density contrasts δ(z) and the velocity contrast ν(z) with cosmic-time. We notice that, the second-order evolution equation (124) is an opened system. While, the evolution equation for velocity Eq. (125) is a closed system, and easy to construct the exact solution of the velocity contrast. Then the exact solutions of the velocity contrast yields as The integration constant c 1 and c 2 can be determined by the imposing initial conditions for plotting. Those constants are worthy to find the exact solution for the density contrasts, and we present in the following as Consequently, the second-order evolution equation of the matter density equation (124) becomes a closed system by substituting the solutions of the velocity contrast Eq. (128) and it's first order derivative. Then, the exact solution is given as where c 3 and c 4 are the integration constants and they are given as In the following figures we present growth the matter density and velocity fluctuations with cosmological redshift in the GR approach. We set the initial conditions at V in = V (z in 1100) = 10 −5 andV in =V (z in = 1100) = 0 and ∆ in = ∆ m (z in 1100) = 10 −5 and∆ in =∆ m (z in = 1100) = 0. For the case when f (T ) = T , δ(z) = δ GR (z) which coincides with TEGR and the results are exactly the same as GR. From this  plots, we depict clearly the contribution of dust component of the universe for the fluctuations of matter density and velocity are growing with decreasing red shift. In the following section we will consider two paradigmatic f (T ) gravity models to clearly see the contributions of both (dust and torsion) fluids for the growth of the fluctuations as well.

The growth of matter and velocity-density contrasts in f (T ) gravity models
For more simplicity, we first find the following parameters in f (T ) gravity models as Consequently, the evolution Eqs. (108)-(116) for f (T ) gravity model are given as (146)

f (T ) power-law model
In this subsection, we consider the paradigmatic power law model which is considered to be the simplest and compatible with the cosmic acceleration for n > 1.5 [65], it is given as where µand n are dimensionless constants, and in background cosmology we have T = −6H 2 and T 0 = −6H 2 0 is the present value of the torsion scalar. For n = 1 this model reduces to the GR limit. We assume the scale factor a(t) of the form where m = 2 3(1+w) is a positive constant with w = 0 and normalized coefficient a 0 and t 0 . The Hubble parameter is quantities X ,Ω m and Y as defined in Eq. (31) become To make sure that our evolution equations are dimensionless, we redefine the following normalized quantities as: Case I: Solving the whole system for the power-law f (T ) model In this sub-subsection, we present the numerical results of the velocity and matter density fluctuations for power-law f (T ) gravity. We start by solving the whole system of the first-order evolution equations (138) -(143). We have evaluated the numerical solutions simultaneously to analyze the density fluctuations. We set the initial conditions at V in = V (z in 1100) = 0, F in = F(z in 1100) = 0 , B in = B(z in 1100) = 0, A in = A(z in 1100) = 0, Z in = Z(z in 1100) = 0 and ∆ in = ∆ m (z in 1100) = 10 −5 . The numerical results are presented in the following Figs. 3 -4. For the case of n = 1, the numerical results of GR are recovered. We have noticed the growth of the density contrast for values of n > 1.5, see Fig. 4 and what we noticed here is that only at this specific choice of the initial conditions, the k dependence does not make any difference in the behaviour or the amplitude of the density contrast for long-and short wavelength (different values of k) . In the following, we will solve the whole system of the second-order perturbation equations, we have evaluated the numerical solutions from Eqs. (144) -(146) simultaneously to analyze the density and velocity fluctuations with redshift. For further analysis, we study the growth of the matter density fluctuations with redshift and present the numerical results for different initial conditions as presented in [66], to see how sensitive the results are to change the initial conditions:  As we mentioned in the earlier and as presented in [65] this model has an observational valid only for the range of n > 1.5 for the background universe, However, we have tried to study the behaviour of the velocity and the density contrast for different n ranges, i.e., for n < 1, n closes to GR and n > 1.5. For instance, the numerical results for set I of the initial conditions, we have noticed that the velocity and density contrast are decaying for values of 0.5 < n ≤ 0.99 and 1 < n ≤ 2, and the decaying is faster with increasing the values of n. For values of n > 2 we notice the growth of the velocity and the density contrast which is realistic compared to the observational expectations. The behaviour of the growth of density contrast is summarised in Table. 1. For set II of the initial conditions, we observe the growth of the velocity and density contrast for values of n > 1.5, while we have to note that the amplitudes of the velocity and density contrast are decreasing with increasing the values of n, we also have to point out that we have noticed orders-of-magnitude deviations from limiting general relativistic results. We present the behaviour of the density fluctuations in Table. 2 We evaluated Table 1: The behavior of δ(z) for set I.
Range of n Density contarst Expected Remark 0.5 < n ≤ 0.99 decreasing decreasing realistic 1 ≤ n ≥ 2 increasing increasing realistic the numerical results of the whole system by using the initial conditions sets I and II. Based on these results we can conclude that the growth of the velocity and density contrast is very sensitive to the values of n and the initial conditions. The results obtained by using conditions I are in complete disagreement with the theoretical and While ψ = − 3 2 and Y = −2(n − 1). By choosing set I of the initial conditions, We notice the growth of the density and velocity fluctuations increases with n but they are highly nonlinear compared to the GR results and the GR results are recovered as shown in Figs. 1 and 2.

Case II: Quasi-static approximations
In this sub-subsection, we present the results of the velocity and density contrast in the quasi-static approximations. We notice that Eq. (119) is an opened system and Eq. (118) is a closed system. We apply the same technique as GR limits and we find first the exact solution of velocity contrast and then the density contrast. So, the exact solution of Eq. (118) is given as where where η 1 = Hη and ζ 1 = 1 H 2 ζ. After we computed the integration constants c 5 and c 6 by imposing the initial conditions, the exact solution for density contrast Eq. (119) is The background quantities X ,Ω m and Y as defined in Eq. (31) becomẽ .

(157)
Case I: Solving the whole system for the generalised f (T ) model Based on the above definitions, the numerical results of this system of the first-order equations (138) -(143) are presented in such way that we set our free parameter Ω m = 0.32 based on the observational expectations. We set the initial conditions at As we mentioned earlier only at this specific choice of the initial conditions, the k dependence does not make any difference in the behaviour or the amplitude of the density contrast. The numerical results of this system of Eqs. (138) -(143) are presented in Figs. 7 -8, For the case of n = 0, the numerical results of GR are recovered. We have noticed the growth of the density contrast for values of n close to 0 for different values of k and as we mentioned before the k dependence does not make any difference in the behaviour or the amplitude of the density contrast only at this specific choice of the initial conditions.  In the following we will solve the whole system of the second-order perturbation equations, evaluating the numerical solutions from Eqs. (144) -(146) simultaneously to analyse the density and velocity fluctuations as a function of redshift. For further analysis, we evaluated the numerical results for different initial conditions as we did for the power law model. The numerical results of this system of Eqs. (144) -146 for set I of the initial conditions , the GR results are recovered as in Figs. 1 and 2 for n = 0 and we noticed the growth of the matter density and velocity for very small values of n but the results are still highly nonlinear compared to the GR results.
As we mentioned this model is tested to be valid for very small values of n 1, we tried set II of the initial conditions and we noticed the decay of the velocity and the density fluctuations for values of n 1. We can conclude based on the results obtained that this model is viable for our study only if we choose the initial conditions for the first derivatives of the variables ∆ m , V and F to be zero at the initial red-shift. In other words, set II of the initial conditions does not provide good results for this model. We also found that the exact solutions in the quasi-static approximations for this model could not be evaluated.

Conclusions
This work presented a detailed analysis of scalar cosmological perturbations in the f (T ) gravity theory using the 1 + 3 covariant gauge-invariant approach. We explored the integrability conditions of the so-called quasi-Newtonian cosmological models in the context of f (T ) gravity, a first such study to the best of our knowledge. We showed that for such cosmological models to exist, they must satisfy certain integrability conditions on the generalised Einstein field equations. The two integrability conditions derived and presented here allow us to describe a consistent evolution of the linearised field equations of quasi-Newtonian universes. We defined the gauge-invariant variables and derived the corresponding evolution equations. We derived the complete set of the firstand the second-order evolution equations of these perturbations These results agree with the GR results when f (T ) = T . We studied the behaviour of matter energy density perturbations with redshift for different ranges of the model-defining parameter n by considering two of the f (T ) gravity models. We employed the so-called quasistatic approximation technique on small scales. Some of the specific highlights of this work are as follows: in the first model f (T ) = µT 0 ( T T 0 ) n , we presented the ranges of values of n for which the perturbation amplitudes δ(z) grow or decay. For instance, the numerical solution of the first-order perturbation equations shows the growth of the density contrast for values of n > 1.5 only at one particular choice of initial conditions and that the wave number k has no effect on the behaviour of the density contrast (see Fig. 4). We also found numerical solution of the whole system of second-order equations and for different sets of initial conditions. For instance, for set I we observed unrealistic behaviour of the density perturbations for 0.5 ≤ n ≤ 0.99 and for 1 < n ≤ 2, but for values of n > 2, we found that the growth of the density contrast was consistent with observational expectations. The results obtained in the first model generalised existing GR results as it is shown in Figs. 1 and 2. We also noticed that the behaviour of the density fluctuations is highly dependent on our choice of the initial conditions. For instance, for set II of the initial conditions, we found results which are highly non-linear compared to expected GR results. The results obtained using the first set of initial conditions (Set I) are in complete disagreement with the theoretical and observational expectations. On the other hand, applying the quasi-static approximation to his model, the growth of the density perturbations can be explored for the range of n > 1.5 for Set I of the initial conditions. By comparing the complete system with the quasi-static approximation, for the same values of n and the same choice of the initial conditions, we found that there is a big difference in the amplitudes of the density fluctuations. In conclusion, the quasi-static approximation appears to be not applicable for this model.
In the second case: f (T ) = T + µT 0 (T /T 0 ) n , with the value of n constrained by solar system tests to be in the range n 1, we obtained a growth in the velocity perturbations and the density contrast δ(z) for 0 ≤ n ≤ 0.4. We considered a range of possible initial conditions, and found that this model is only viable when the first derivatives of the variables ∆ m , V and F are equal to zero. In conclusion, we have developed a covariant framework for studying quasi-Newtonian cosmologies in the f (T ) gravity theory. Future studies in this direction using more realistic (i.e., in terms of cosmological and astrophysical viability) models of the theory using actual astronomical data might shed more light on whether such a framework can be a viable cosmological alternative to late-time cosmology.