{'fraction_non_alphanumeric': 0.034536053504187036, 'fraction_numerical': 0.015343112459956162, 'mean_word_length': 4.322614418187257, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 11, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
General Perfect Fluid Perturbations of Homogeneous and Orthogonal Locally Rotationally Symmetric Class II Cosmologies

First order perturbations of homogeneous and hypersurface orthogonal LRS (Locally Rotationally Symmetric) class II cosmologies with a cosmological constant are considered in the framework of the 1+1+2 covariant decomposition of spacetime. The perturbations, which are of perfect fluid type, include general scalar, vector and tensor modes and extend some previous works where vorticity perturbations were excluded. A harmonic decomposition is performed and the field equations are then reduced to a set of eight evolution equations for eight harmonic coefficients, representing perturbations in density, shear, vorticity and the Weyl tensor, in terms of which all other variables can be expressed algebraically. This system decouples into two sub-systems, one for five and one for three coefficients. As previously known, vorticity perturbations cannot be generated to any order in a barytopic perfect fluid. Hence the time development of existing first order vorticity perturbations are seen to be completely determined by the background. However, an already existing vorticity will act as source terms in the evolution equations for the other quantities. In the high frequency approximation the four independent Weyl tensor harmonics evolve as gravitational waves on the anisotropic background in the same manner as in the case without vorticity, whereas vorticity gives a first order disturbance of sonic waves.

First order perturbations of homogeneous and hypersurface orthogonal LRS (Locally Rotationally Symmetric) class II cosmologies with a cosmological constant are considered in the framework of the 1+1+2 covariant decomposition of spacetime. The perturbations, which are of perfect fluid type, include general scalar, vector and tensor modes and extend some previous works where vorticity perturbations were excluded. A harmonic decomposition is performed and the field equations are then reduced to a set of eight evolution equations for eight harmonic coefficients, representing perturbations in density, shear, vorticity and the Weyl tensor, in terms of which all other variables can be expressed algebraically. This system decouples into two sub-systems, one for five and one for three coefficients. As previously known, vorticity perturbations cannot be generated to any order in a barytopic perfect fluid. Hence the time development of existing first order vorticity perturbations are seen to be completely determined by the background. However, an already existing vorticity will act as source terms in the evolution equations for the other quantities. In the high frequency approximation the four independent Weyl tensor harmonics evolve as gravitational waves on the anisotropic background in the same manner as in the case without vorticity, whereas vorticity gives a first order disturbance of sonic waves.

I. INTRODUCTION
Present observations of cosmological parameters, like the large scale structure, anisotropies of the background radiation (CMB) and the cosmological redshift are well described by an expanding isotropic and homogeneous background metric, where the time development is sourced by a cosmological constant and a fluid consisting of baryonic matter, radiation and dark matter, the ΛCDM model [1,2]. However, there are deviations. For example, the predicted amplitudes in the power spectrum of the CMB differ from the observed for large angles [3][4][5][6], even though the uncertainties are high for these angles, and the upper limits of a possible anisotropy in the cosmological redshift are rather big [7][8][9]. Hence it is of interest to investigate the properties of other models, like, e.g. the evolution of perturbations on anistropic backgrounds. If the background has one spatial direction of anisotropy, so that it is Locally Rotationally Symmetric (LRS), a gauge invariant perturbation theory built on the 1+1+2 covariant split of spacetime [10][11][12], which generalizes the 1+3 covariant split [13][14][15][16][17], is suitable for the study. For other perturbative methods used in general relatvity see e.g. [18][19][20][21][22] and for some previous works on perturbations in anisotropic universes see e.g. [23][24][25][26][27][28].
In some earlier papers [29][30][31] we studied vorticity free perturbations on a class of anisotropic cosmological models, given by the homogeneous and orthogonal locally rotationally symmetric (LRS) spacetimes of class II [32,33]. The restriction to zero vorticity was partly due to the significant simplifications this implies for the commutator relations, but was also motivated by that vorticity cannot be generated in a perfect fluid with barotropic equation of state, cf. e.g. [14,34,35] and that an exist-ing vorticity in a homogeneous and isotropic universe decays with time [21]. However, the vorticity acts as source terms for the other quantities and hence can give imprints on other measurable quantities. For works on vorticity perturbations on Robertson-Walker backgrounds see e.g. [21,34,[36][37][38][39].
In this paper we present a general treatment of a class of perturbed LRS cosmological backgrounds, with nonvanishing vorticity on the perturbed model. We start by assuming a homogeneous and hypersurface orthogonal cosmological background model which is of LRS class II, characterized by that the vorticity, magnetic part of the Weyl tensor, and twist of the 2-sheets all vanish [32]. The matter content of the model is assumed to be described by a perfect fluid, and a nonvanishing cosmological constant. The choice of an anisotropic background gives us a preferred spatial direction, which motivates our choice to utilize the 1+1+2 covariant formalism [10], in performing the perturbation. This formalism entails first splitting spacetime in the 1+3 covariant formalism [13,14], with the average velocity of matter being the preferred direction for each point. After this split we are left with 3-surfaces in spacetime which we further decompose, this time with the anisotropy as the preferred direction, and this gives us all relevant quantities in the 1+1+2 covariant formalism.
In order for the perturbation to be gauge invariant we define all our variables in terms of quantities that vanish on the background. Due to the Stewart-Walker lemma [40] these variables are ensured to be gauge invariant. With the background and perturbation method properly defined, we proceed to linearize the propagation equations and constraints given in the 1+1+2 formalism, and then decompose the gauge invariant variables using a harmonic decomposition. Our initial system, consisting of scalar, vector, and tensor equations, is then reduced to a system containing only scalar time evolution equations and constraints. This system is then solved and the result is a set of eight harmonic coefficients, out of which all other variables can be expressed, and their associated time evolution equations.
The final set of harmonic coefficients decouples into an even and odd sector, containing five respectively three variables. The vorticity appears in both these sectors, with one degree of freedom in each, and its time evolution equations completely decouples from the other variables. Assuming a linear equation of state, we get an analytical solution for the vorticity, which depends on time through the expansion coefficients, and on the comoving wavenumbers through integration constants dependent on the initial conditions. The consistency of this solutions is checked specifically through the fall off rate for the case of dust and radiation, and the conservation of angular momentum.
Even though the evolution of the vorticity is decoupled from the other variables, the vorticity appears as a source term in the time evolution of five of the remaining six harmonic coefficients. We move on to examine these evolution equations in the geometrical optics approximation by deriving second order wavelike equations for the harmonic coefficients in the limit of large harmonic numbers k. On keeping terms to the two highest orders in k, damped shear waves are then found to propagate with the speed of sound and to be sourced by vorticity. The evolution equations for the four components of the Weyl tensor are to this order unaffected by the vorticity and represent gravitational waves.
The paper is organized as follows: In section II the theory behind the 1+3 and 1+1+2 covariant splits of spacetime are briefly summarized and in section III the cosmological backgrounds are described. The perturbative method is described in section IV and harmonic decompositions of the first order quantities are shown. In section V the two reduced systems of evolution equations are given and in section VI the equations for the vorticity perturbations are solved and discussed. Finally, in section VII the high frequency approximation is considered and then the results are summarized in section VIII. The commutation relations between the different differential operators are given in appendix A and properties of the harmonics in appendix B. The harmonic expansion of the linearized equations are presented in appendix C and then 27 harmonic coefficients are solved for algebraically in terms of the 8 remaining in appendix D.
II. THE 1 + 3 AND 1 + 1 + 2 COVARIANT SPLITS OF SPACETIME We will in this paper utilize the 1 + 1 + 2 covariant split of spacetime, since we have two preferred directions on the background. This method builds on the 1 + 3 covariant split of spacetime, and in this section we will provide a brief summary of both these formalisms. More details about the 1 + 3 formalism can be found in [13,14], and about the 1 + 1 + 2 formalism in [10,11].
We start by looking at the 1 + 3 covariant split of spacetime, which is a suitable method when there exists a preferred unit timelike vector u a . We use this vector u a to define a projection tensor h b a ≡ g b a + u a u b , where g ab is the 4-metric, which projects onto the 3-surfaces orthogonal to u a . This projection tensor allows us to covariantly decompose vectors and tensors into spatial and timelike parts. Using u a and h ab we can define the covariant time derivative and projected spatial derivative asṪ c···d a···b ≡ u e ∇ e T c···d a···b (1) and respectively.
The independent variables of relevance are now given by the decomposed form of the Ricci tensor, the Weyl tensor and the covariant derivative of u a . The Ricci tensor can be expressed, using Einstein's field equations, through the energy momentum tensor which is decomposed as T ab = µu a u b + ph ab for a perfect fluid. Here, µ ≡ T ab u a u b is the energy density and p ≡ 1 3 T ab h ab is the isotropic pressure of the fluid. The Weyl tensor is decomposed into an electric and magnetic part, E ab ≡ C acbd u c u d and H ab ≡ 1 2 ε ade C de bc u c , respectively, where we have introduced the volume element on the 3- Lastly, the covariant derivative of u a is decomposed as where we have the acceleration A a ≡ u b ∇ b u a , rate of volume expansion Θ ≡ D a u a , rate of shear σ ab ≡ D a u b and vorticity ω ab ≡ D [a u b] . Here the angular brackets denotes the orthogonally projected tracefree symmetric part of a tensor, defined as T ab ≡ h c (a h d b) − 1 3 h ab h cd T cd , and the square brackets denotes the antisymmetric part of a tensor. We also further define a vorticity vector as ω a ≡ 1 2 ε abc ω bc . The evolution in the direction of u a and the constraints for the independent variables introduced above are now given by the Ricci identity for u a , the once and twice contracted Bianchi identities and the commutation relations of the new projected derivatives, and can be found in e.g. [14].
We can further decompose our physical quantities in what is known as the 1+1+2 covariant formalism. This uses a preferred spatial direction, given in terms of a unit vector n a which is orthogonal to u a , i.e. it resides on the 3-surfaces. As previously, we introduce a projection tensor N b a = h b a − n a n b which will project vectors and tensors onto the 2-sheets orthogonal to n a . This allows us to decompose any spatial vector into a scalar part along n a and a vector part on the 2-sheet, and any spatial tensor into a scalar part along n a , a vector part on the 2sheet and a projected, trace-free and symmetric 2-tensor. We also have the derivative along n a and the projected derivative on the 2-surfaces defined as T c···d a···b = n e D e T c···d a···b (4) and respectively.
Decomposing the kinematic quantities given in the 1 + 3 formalism gives us A a = An a + A a , ω a = Ωn a + Ω a and σ ab = Σ n a n b − 1 2 N ab +2Σ (a n b) +Σ ab . The electric and magnetic part of the Weyl tensor, E ab and H ab , are decomposed in the same way as the shear, in terms of the variables E, E a , E ab and H, H a , H ab . We also get new kinematic quantities by decomposing D a n b andṅ a as D a n b = n a a b + 1 2 φN ab + ξε ab + ζ ab andṅ a = Bu a + α a , where we have the acceleration of n a , a a ≡ n c D c n a , the 2-sheet expansion φ ≡ δ a n a , the twisting of the 2sheets ξ ≡ 1 2 ε ab δ a n b , the shear of n a , ζ ab ≡ δ {a n b} , and α a ≡ N b aṅb . It follows, since u a n a = 0, that B = A. In the above decomposition we have introduced the volume element on the 2-sheets ε ab ≡ ε abc n c , and also the shorthand notation of curly brackets to denote the projected, trace-free and symmetric part of a 2-tensor. We will also introduce the notation of a bar over a vector index to denote that index to be projected onto the 2-sheets, i.e. ψ a ≡ N b a ψ b . The Ricci identities for u a and n a and the Bianchi identities can now be decomposed into evolution and propagation equations in the direction of u a respectively n a , and constraint equations. For the full set of these equations, see e.g. [10]. We also get commutation relations between the operators˙,ˆand δ a which can be found in appendix A.

III. BACKGROUND SPACETIMES
The properties of the background spacetimes are here breifly summarized. For more details the reader is referred to [30,32,33].
As backgrounds we will use the class of homogeneous, hypersurface orthogonal and locally rotational symmetric (LRS) perfect fluid solutions to Einstein's equations with vanishing magnetic part of the Weyl tensor. For technical reasons we exclude the hyperbolic and closed Friedmann models, which are different from the others both in form of line element and harmonic decomposition [30]. In the classification of [32] the used backgrounds belong to LRS class II.
These spacetimes have two prefered directions, given by the 4-velocity of the fluid, u a , and the direction of anistropy, n a . The 2-sheets perpendicular to n a are maximally symmetric with 2D curvature scalar R = 2K/a 2 2 , where a 2 (t) is the radius of curvature and K = ±1 or 0 for spheres, pseudo-spheres or planes, respectively. In terms of the quantities of the 1+1+2 formalism introduced in section II, these spacetimes are characterized by the following nonzero quantities On assuming a nonzero cosmological constant Λ, they satisfy the following set of evolution equationṡ and E is given algebraically as Hence a closed system is obtained if an equation of state p = p(µ) is given. The 2D curvature scalar may be expressed in terms of the quantities in S (0) as The line-element can for this class of metrics be written as in terms of the two scale factors a 1 (t) and a 2 (t). For spheres (K = 1) f 1 (ϑ) = sin 2 ϑ, for pseudo-spheres (K = −1) f −1 (ϑ) = sinh 2 ϑ and for planes (K = 0) f 0 (ϑ) = 1 (or f 0 (ϑ) = ϑ 2 ). The expansion and shear are then given by In some previous papers we have studied perturbations on the backgrounds given in section III [29,30]. Here these studies will be extended to the most general first order perturbations, consistent with a perfect fluid, by including also vorticity of the fluid. To avoid the gauge problem in perturbation theory, i.e. the problem of identifying the perturbed spacetime with the background, we will use an approach based on the covariant split of spacetime and quantities which are zero on the background.
As variabels we will use the scalars which are nonzero on the background, i.e. G = {µ, p, E, Θ, Σ}, and the remaining 1+1+2 covariantly defined quantites from section II, which are zero on the background. The latter are gauge-invariant due to the Stewart-Walker lemma [40].
To get a consistent first order system we need also expand the quantities G as can in a convenient way be represented by the gradients of the G, i.e. by δ a G andĜ. These are also zero on the background for which the quantities only are functions of time. As will be shown later in section IV C, where a harmonic decomposition of the first order quantities is done, theĜ derivatives may be expressed in terms of the derivatives Z a ≡ δ a G. Hence, instead of G (1) , we introduce the following first order quantities Z a which vanish on the background. The complete set of first order quantities which vanish on the background is then given by The frame is partly locked by using the 4-velocity of the fluid as the preferred timelike vector u a also for the perturbed spacetime. The direction n a of anistropy is not well-defined in the perturbed spacetime and we will later restrict the frame, see section V. On imposing the Ricci identities for the vectors u a and n a and the Bianchi identities together with the commutator relations (see appendix A), evolution equations along u a , propagation equations along n a and constraints are obtained for the covariant quantities. The exact equations are found in [10]. The linearized equations for the first order quantities are given below.

A. Linearized Equations
The evolution equations arė The equations containing a mixture of evolution and propagation contributions are The equations containing only propagation contributions are Lastly, the constraints are

B. Harmonic Decomposition
To transform the linearized partial differential equations from the previous section into ordinary differential equations in time we perform a harmonic decomposition of the covariant quantities. The procedure is here shortly described and for more details the reader is referred to [11,30,41,42]. Different types of harmonics used in relativity can, for example, be found in [43][44][45][46].
Scalars are written as where the harmonic coefficients Ψ S k k ⊥ are functions of time. Here the eigenfunctions P k satisfy where the constants k are dimensionless comoving wave numbers in the n a direction and the physical wavenumbers are given by k /a 1 in terms of the scale factor a 1 = a 1 (t). In terms of the metric (14) the P k can be represented by e ik z . The eigenfunctions Q k ⊥ satisfy [10]: where a 2 is the scale factor of the 2-sheets, and k ⊥ are the dimensionless comoving wavenumbers along the 2sheets. When R > 0 the 2-sheets are spheres and the eigenfunctions can be represented by the spherical harmonics Y m l with k 2 ⊥ = l(l + 1), l = 0, 1, 2, .... Due to the background symmetry the m-values will not appear explicitly in the equations. When R ≤ 0, and the 2-sheets are open, the k ⊥ take continuous values. For continuous values the sums go over into integrals.
Vectors Ψ a are expanded in terms of the even and odd vector harmonics [11,41,42] as Note the factor a 2 appearing in the definitions of the vector harmonics in (55), giving the harmonic coefficients Ψ V k k ⊥ and Ψ V k k ⊥ the right physical dimensions. Similarly, a tensor Ψ ab can be expanded in terms of the even and odd tensor harmonics as Some properties of the vector and tensor harmonics are found in appendix B.

C. Relations between harmonic coefficents
As mentioned above, the hat derivatives,Ĝ, of the elements G = {µ, p, E, Θ, Σ} may be expressed in terms of Z a ≡ δ a G. ExpandĜ and Z a as and respectively. From the even part of the commutation relation (A3) one then obtains Hence allĜ are given by (59) withG S k k ⊥ from (61). Similarly, with the commutation relation (A4) and the properties of the harmonics in appendix B, one can show that the odd coefficients of Z a satisfy The degrees of freedom of the vorticity can also be decreased by substitution of the harmonic expansions of Ω and Ω a in equation (41), which gives In the following sections and the appendicies we will drop the indicies k k ⊥ on harmonic coefficient, since their meaning will be obvious due to the superscripts S, V and T .

V. EVOLUTION EQUATIONS
By decomposing the system of linearized equations given in section IV into harmonics we get a new system of scalar evolution equations and constraints, which can be found in appendix C. The new system contains a total of 30 evolution equations, 28 constraints and 37 scalar variables, and it decouples into an even and odd sector.
Before solving the new system of equations we can simplify it further. First, we note that we still have the freedom to partially fix the frame on the perturbed model, in choosing a direction for n a . Based on the reasoning in [29] we choose n a such that a a = 0, which will simplify the system significantly. Also, we choose to look specifically on barotropic matter perturbations, which means that p = p(µ), and hence that p V = c 2 s µ V , where c 2 s =ṗ/μ is the square of the matter speed of sound.
With the choices above the even sector can now be solved in terms of the five free variables and the odd sector in terms of the following three free variables, Note that for vectors and tensors which are odd by definition, like the magnetic part of the Weyl tensor, the rôles of quantities without and with an overbar are interchanged. Hence, H T belongs to the odd sector, whereas H T belongs to the even sector.
It is of course possible to choose eight different variables than the ones above, though our choice is motivated by the results in [29,30], where E and H represents gravitational waves. The remaining even variables, can now all be expressed in terms of the previously introduced free variables. These expressions are found in appendix D, and their consistency has been checked against the propagation equations found in appendix C. The equations of interest are now the evolution equations for the variables given in (D1) and (D2). For the even sector we get the following equations, Above, we have introduced the additional notationsk For the odd parity sector we get the following evolution equations,Ω where we have introduced Interestingly enough, the evolution of both the even and odd part of the vorticity completely decouples from the other variables, and hence we will treat it separately in the next section. The evolution equations for the remaining variables are consistent with previous research, see e.g. [29,30], with the addition of the vorticity appearing as a source term for all variables except Σ T . These equations are still quite cumbersome to solve analytically though, and so we will proceed to examine them in the high frequency limit in section VII.

VI. VORTICITY
As is well known both in the non-relativistic, cf. e.g. [35], and relativistic cases, vorticity cannot be generated in a perfect fluid with barotropic equation of state p = p(µ), see e.g. [14,34] for a proof in general relativity. Their argument is in short as follows: in the 1+3 covariant split the evolution equation for the vorticity vector ω a is given by [14]. Using the twice contracted Bianchi identities D a p + (µ + p)A a = 0 andμ + Θ(µ + p) and the commutator relation where Ψ is an arbitrary scalar, it follows that the last term can be rewritten as c 2 s Θω a when p = p(µ), so that Hence there are no source terms, and since this equation is exact, vorticity cannot be generated to any order in a barotropic perfect fluid. For generation mechanisms due to entropy flow see e.g. [36,37] or from N-body simulations see [47]. However, as seen from section V, an already existing vorticity will act as source terms in the evolution equations for the density and electric and magnetic parts of the Weyl tensor, and may hence influence the growth of structures and give imprints on sonic and gravitational waves. For high frequency gravitational waves the effects of vorticity will be negligible though, as will be seen in section VII.

A. Evolution of vorticity
On projecting (87) along n a and onto the 2-sheets with N a b , using that ω a = Ωn a + Ω a , we obtaiṅ anḋ 89) respectively. To first order we then havė These are in agreement with the scalar equations (78) and (66) for the harmonic coefficents Ω S and Ω respectively. By assuming a linear equation of state, p = (γ − 1)µ, and substituting equations (15) and (16) for Θ and Σ in terms of the scale factors a 1 and a 2 , one can now easily solve for either Ω and (Ω a Ω a ) 1/2 or the corrresponding harmonic coefficients Ω S and Ω V 2 . The harmonic coefficients are given by and where C Ω = C Ωk k ⊥ and C Ω = C Ωk k ⊥ are functions of the harmonic numbers k and k ⊥ , given by the initial conditions. For the dust case (γ = 1) and radiation (γ = 4/3) we then get and 1 A scalar equation for Ω a Ωa can be obtained by multiplying (91) with Ωa giving (Ω a Ωa)˙= − Σ + 4 3 − 2c 2 s Θ Ω a Ωa (for a scalarΨ reduces to ∂Ψ/∂t in the given background coordinates (14)). 2 Due to the definition of the vector harmonics, see (55), the harmonic coefficients of the vorticity have the right dimension and time dependence.
respectively. To compare with earlier results for the flat Friedmann universes [21] let us first consider the isotropic case with a 1 = a 2 ≡ a. Then Ω S = C Ω a 3γ−5 and Ω V = C Ω a 3γ−5 (98) so that for dust and radiation the vorticity goes as a −2 and a −1 , in agreement with the results in [21]. As a check of the result for the vorticity, we here consider the angular momentum to first order. Consider a small comoving volume element with uniform angular velocity, which is proportional to the vorticity. Since it is a first order quantity, the other quantities can be treated to zeroth order. Note also that to this order the Lorentz factor is one. Due to the background isotropy around n a , the component of the angular momentum along n a should be conserved in the general case, as all its components should in the fully isotropic case.
From equations (9) and (15) with p = (γ − 1)µ the equation is obtained for some C µ . Volume increases as V ∼ a 1 a 2 2 and radius perpendicular to the n a axis as R ∼ a 2 . Hence the time dependence of the component of the angular moment parallel to n a , L , of a comoving element with massM goes as i.e., it is independent of time. For rotation around an axis in the xy-plane perpendicular to n a , we may choose this axis to be, say, the x-axis due to the background isotropy. Then the radius perpendicular to the axis of rotation will go as R ∼ a 1n + a 2ŷ (modulo non time-dependent coefficents in front of the two different directions), so that the time dependence of R 2 ∼ a 2 1 + a 2 2 . The time dependence of the component of the angular moment perpendicular to n a , L ⊥ , then goes as and hence L ⊥ is conserved only if a 2 ∼ a 1 , i.e. if the background is isotropic.

B. Vorticity as source terms
From equations (78), (79) and (80) for the odd parity sector we can derive the following inhomogeneous wave equations with damping for E and Ḧ andḦ where Hence vorticity, which as seen is independent of the other quantities, acts as a source for the E and H perturbations. However, as will be seen in section VII, for large wave numbers, k , k ⊥ , the source terms become negligible. Similarly we can from the even parity equations (66), (67), (68), (69) and (70) derive wavelike equations sourced by the density gradient µ V , its time derivativė µ V and the vorticity component Ω V for the shear Σ T and the components E T and H T of the Weyl tensor. Note however that to obtain pure equations in Σ T , E T or H T , sourced only by the vorticity, one needs to go to higher order differential equations than second order.
For the shear one gets where for the magnetic part of the Weyl tensor where q H0 = 3 2 and for the electric part of the Weyl tensor where with the definitions Their behaviour for large wave numbers will be considered in section VII.

VII. GEOMETRICAL OPTICS APPROXIMATION
We now turn our attention to the high frequency limit, also known as the geometrical optics approximation. Following the definitions given in [48,49], our current variables can in the limit be expressed as where the physical wave numbers along the direction of anisotropy are given by k /a 1 and along the 2-surfaces by k ⊥ /a 2 . Hence we define the total physical wave number k as Our results will now be compared to previous results, see [29] and [30], where the only differences should be related to the impact of Ω S and Ω V . As we saw in the previous section though, the evolution of the vorticity completely decouples from the other variables, independent of the magnitude of the frequency. However, we will see that the vorticity appears as a first order source term in the evolution equation for the shear. We now rewrite the second order equations from section VI B in the form of the general wave equation for a time dependent scalar Ẍ X + 2ζωẊ + ω 2 X = Z where ζ represents the damping parameter, ω represents the undamped angular frequency and Z is some time dependent factor which acts as a source term, producing forced oscillations. On assuming that the change of the scale factors are negligible over one period, we have that the real angular frequency is obtained from ω 1 − ζ 2 with |ζ| < 1, which means that the propagation speed of the wave is given by ω 1 − ζ 2 /k. First we consider the even sector, given by the equations (67), (68), (69) and (70). By keeping only the highest order terms in k they produce undamped wave equations, where the density waves µ V and the shear waves Σ T progagate with the speed of sound c s , and the Weyl tensor components E T and H T propagate as gravitational waves with the speed of light. To next order in k also damping and source terms are introduced. The equations (118) and (113) for E T and H T then reduce tö where Hence, H T and E T represent damped gravitational waves which to zeroth and first order in 1/k propagates at the speed of light. This can be seen as follows: From equation (131) for E T we can read off the undamped angular frequency ω to be k. Hence the propagation velocity becomes 1 − ζ 2 and ζ will be given by q E1 /2k. Expansion of this gives that ζ is of order 1/k, so that the propagation velocity is propagate with the speed of light up to first order. Note that µ V , due to its definition as a gradient, in itself is one order higher in k than the other quantities. This can also be seen from equations (67) and (68). The leading coefficient on the right hand side in (67) goes as k 2 , whereas the leading coefficient in (68) goes as k 0 . On noting that each time derivative picks out a factor of order k, we then for consistency need µ V to be one order higher in k than Σ T . Hence equation (131) for the electric part of the Weyl tensor E T is found to be sourced by the density gradient µ V to first order. 3 The shear equation (110) becomes where s Σµ is given by (112), and the evolution equation (67) for the density gradient becomeṡ on keeping terms to the two highest orders in k. Substitution of (136) in (137) gives The density gradient µ V on the right hand side can be eliminated by using (68), which to leading order is given by Substitution of (138) into (137) gives where H T and Ω V are determined from (132) and (66) respectively. Hence the shear waves propagate as damped 4 sound waves with Ω V and H T appearing as source terms with coefficients of order k, i.e. they are corrections of the same order as the damping. The last variable µ V represent density waves which also propagate at the speed of sound, but π/2 out of phase relative to the shear waves, as can be seen from (138).
In the same manner the evolution equations for the odd parity sector, equation (102) and (103), can be rewritten asË where Since there is no appearance of the vorticity in these equations the result is exactly the same as the one obtained in [29,30] where the vorticity vanishes on the perturbed model. Hence, both E T and H T represents decoupled gravitational waves, which will to first order propagate at the speed of light.

VIII. CONCLUSIONS
In this paper we have presented a general treatment of perturbations of orthogonal LRS class II cosmological backgrounds with non-vanishing vorticity. We have utilized the 1+1+2 covariant decomposition of spacetime and a harmonic decomposition of the resulting quantities, and the result is given as time evolution equations for the following eight harmonic coefficients, which decouple into an even sector and and odd sector The remaining 27 harmonic coefficients can be desribed by the eight given above, which then completely describes the full model in a gauge-invariant way, utilizing the Stewart-Walker lemma. We found that the evolution of the vorticity decouples from the other variables, as is expected, and given a linear equation of state it can be solved for analytically. The vorticity also appears as a source term for the other free variables. In the high frequency limit the free variables are found to represent graviational, density and shear waves, and the vorticity source terms disappears for all variables except the shear and density waves.
A natural extension to the work presented in this article would be to introduce vorticity creating mechanisms by considering, e.g., imperfect fluid perturbations with viscosity and heat flow. This would imply the introduction of more terms in the energy momentum tensor T ab , one symmetric and trace-free 3-tensor π ab and one spacelike vector q a , together with coefficients of dynamic and bulk viscosities and of heat conductivity. For the exact equations with a general energy momentum tensor in the 1+1+2 formalism see [10], and for the theory of relativistic imperfect fluids [50][51][52]. Since a heat flow alone cannot generate vorticity to first order in perturbation theory, [36,37], it would then be of interest to consider also second order perturbations. For works on second order perturbations in covariant perturbation theory see, e.g. [53].

Appendix A: Commutation Relations
For a zeroth-order scalar field Ψ on orthogonal and homogenous LRS class II backgrounds, with φ = 0, the following first order commutation relations hold: For a first-order 2-vector Ψ a the following commutation relations, to first order, hold: and for a first-order trace-free and symmetric 2-tensor Ψ ab it holds that:

Appendix B: Harmonics
For both the vector and tensor spherical harmonics there exists orthogonality, algebraic and differential relations that are essential when decomposing the linearized 1+1+2 equations in terms of harmonic coefficients. For the vector spherical harmonics, the orthogonality relations are the algebraic relations are and the differential relations arė For the tensor spherical harmonics the orthogonality relations are the algebraic relations are and the differential relations arė for the 2-vectors arė and for the 2-tensors arė The constraints for the scalars are for the 2-vectors are ik and for the 2-tensors are

Odd Parity
The evolution equations for the scalars are for the 2-vectors areΩ and for the 2-tensors arė ik a 1 a 2 α V = − 1 a 2

Odd parity
For the odd harmonic coefficients the following equations hold,