Induced cosmological anisotropy by a gauge-gravity interaction

We present a simple model which generates cosmological anisotropies on top of standard FLRW geometry. This is in some sense reminiscent of the mean field approximation, where the mean field cosmological model under consideration would be the standard FLRW, and the anisotropy is a small perturbative correction on top of it. Using a supergravity-inspired model, we confirm that the stable fixed point of our model corresponds to standard FLRW cosmology. We use a Bianchi VII$_0$-type model supplemented with a scalar and $U(1)$ gauge fields, and we show that the anisotropies of the geometry are generated by the non-trivial interaction between the gravity sector and the $U(1)$ gauge sector. Studying the attractor flow, we show that the anisotropies are present at early times (high redshift) and decay asymptotically to an FLRW attractor fixed point. With such a mechanism, observations of non-isotropy are not contradictory to FLRW geometry or indeed the $\Lambda$CDM model. Such models could in principle shed some insights on the present cosmological tensions.


Introduction
One of the most successful cosmological models based on General relativity is the base Lambda Cold Dark Matter (ΛCDM) model.This tremendously well established model of cosmology assumes a flat universe, cold dark matter (CDM) and a positive cosmological constant, and is the simplest cosmological model which is fairly in good agreement with current observations.As the current de-facto standard model of cosmology, the spacetime geometry in ΛCDM is that of the homogeneous and isotropic Friedmann-Lemaître-Robertson-Walker (FLRW), where the only inhomogeneities allowed are those of small perturbations, which are actually the sources of some of the most important cosmological observables.Observations of the Cosmic Microwave Background (CMB) [1], Baryon Acoustic Oscillations (BAO) [2], and Large Scale Structure [3] have for a long time been satisfactory proof that the Universe is evolving very closely along the predictions of the ΛCDM model.
Although the standard cosmological model has been a resounding success, there exist several problems which emerge when confronting the model with data.One of the most topical is the Hubble tension, the discrepancy between the value of the Hubble constant H 0 when measured in the local Universe versus with CMB observations, a tension which is currently reported at 5σ [4].This is just one of a slew of "cosmic tensions" persistent within the ΛCDM paradigm, an overview of which can be found in [5].These cosmic tensions are not the only threats to ΛCDM: different types of anomalous anisotropies have been reported both in the early and late Universe, such as quadrupole-octopole alignment in the CMB [6,7], anomalous bulk flow [8][9][10] 1 , radio-galaxy dipoles [7,11], and possible variations in the fine-structure constant [12].Also, recent hints of cosmic birefringence, the rotation of the polarisation plane of CMB photons, were reported at over 3σ in the Planck EB power spectrum [13][14][15].This is in sharp contrast to the ΛCDM prediction (no birefringence) and would have profound implications for fundamental physics if confirmed.It seems clear that the ΛCDM model may need to be revised.
In [16] the authors considered a model which is closely related to the case studied in our present work.Essentially, the model in this paper is a special case of the phenomenological model considered in [16], where the coupling between the gauge field and scalar field is taken to be minimal.The authors of [16] study an inflationary scenario with a vector field coupled to an inflaton field and show that the inflationary universe is endowed with spatial anisotropy for a wide range of coupling functions f (ϕ), where ϕ is the inflaton field; importantly, the gauge-field ansatz considered in [16] is gauge inequivalent to the gauge field considered here, and we consider the evolution of the universe without specifying the inflationary scenario.The authors in [16] focuses on the early universe evolution of the gauge fields and the inflaton fields, whereas we are more interested in the full history of the cosmological evolution after inflation. 2n this paper, we introduce an abelian version of the chromo-natural models discussed in [17][18][19][20] 3 .This type of model has been shown to arise naturally in N = 4 supergravity, and has been used to study spacetimevarying couplings as discussed in [21] and others.We begin our analysis in a very general way by using the Bianchi I spacetime before specialising to Bianchi VII 0 in Section 3, and by employing a perturbative scheme, we show that the model contains ΛCDM at the zeroth order, and that FLRW geometry is a stable point in the attractor flow.As such, there is no contradiction between the observed cosmological tensions and anisotropies and the ΛCDM model.
This paper is organized as follows: in Section 2 we introduce the model and the theoretical details; in Section 3 we discuss the covariant equations of motion and their perturbative expansions; Section 4 contains the numerical solutions, where we also present our main results; in Section 5 we present the the dark energy equation of state generated by the gauge field and anisotropies; in Section 6 we compare the model to ΛCDM using lowredshift data, and we conclude in Section 7. Appendix A contains a short treatment of the general Bianchi classification; in Appendix B we present the Killing symmetry of 1-form fields and the 2-form fluxes; Appendix C and Appendix D contains the metric gauge choice and our procedure for generating initial conditions, respectively.Finally, we present the relevant Einstein equations and the perturbative expansions in Appendix E.
We use c = ℏ = κ = 8πG = 1 and the metric signature (− + ++) throughout the paper.When studying the behaviour of the model, we focus on the time after recombination, i.e. redshift z < 1100 and thus focus on the matter and Λ dominated eras.

Gauge-Axion model
In this section, we focus on the bosonic part of a supergravity-inspired model with the action where κ = 8πG (which we set to unity from now on), R is the Ricci scalar, Λ is the cosmological constant, ϕ is the pseudoscalar axion field, Θ is the axion decay constant, and L PF is the canonical Lagrange density for a perfect fluid containing baryonic matter, dark matter, and radiation.Here, is the field-strength tensor for the gauge field A µ and F µν = 1 2 ϵ µναβ F αβ is its dual where ϵ µναβ is Levi-Civita tensor.The new field ϕ can be thought of as a candidate for axionic dark matter and/or dark energy.The gauge-axion Lagrangian considered in this work is very general, which can encompass a very general class of Bianchi models; viz, Bianchi type I.In what follows, we consider the gauge field as some dark sector component, but since it is massless, it may in principle be thought of as the photon; in this case, our solutions will be further constrained by e.g.primordial magnetic fields [22].
We note here that a stringent supergravity model would not allow us to have any explicit cosmological constant term in the action.However for the present paper where we mostly study an effective cosmological model, such constraints coming from supergravity can be relaxed and we present our action with explicit cosmological constant term.
In the rest of the paper we will mostly focus on the abelian U (1) gauge field A µ , which together with the ansatz chosen makes all contributions from the symmetry-breaking term (∝ Θ) vanish 4 .We note here that for the most general gauge field and metric ansatz, i.e full dependence on the time and the spatial coordinates, the symmetry-breaking term does not vanish and has non-trivial contributions which we defer for future study.The model that we use in the present analysis can also be considered as minimally coupled Quintessence with electromagnetic fields [23][24][25][26].In minimally coupled Quintessence models the Quintessence (scalar) field couples to the Maxwell term 5 , which is in contrast to the gauge-axion model where the pseudoscalar axion couples to the CP -violating Θ term.
It is also worthwhile to note that our analysis can be extended to nonabelian sectors, viz.SU (2) or SU (3) gauge groups [27][28][29], which, when coupled to the axion field would encode a QCD axion, which is among one of the most compelling candidates for physics beyond the standard model (BSM).This axion solves the strong CP problem [30,31] and is potentially a natural candidate for cold dark matter [32,33].In string theory, a similar spectrum of particles dubbed axion-like particles (ALPs) can be identified as ultralight dark matter with a broad mass range and interesting cosmological consequences [34][35][36].In general, the abundance of axion-like dark matter is determined by the axion mass term and the coupling of the axion to the gauge sector, i.e the decay constant, which depends on the cosmological epoch when the Peccei-Quinn (PQ) symmetry breaking takes place [37,38].In the non-abelian case, the axion (pseudo-scalar) and gauge term coupling do not vanish and is being contributed by the Chern-Simons type terms.In this case the Einstein equations will have one extra term proportional to the structure constants f a bc (which vanish in the U (1) limit).The equations of motion derived from Eq. ( 1) are given below.The Einstein equations where we add the stress-energy tensor for a perfect fluid, T PF µν .We have simplified Eq. ( 2) by including the deviation from the base ΛCDM in T AN µν , which we call the anisotropic stress-energy tensor; it takes the form Equations of motion for ϕ and A µ We choose as our starting point the Bianchi I metric, which we parametrize as where α(t) and β i (t) are the isotropic and anisotropic scale factors, respectively (for details, see Appendix A).The factor two in the exponentials has been introduced so that the isotropic scale factor matches its FLRW equivalent, i.e. a(t) = exp(α(t)), and ȧ/a = α.We also adopt the temporal gauge for the gauge fields and write In Appendix B we explicitly show that the 1-form gauge field is invariant under the Killing symmetry of the metric (6), which allows us to expand the 1-form field as follows where e i are the spatial triads, which take the following form6 , (δ i is the Kronecker delta) With the Bianchi I metric (6) with R 3 symmetry, we can write the gauge field as which allows us to rewrite the 1-form fields in terms of some scalar functions which we call ψ i (t), α(t) and β i (t).In the following section we proceed by writing the most general coupled differential equations for the metric ansatz (6) and the 1-form fields (7).In the rest of the paper we will focus only on the U (1) 1-form field strength, and it can be shown that the symmetry-breaking term proportional to Θ vanishes identically for the abelian sector.The most general solution for non-abelian 1-form field strength will be discussed in the forthcoming paper [39].

Equations of motion and their solutions
We substitute the metric (6) into the equations of motion ( 2), (4), and (5), and explicitly write out the results for each index value; after some simplification, we can write the scalar-field equation as With our gauge choice (temporal gauge), the temporal component of the gauge-field equation vanishes, and we can write the spatial components as We write out all the components of the Einstein equations (2) in a similar manner; these are somewhat lengthy, and we show the first Friedmann equation (µ = ν = 0 component) here (the rest can be found in Appendix E) In the rest of this paper we incorporate the contribution from the cosmological constant Λ into the stress-energy tensor for the perfect fluid as follows and we work only with T PF µν (without tilde) from now on.
The stress-energy tensor T PF µν for the perfect fluid is given by 7 where ρ = i ρ i is the energy density, p = i w i ρ i is the pressure, and w is the equation of state parameter, which takes the values w = −1 for the cosmological constant, w = 1/3 for radiation, w = 0 for baryonic matter, and w = −1/3 for curvature.Taking the flat (zero spatial curvature) case, the components of T PF µν for the homogeneous and isotropic (zeroth order) limit reads as follows where we have denoted the zeroth-order part with a superscript (0).The full order can be found in Appendix E.
To simplify the equations of motion, we rewrite the components of the gauge field (10) by introducing two new scalar fields, σ(t) and γ(t) and redefine the ψ i 's as which will be useful when reducing the solutions to the homogeneous and isotropic (FLRW) limit 8 .Given these redefinitions, it is easy to see that the 7 The stress-energy tensor is given by for a boosted fluid.In this paper we consider a fluid four velocity given by u µ = (1, 0, 0, 0), with the normalization u • u = −1.Note that the velocity field does not receive any corrections from the non-trivial metric evolution. 8The number of degrees of freedom is the same.
isotropic condition is The metric as written in Eq. ( 6) has had its symmetries broken down to R × R × R, which is equivalent to the Bianchi I spacetime; in order to restore which sets the components of the gauge field to A 2 = A 3 .This choice brings us to the final metric which we use in the rest of this paper as which is equivalent to Bianchi VII 0 .The symmetries of this metric encapsulates the idea that the universe has a kind of preferred direction or symmetry axis, along which the cosmic expansion evolves differently.

Perturbative Analysis
The equations of motion in Section 3 have now been reduced to a system of coupled second-order scalar differential equations.At the zeroth order, the universe evolves in a isotropic and homogeneous space-time, and the first order contribution of the gauge-field driven anisotropy is small.The intergalactic gauge field decays away rapidly [40] in the late-time evolution, and since this gauge field is driving the anisotropies, we can therefore study them perturbatively.
In order to obtain numerical solutions, we use a perturbative approach and employ the following scheme: • Expand all scalar degrees of freedom ζ = {α, β i , ϕ, ψ, σ} in a perturbative series around their equilibrium fixed points (homogeneous and isotropic fixed point) and retain only the linear order in perturbations where ϵ is a book-keeping device for perturbative order.
• Plug the zeroth-order solutions back into the equations, where they act as seed solutions for first order.
By introducing the perturbative parameter ϵ, we explicitly note that the metric anisotropies are small, but we have not linearised the new metric functions β i .Following the above scheme we write out the perturbative expansions around the homogeneous and isotropic fixed points as where we have used the remaining gauge freedom in the metric to set α (1) (t) = 0 (For details, see Appendix C).We have also set σ (0) (t) = 1 and β (0) i (t) = 0, since this represents the homogeneous and isotropic zeroth-order background; moreover, we set γ(t) = 0 to restore the planar SO(2) × R symmetry.
The perfect fluid evolves according to the continuity equation, which in the ΛCDM case reads ρ + 3H(1 + w)ρ = 0.This equation changes due to the present non-trivial Bianchi V II h geometry [41,42].The implications and perturbative corrections to the contiuity equation and T PF µν are presented in Appendix E.
From now on, expressions of order ϵ will always be enclosed in square brackets.

Zeroth order
As a first consistency check, we start with the zeroth-order vacuum equations, where we set (T PF µν = 0) and ϕ (0) (t) = ψ (0) (t) = 0, leaving us with a system of equations which is the flat-space vacuum.In this case, the system we need to solve is the two Friedmann equations, which read The allowed solution of the above equation gives a constant solution for the α (0) .With the identification of α (0) (t) = log a(t) gives the physical scale factor, and which reduces to the familiar solution for a static Universe.
Adding now a radiation term in the stress-energy tensor, the Friedmann equations read which solves as and the corresponding scale factor reads a(t) = (2H 0 Ω 0 r ) 1/2 √ t, which is consistent with standard FLRW evolution.We now turn our attention to the more general case when ϕ (0) and ψ (0) are non-zero.Here, the dynamical variables are ϕ (0) , ϕ (0) and ψ (0) , and we have the following three equations for the scalar, gauge field, and Einstein parts, respectively By examining the full set of equations in Appendix E, we notice that all terms containing σ (1) or β i , i.e. the anisotropic variables, are proportional to ψ (0) or its time derivative.This influences our choice of initial conditions in the numerical solutions: if we simply choose ψ (0) (0) = const.and ψ(0) (0) = 0, we obtain a solution proportional to a constant ψ (0) , and this can simply be gauged away.In order to obtain a meaningful solution, we therefore have to implement a non-zero ψ(0) as our initial condition.A description of our method for choosing consistent initial conditions can be found in Appendix D.

First order
The scalar equation (11) reads where the factor 2 on β (1) 2 comes from β 2 = β 3 .The expressions for the gauge field and Einstein equations are rather lengthy, and we will only display the zeroth order in this section, including the full equations in Appendix E.
For the gauge field in Eq. ( 12), the zeroth-order expressions are identical µ = 1, 2, 3, but the equations differ at first order, and due to the symmetries, the µ = 2 and µ = 3 components are equal.Keeping to our choice of a positive sign for σ (0) = +1, all the spatial components are identical (at zeroth order), and read The first Friedmann equation (µ = ν = 0 component of the Einstein equations) read T PF,(0) 00 The spatial diagonal components (µ = ν = i) are identical at zeroth order and read T PF,(0) 11 = −e 2α (0) 2α (0 We choose a simple ϕ 4 -type potential for V (ϕ) as where V 0 is a constant, and we expand V (ϕ) and its derivatives; the potential reads In order for the kinetic term to not dominate over the potential at all times, we have set the value of the constant, V 0 = 10 −3 in our numerical computation.

Numerical solutions
We solve the full system of coupled differential equations for scalar, gauge field, and Einstein parts order-by-order and present the relevant solutions here; the full equations can be found in Appendix E. When generating these solutions we fix the background FLRW cosmology to the parameter set The qualitative behaviour of these solutions indicate that the field content ϕ(t) and A µ (t) have considerable contribution in the early Universe before decaying exponentially, and eventually flowing to the homogeneous and isotropic attractor fixed point, which exactly corresponds to FLRW.The initial conditions for all the vari-Zeroth order Table 1: Boundary conditions used in the numerical solutions, defined at t f = 20 Gyr.
ables are in general coupled, and need to satisfy the equations of motion; therefore, the conditions shown in Table 1 are the ones we choose as "primary", whilst the rest are derived.In Appendix D we present our method for finding the rest of the boundary conditions from the Einstein equations in a consistent way.
From the zeroth-order equations we can solve the isotropic part of the scale factor α (0) from the zeroth-order Einstein equations.Here we have imposed boundary condition at the isotropic fixed point and solved the evolution of the Einstein equations.The evolution of the zeroth order scalar and the gauge fields, ϕ (0) and ψ (0) respectively.
The second order differential equations governing the evolution of the Einstein equations, 1-form gauge fields and the scalars are roughly damped harmonic oscillators, the solutions of which contain both growing and decaying modes; however, to be consistent with observations of the late-time universe, the evolution should settle down to homogeneous and isotropic solutions, viz.FLRW universe.In order to keep consistency with the cosmic no-hair theorem (the scalar/hairy solution should decay at late times) we have imposed the boundary condition at (t ∼ 20 Gyr); the evolution at early times is governed by the Einstein equations.
In our numerical solutions we retain the decaying solutions.

Numerical results:
• In Figure 1 we present the solution of the isotropic scale factor.Our result at current epoch, viz.t 0 = H −1 0 = 13.7 Gyr, in good agreement with the results in [43].The isotropic scale factor has been plotted against the scale factor of ΛCDM (which has been normalized to unity at the present time.The deviation from the ΛCDM value can be attributed to the scalar and gauge fields in the present model under study. Next we focus on the deceleration parameter, which for ΛCDM is canonically defined in terms of the scale factor (a(t)) as In Figure 2 we compare the deceleration parameter for the model under consideration with ΛCDM, and we notice that the present model has marginally faster expansion (q more negative), with the difference being most pronounced between t = 3 − 10 Gyr.This faster expansion is expected to play a crucial role in alleviating H 0 tension in this model.
• In Figure 3 we present the solution for the scalar fields.The scalar field profile starts with a non-zero divergent nature in the early universe, before rapidly decaying and finally saturating to zero at very late asymptotic times.This axion-like particle can be attributed to the scalar dark sector contributing to either dark energy (and/or dark matter).In the following section 5 we examine the energy equation of state, which confirms our observations here.We also show the evolution of the equation of state for the scalar field ϕ in Figure 4, which can be seen to exhibit kination behaviour for most of cosmic history, only decreasing in value slightly at very early times.
• In Figure 5 and 6 we show the behaviour of the fields ψ and σ, both of which take on very small values, even at early times, before flowing to the attractor fixed point asymptotically, which is consistent with our construction.Essentially there will be no residual gauge fields in the future and only residual gauge-field contributions would survive to the present epoch ∼ 13.7 Gyr; this is consistent with present observations.
One crucial point at this juncture is to bear in mind the overall picture: the backreaction from the U (1) gauge fields are generating the anisotropies in the early Universe, and the anisotropies settle down to their fixed-point values as the gauge field saturates to the attractor fixed points.
• The zeroth-order solutions of the Friedmann equations dictate the isotropic evolution of the universe, which is the base ΛCDM; however, we notice that there is some deviation due to the residual presence of the scalar and gauge-field contributions, where the contribution from the anisotropic parameters appear as perturbative corrections.
The anisotropic contributions to the metric, β 1 and β 2 , are suppressed by order 10 −6 as compared to the isotropic scale factor, which is in agreement with the observational constraints where the anisotropy in the universe is comparatively very small as compared to the isotropic scale factor.In Figure 7 we show the evolution of the anisotropic scale factors exp(β 1 ) and exp(β 2 ), which flow towards the stable fixed point at late times, exactly the isotropic limit (Note that β 1 should be further suppressed by ϵ), in keeping with observational results.The apparent mirror similarity in Figure 7 is a consequence of the coupled nature of the equations of motion, where we are only able to choose three out of the four initial conditions related to the β (1) i 's (as seen in Table 1), and the fourth condition is then imposed for self-consistency (as shown in Appendix D), which selects the depicted solutions for the anisotropic scale factors.We also present the total anisotropic scale factor, which is the exponential sum of the β (1) i 's, where we clearly see that it saturates to unity at late times, since the anisotropies decay; Figure 8 depicts this behaviour, clearly showing the return of homegeneity and isotropy at late times.
• In order to quantify the evolution of the anisotropic degrees of freedom, we define the average Hubble parameter H as follows In Figure 9 and 10 we show the full contribution of the anisotropy to the Hubble parameter compared to base ΛCDM.From these two plots we can see that the average Hubble parameter H is slightly smaller than its ΛCDM counterpart at all times, but that this difference is larger at early times.We also see that when compared to the isotropic limit of the present model (Figure 10), the effects of the anisotropies are on the order of ≤ 10 −7 throughout the history of the universe, though divergent as very early times 9 .
Using the isotropic Hubble parameter (Eq.( 32) for β (1) i → 0) we can construct the time-dependent energy densities for matter and Λ as Ω X (t) ≡ ρ X /ρ c , where ρ c is the critical density.Using these quantities we can establish the relative contributions of matter and Λ to the total energy budget of the Universe across cosmic history.We also form the analogue of the energy densities when taking anisotropic evolution into account the average Hubble parameter and scale factor.We plot these quantities in Figure 11 (where quantities formed with the average quantities are denoted with an overbar).Due to the attractive nature of the potential, we see a generally lower values at early (late) times for matter (Λ), which causes the deviation in the deceleration parameter seen in Figure 2.
The effects of the anisotropic variables on cosmic evolution may be important when studying the H 0 tension and other cosmological puzzles, but a detailed treatment of observational signatures consistent with the observational signatures of the H 0 tension requires some more exploration and lies beyond the scope of this paper, although we give some brief comments below.
We end this section with some plausible implications of our axion-anisotropic cosmological model on the resolution of the present cosmological tensions.A naive observation from the solution of the average Hubble parameter from Figure 9 indicates that the value of Hubble parameter is lower than in the base ΛCDM model, especially at very early times.A natural question to ask at this juncture is: Can the Hubble tension be resolved in the presence of some extra degrees of freedom on top of standard FLRW cosmology?Let us briefly present the possibility of the model under consideration in resolving one the specific cosmological tension; viz, the H 0 tension.For an efficient resolution of the H 0 tension in the context of any effective-field theory approach, the predicted Hubble parameter should be large (∼ 73 km s −1 Mpc −1 ) compared to the standard prediction from the astrophysical models of ΛCDM.A quick comparison of the Hubble parameter for the model under consideration and that of ΛCDM in Figure 10 indicates that the Hubble parameter of ΛCDM should be higher; a naive conclusion would be that the model presented in this paper not efficient in resolving the H 0 tension effectively.Some plausible explanations for this include • In Figure 12, we observe that the dominant contribution to the dark energy induced by the anisotropic matter sector is controlled by the scalar fields; effectively, the kinetic terms of the scalar fields are dominant (which is why the energy equation of state saturates to unity) and the contribution of the gauge fields are negligible.This gives a possible explanation: as the Universe starts to expand under the gravitational force, the scalar fields tries to counterbalance the expansion; thus, there is small dip in the Hubble parameter compared to ΛCDM.
• This is in general true for any EFT which has dominant contribution from the bosonic sectors.
In [27] the authors showed that a rolling axion coupled to a non-Abelian gauge field has the potential to provide a viable solution to the Hubble tension.The pertinent point made in [27] is that the axion fields coupled to non-abelian gauge fields provides some additional friction term (thermal friction) to the gravity system, and thus have a potential solution to stabilize the Hubble tension.

Anisotropic dark energy
From our construction it is worthwhile to investigate the anisotropic contribution to the energy equation of state.We can write the anisotropic stress-energy tensor (3) in the standard form as In the particular case of homogeneous and isotropic cosmological models, we can assume an equation of state of the form Figure 1: The isotropic scale factor a(t) = e α 0 (t) compared with the ΛCDM model.
and in the presence of anisotropic matter sources and geometry, the total pressure and the total energy density can similarly be split into isotropic and anisotropic parts from which we can determine the effective equation of state parameter w t for the cosmic fluid, as was also noted in [41,42,[44][45][46].Note that we show in Appendix E that the perfect-fluid part also receives corrections at order ϵ; these contributions are coupled to the anisotropic degrees of freedom, and we count them as part of ρ AN 1 and P AN i .In Figure 13 we show the evolution of w t as a function of time, and we observe that it stays negative throughout all of cosmic history, and is close to, but always lower than, the ΛCDM model.From the point of view of the perfect fluid, the negative values of the equation of state parameter are to be expected, since we neglect the radiation term w r = 1/3, and ω ≤ 0 for both matter and cosmological constant.
It is also interesting to examine the contribution to w t from the anisotropic variables.First of all, by examining the anisotropic energy density ρ AN 1 in Eq. (35) and comparing it to the perfect fluid, we see that ρ PF dominates, and the anisotropic parts make up on the order of 10 −5 of the total energy budget of the system.Moreover, when examining the equation of state for the anisotropic contribution (which we may call w AN ), we see that up to a This corresponds to a stiff matter fluid, which has been studied in the context of both classical and quantum cosmology in [47] and [48] and others.Specifically, it was found in [47] that a stiff fluid may lead to a bouncing solution of the Wheeler-de-Witt equation.

Exploring the parameter space
Given the predictions our model makes, it is interesting to compare it to some available data.In this Section, we perform a post-fit analysis10 using late-time cosmological data at the background level (using only distance measures).In this paper, we have considered the case of vanishing radiation density (Ω 0 r = 0), and our model should therefore provide its best fit at low redshift; late-time data should therefore be sufficient to gain some insight of the overall fit of the model.To accomplish this, we employ a combination of two robust local-Universe datasets as described below.
We use the Pantheon+ catalogue of Type Ia supernovae (SNeIa) with SH0ES Cepheid host calibrators [49,50], which is a set of 1701 light curves Figure 3: The behaviour of the full scalar field ϕ(t) = ϕ (0) + ϵϕ (1) .and 1550 resolved SNeIa in the redshift range 0.001 < z < 2.26.The inclusion of Cepheids with known distances provides a robust calibration of the SneIa lightcurves and breaks the degeneracy between H 0 and Ω 0 m .In order to compare the model with this dataset, we construct the theoretical distance modulus as µ th (z hel , θ) = 25 + 5 log(d L (z hel , θ)), (36) where z hel is the redshift in the heliocentric frame, θ is a vector containing the model parameters, and d L is the luminosity distance (for full definitions of the distance measures, see for example Appendix A of [51]).On the data side, the observed distance modulus reads µ data = m − M , where m is the standardised apparent magnitude in the blue band, and M is a fiducial absolute magnitude calibrated using the Cepheid host distances.In order to compare the model with ΛCDM, we also compute the χ 2 values for our model.We form the measure ∆µ depending on whether the SNeIa data points has an associated Cepheid host as with the corresponding χ 2 measure being where ∆D = D theory − D data , and C tot is a covariance matrix containing statistical and systematic uncertainties for both the SNeIa and Cepheids.We also include measurements of the Hubble parameter from Passivelyevolving Early-Type Galaxies (ETG), which have an old stellar population and thus a low star-formation rate.It is possible to reliably trace the spectral properties of ETG's along cosmic time (independent of the cosmological model), making ETG's a standardisable clock (they are also known as Cosmic Chronometers (CC)).For this purpose, we use a sample in the range 0 < z < 1.97 [52][53][54].In order to construct the χ 2 for the CC's, we follow the same prescription as above and write where ∆H = H theory − H data , and C CC is a covariance matric containing statistical, sample-contamination, model dependence, and stellar metallicity uncertainties 11 .We investigate the fit of our model to these two data sets using three sets of parameter values.Since we are using late-time data and are considering a flat Universe, the free parameters are {H 0 , Ω 0 m }.Since we are not performing a Bayesian likelihood analysis at this stage, we will fix the fiducial absolute magnitude M = −19.5, which is close to the canonical value.In this analysis, The behaviour of ψ(t) = ψ (0) + ϵψ (1) .
we use the average Hubble parameter (32) in the definition of the distance measures.Here, we are varying only the standard cosmological parameters and do not consider the contribution to the energy densities of the scalar and gauge field, which are fixed by initial condition.As such, we are likely overestimating the value of Ω 0 m which in principle obtains contributions from the new scalar (depending on the equation of state), but this approach is enough to give an indication of the overall fit at late times.Since our numerical results in the previous sections indicate that the anisotropic effects are small at late times (low redshift), we pick the parameter values to lie close to, but slightly deviating from the ΛCDM values.As such, we are able to estimate the deviations induced by our model as well as its sensitivity to the parameter values. 12Figures 14 show the Pantheon+SH0ES data as a function of redshift (where h = H 0 /(100kms −1 Mpc −1 ) is the dimensionless Hubble parameter).We observe that in general, raising the Hubble parameter from the CMB value (h = 0.68) to the local Universe value (h = 0.73) provides a better fit to the data.The same trend can be seen in Figure 15 showing the data from Cosmic Chronometers, where the combination h = 0.73, Ω 0 m = 0.27 provides the best fit to the data points.This also has implications for the age of the Universe as chosen in our analysis.The behaviour of σ(t) = ϵσ (1) .
In order to compare with the ΛCDM model, we compute the χ 2 statistic as described above, and we show the values in Table 2.We find that ΛCDM, (which for the present parameter values gives a χ 2 ∼ 10 3 , the minimum χ 2 being 1600) fits the Pantheon+ SH0ES data in Figure 14 significantly better than our model (where we find χ 2 ∼ 10 5 ), the difference in χ 2 being around 10 5 .We also note that the best-fit parameter set (out of those considered) is not the same for our model as for ΛCDM ; in fact, the best fit to our model (χ 2 = 2.57 • 10 5 for {h = 0.73, Ω 0 m = 0.27}) is the worst fit for ΛCDM (χ 2 = 3.68 • 10 3 ).The conclusion of this simple comparison is that even though ΛCDM provides a better fit overall, the best-fit parameters are likely to be different in our model.
For the CC data, the difference is much smaller, and for one parameter combination, the χ 2 differ by a factor of 4 as compared to ΛCDM .Table 2 shows the χ 2 values for the set of parameters chosen.For this dataset, our model shows very low sensitivity to the parameter set chosen, the best and worst χ 2 differing only by 460.52 − 460.41 = 0.11.In contrast, ΛCDM shows greater variability, with the χ 2 ranging from 14.84 (best), to 187.29 (worst).Interestingly, the parameter set giving the lowest χ 2 (h = 0.7, Ω 0 m = 0.3) is the same for both models, which is in contrast to the Pantheon+ SH0ES case described above.Overall, the fit to Pantheon+SH0ES data is worse for all parameter choices (not including the naturally higher χ 2 due to the number of data points being higher), although it is possible that our model prefers  a different value of M .A complete Bayesian inference analysis of the model (through which we will be able to place error bars and significance levels on the model parameters) is challenging, and will be presented in a forthcoming paper.

Discussion & Conclusions
In this paper we introduce an axion-electrodynamics model for the purpose of generation of cosmological anisotropies.Working with abelian gauge fields, we choose the components of the gauge field A µ to be aligned with the Killing vectors of the Bianchi VII 0 metric, and we show that the field content satisfies the same isometries as Bianchi VII 0 .We solve the resulting equations of motion numerically using a perturbative scheme where the zeroth order is the homogeneous isotropic limit; in this way, we obtain the canonical ΛCDM solutions at zeroth order, with anisotropic contributions appearing at first order.Thanks to the parametrisation of the gauge field, we obtain solutions to the anisotropic scale factors β (1) i which are driven by the evolution of the gauge-field A µ , and by constructing the average Hubble parameter H, we see that the deviation from ΛCDM is largest in the early universe, before relaxing down to the asymptotic ΛCDM fixed point.The magnitude of H is always smaller than H ΛCDM , and a negative slope at all times, which may have implications for the Hubble tension.Simultaneously, the isotropic scale factor exhibits approximately standard ΛCDM evolution throughout the history of the Universe, although the amplitude is consistently higher.Our solutions for the anisotropic scale factors exp(β 1 ) and exp(β 2 ) are very similar in amplitude, but not identical; this is a desirable feature, since cosmological anisotropies are expected to be small, and by evaluating exp(β (1) 1 ) and exp(β 2 ) at the present time (t 0 = 1/H 0 ), we find that the anisotropic expansion is on the order of 10 −7 − 10 −8 ; by examining H in Figure 9, we see that a large part of the anisotropies have decayed away at t = 5 Gyr.The scalar field ϕ exhibits steep falloff in the early Universe and settles down to a small constant at late times, and we find similar behaviour in ψ and σ, which parametrize the gauge field.A related model was studied in [16] and similar results were found, but as discussed in the Introduction, this is gauge-inequivalent to our model.
Taken together, these results indicate that most non-trivial effects will be contained to the early universe.Whilst this does safeguard late-time evolution against large anisotropic effects, this is not necessarily desirable, since early-Universe processes (inflation, BBN, recombination etc) are very sensitive to the field content and initial conditions; in particular, early-Universe observables such as the sound horizon may be modified in the presence of anisotropies, in an analogous way to that of early dark energy [55].However, this lies beyond the scope of the present work.For studies regarding anisotropies in the inflationary era, see for example [16,[56][57][58][59].
In Appendix E we find that the perfect-fluid part of the total stressenergy tensor receives anisotropic corrections perturbatively, both in the energy density and in the pressure.We also find off-diagonal components to the stress-energy tensor, which act as constraint equations, as was also studied in [60].The anisotropic part of the energy density has been studied as anisotropic dark energy, for example in [44] and [45], although at the background level.There are also interesting connections to the quadrupole anomaly in the CMB [61].
The most important result of this work is the generation of cosmological anisotropies; we have shown that it is possible to find solutions which closely resemble those of ΛCDM at zeroth order, whilst containing a small degree of anisotropic correction at order ϵ.An important note is that we are likely overestimating the magnitude of the dark-energy density Ω Λ : since the extra field content {ϕ(t), ψ(t), σ(t), β 1 (t), β 2 (t) can be interpreted as dynamical dark energy, the total dark-energy density should read Ω DE = Ω Λ + Ω ϕ + . .., but because of the small scales of the anisotropies and the field ϕ(t), this would be a very small correction 13 .
The observational status of cosmological anisotropy is rapidly evolving, with some groups claiming very strong results, such as anisotropic acceleration (anomalous bulk flow) in the direction of the CMB dipole at 3.9σ significance [63] and a 3σ hemispherical power asymmetry in the Hubble constant, also aligned with the CMB dipole14 [65].Together with probes such as fine structure-constant variation and preferred directions in the CMB results in compelling evidence that the cosmological standard model needs revision, and we have provided a mechanism through which such preferred directions can arise dynamically from a well-motivated field theory.This is of course not the only model which can generate cosmological anisotropies; in particular, models exhibiting spacetime-symmetry breaking are known to contain preferred directions in the form of timelike vector fields.For example Hořava-Lifshitz gravity [66] Einstein-Aether theory [67], and bumblebee gravity [68], all of which have received significant attention in recent years, contain preferred frames of reference.On the other hand, spacetime-symmetry breaking in gravity has been tightly constrained using the Standard-Model Extension effective field theory, restricting the available parameter space for all spacetime-symmetry breaking models [69].Our construction has the advantage of keeping these well-tested spacetime symmetries intact, and instead postulating the existence of the fields ψ(t) and A µ (t), and in this sense, it can be considered a scalar-vector model.In this paper, we presented an estimate of the fit of the average Hubble parameter of our model with late-time cosmological data for certain parameter values and contrasted it with the fit using ΛCDM.Overall, the standard cosmological model is a better fit to late-time probes.
It is worthwhile to mention [70], which have partial overlap, and of course are compatible with some of the results and statements presented in this paper.However it is important to note that the authors of [70] considered so called flowing dark-energy cosmology, where the tilt parameter is non-zero at late times (for details see Section V of [70]).This is in contrast to the model in the present paper, where all anisotropies decay to a homogeneous and isotropic fixed point, in keeping with the cosmic no-hair theorem [71].Also, even though FLRW is stable in our setup, the possibility of a tilt instability in the FLRW geometry which could potentially evade detection through the cosmic no-hair theorem was raised in [72].
A natural extensions and applications of this work would be to consider an SU (2) gauge field, as was done in the context of cosmic birefringence in [18], as well as computing imprints of anisotropy on the CMB, by introducing angular dependence of the metric functions.All of these applications, as well as parameter constraints on the present model by means of cosmological data and a Markov-Chain Monte Carlo (MCMC) algorithm are forthcoming  where α(t) and β i (t) are the isotropic and anisotropic scale factors, respectively.The factor two has been introduced so that the isotropic scale factor matches its FLRW equivalent, i.e. a(t) = exp(α(t)), and ȧ/a = α.
In this Appendix we have presented a general treatment of the Bianchi type, but in the rest of the paper (and Appendix B) we specify our geometry to Bianchi type VII 0 , which is suitable for our purposes and keeps the equations tractable.In a nutshell, Bianchi-type geometries are classified by their Killing vector fields, and in our present work we only consider Bianchi Type VII 0 .For an exhaustive Bianchi classification, we refer the reader to [76].We have three Killing vectors associated with the (B.1), The Killing vectors satisfies the following condition . Here we will use a convenient notation for ω i ∧ ω j just by ω i ω j with the property that ω i ω j = −ω j ω i .
Let us write the 2-form fluxes in the most general form as F = (f 1 dtdx 1 + f 2 dtdx 2 + f 3 dtdx 3 ) + (g 1 dx 1 dx 2 + g 2 dx 2 dx 3 + g 3 dx 3 dx 1 ) , (B.2) where the f i 's and g i 's can be arbitrary functions of (t, x i ).The equation of motion of the 2-form fields are given by dF = 0,   Where we have defined the volume form as Simultaneously satisfying (B.3) and (B.5) gives us the following relations, which implies that f i and g i are functions of only time t.
A similar analysis can be done for the 1-form gauge fields A i .The 1-form gauge field can be written as where again a 0 and b i can be arbitrary functions of (t, x i ).The fluxes can be computed as F = dA and the equation of motion is trivially satisfied, dF = 0.
Lets us write the Killing equation for the 1-form A and the algebra can be easily worked out Note here that once we consider more general metric ansatz, which depends on angular direction, we have different sets of Killing vectors, and the 1−form field strengths can depend on these variable.We defer this analysis for our forthcoming work [39].

Figure 2 :
Figure2: The deceleration parameter q compared with that of ΛCDM.

Figure 4 :
Figure 4: The behaviour of the equation of state for the scalar field ϕ.

Figure 8 :
Figure 8: The total anisotropic scale factor.

Figure 9 :
Figure 9: The average Hubble parameter in contrast to the pure ΛCDM case represented by H ΛCDM , which is the isotropic part of the Hubble parameter.

Figure 10 :
Figure 10: The average scale factor normalized by the isotropic case.

Figure 11 :
Figure 11:  The evolution of the fully isotropic energy density for matter and Λ (without bar) compared to the corresponding quantities constructed using the average scale factor and Hubble parameter (with bar).

Figure 12 :
Figure 12: The dimensionless number of anisotropic equation of state.

Figure 13 :Figure 14 :
Figure 13: The behaviour of the total equation of state parameter w t compared to that of ΛCDM

Figure 15 :
Figure 15: Hubble parameter measurements from Cosmic Chronometers (orange) with 1σ error bars [52, 53], together with the theoretical prediction from our model for specific choices of parameter values.