On the stability and causality of scalar-vector theories

Various extensions of standard inflationary models have been proposed recently by adding vector fields. Because they are generally motivated by large-scale anomalies, and the possibility of statistical anisotropy of primordial fluctuations, such models require to introduce non-standard couplings between vector fields on the one hand, and either gravity or scalar fields on the other hand. In this article, we study models involving a vector field coupled to a scalar field. We derive restrictive necessary conditions for these models to be both stable (Hamiltonian bounded by below) and causal (hyperbolic equations of motion).


Introduction
Inflationary models including a vector-field sector have been studied following diverse approaches over the past few years. Among the most recent models, a sizable fraction was motivated mainly by the appearance of certain "anomalies" pointed out by observations [1][2][3][4][5][6], later confirmed by Planck 's results [7], which suggested the presence of statistical anisotropies and maybe signals of parity violation in the cosmic microwave background (CMB). Although the statistical significance of these anomalies is still a matter of debate (possible systematic errors, contamination by foregrounds, asymmetric beams, etc.), it is not excluded that they are actually seeded by a source, different than an inflaton scalar field, during the early stages of the Universe. In this context, vector fields arise as suitable, simple and natural candidates to explain the origin of such anomalies as they possess intrinsically a preferred direction. With these motivations and the requirement to generate both the inflationary dynamics and the presence of a detectable level of statistical anisotropy in the CMB within an unified framework, several models involving vector fields have recently been proposed. Their classical dynamics and statistical properties have thus been explored in great details  (for reviews see Refs. [45][46][47]). The determination of cosmological parameters related to the presence of a statistical anisotropy in the CMB can provide valuable information about the mechanisms governing the dynamics of the inflationary Universe, and their possible deviations from the reference single-field model.
In the recent literature, popular models propose to couple scalar and vector fields by modifying the standard kinetic term of the vector as f (φ)F µν F µν [12,20]; or adding a term of the form φF µνF µν (F being the Hodge dual of F ) which couples vectors and "pseudo scalars" or axions [48][49][50][51][52] (see the review [53] for further references); or variants and generalizations of the above ideas including non-Abelian gauge fields [21,22,39,44,45]. These models have been proved to be free from instabilities, in particular they do not possess any longitudinal propagating mode; they also have the virtue of generating a non-diluting amount of statistical anisotropy [20,33] which could leave measurable imprints in the CMB.
Of course, scalar-vector theories offer a very broad set of possibilities, among which the examples mentioned above are somehow the simplest representatives. Apart from Occam's razor, there is a priori no reason to focus on these models specifically, hence one could wonder which subset of all possibilities are worth investigating. This motivates the present article, whose purpose is to identify a class of fundamentally healthy scalar-vector theories, which could then be safely considered candidates for inflationary models, or dark energy, etc. More precisely, we focus on theories involving one scalar field and one (gauge-invariant) vector field, both minimally coupled to spacetime geometry, and we study the necessary conditions for such theories to be stable-their Hamiltonian must be bounded by below-and causaltheir dynamics must be governed by hyperbolic equations of motion. This method has the advantage of being nonperturbative, and thus more general than only studying the behavior of small perturbations about a given background. The application of the healthy models to cosmology, together with the tests of their compatibility with current observations, are beyond the scope of our analysis and left as future projects.
The article is organized as follows. In Sec. 2, we derive the most general form of a theory involving a scalar field and a gauge-invariant vector field, both minimally coupled to gravity, and propose a reasonable restriction motivated by previous works. The Hamiltonian stability of this theory is analyzed in Sec. 3, and then its causality in Sec. 4. Finally, Sec. 5 is dedicated to a summary of the results, followed by a discussion about possible extensions.
2 Building general scalar-vector models

General assumptions
Consider, as a starting point, the most conservative theory in which matter is described by the standard model of particle physics, while being minimally coupled to spacetime geometry governed by general relativity. To this theory, we add two new fields, namely a scalar φ and a vector A µ , which for convenience will be referred to as "dark sector", though it can stand for inflationary models as well as for dark energy models (see, e.g., Ref. [54] for a discussion of the different classes of universality of extensions).
The presence of such new fields potentially offers a huge amount of possibilities, depending on how they couple to standard matter, to spacetime geometry, or simply to each other. Among them, many shall lead to unhealthy theories, typically due to instabilities (e.g., ghosts), or violations of causality. Since obviously we cannot analyse all possible theories, we choose to focus on those satisfying the following three conditions. 1. The fields φ, A µ are uncoupled to standard matter, and minimally coupled to gravity.
In other words, the action of the theory takes the form where S EH is the Einstein-Hilbert action, S SM the action of the standard model of particle physics, and S DS the action of the dark sector. This assumption ensures (a) the non-violation of the equivalence principle, and (b) the constancy of fundamental constants [55,56]. Note that for a scalar field alone, non-minimal couplings to spacetime geometry have been actively studied, e.g. in the context of scalar-tensor theories, and now well understood [57,58]. For a vector field alone, it has been proved that nonminimal coupling generically leads to instabilities [15,19].
2. The action only contains at most order-one derivatives of φ, A µ . This is a sufficient condition to have second-order equations of motion, though not necessary-see for instance Horndeski and Galileon models [59][60][61][62][63].
3. The action is gauge invariant 1 . This restriction is essentially chosen for simplicity. The variety of models breaking gauge invariance is indeed extremely broad, even in the absence of scalar fields (see, e.g., Refs. [64,65]), which would make the analysis performed in the present article much more involved.
The last two asumptions imply that the action of the dark sector reads where g is the determinant of spacetime's metric, and F ≡ dA = (F µν /2) dx µ ∧ dx ν , with F µν = ∂ µ A ν − ∂ ν A µ , is the field-strength two-form associated with the vector field. The latter can only appear through F in L DS , since any other type of term (e.g., A µ A µ or ∂ µ A µ ) would be gauge dependent, which is excluded by assumption 3.

A general class of Lagrangians
Let us construct the most general Lagrangian density for the dark sector, under the assumptions formulated above. As a scalar, L DS can only depend on the scalars that can be constructed from φ, ∂ µ φ, F µν ; in principle, their free indices could be contracted with arbitary tensors-standing for parameters of the theory-and lead to terms of the form However, from a tensorial parameter there generally emerges fundamentally preferred directions in spacetime 2 , that we do not wish in the theories considered here. The only nondynamical tensor escaping from this rule is the Levi-Civita tensor ε µνρσ ≡ − √ −g [µνρσ], where [µνρσ] stands for the permutation symbol, with the convention [0123] = 1. It turns out that any scalar constructed from φ, ∂ µ φ, F µν , g µν , and ε µνρσ can be written as a function of φ, whereF µν ≡ ε µνρσ F ρσ /2 are the components of the Hodge dual F of F ; so that Let us prove this assertion. First, it is clear that the Levi-Civita tensor cannot be involved without being contracted with F µν ; if both indices of the latter are contracted with the former, then it leads toF µν ; if only one index is contracted, then we get So when the Levi-Civita tensor appears, the associated expression can be rewritten in terms of F , whence L DS (φ, ∂ µ φ, F µν ,F µν ). There are two elementary classes of scalars that can be constructed from contractions of ∂ µ φ, F µν ,F µν , namely where Ù F stands either for F or forF . Indeed, since ∂ µ φ has only one index, it always ends a contraction branch, hence if more than two ∂ µ φ are involved in a scalar term, then it can be factorized into chains of the form (2.12). Finally, such Ù F -chains can in general be reduced thanks to the identities 3 14) Indeed, if in a Ù F -chain, an F and anF are contiguous, then we can use Eq. (2.15) to factorize the couple. If there are only F s (or onlyF s) in a chain with strictly more than two Ù F s, then we use Eq. (2.14) to create FF pairs, and so on. The only irreducible chain 4 is therefore F µα F να (or alternativelyF µαF να ), that is, if contracted with the gradient of the scalar field, ∂ µ φF µα F να ∂ ν φ (or alternatively ∂ µ φF µαF να ∂ ν φ = −Z). In this article, we consider Z instead of the untilded contraction, because it will turn out to be more convenient for presenting the results of Sec. 3.
3 These identities can be considered a special case of the following lemma: in a four-dimensional manifold, for any two 2- This can be easily derived by using the contraction of two Levi-Civita tensors ρ . 4 The situation changes, however, in the case of non-Abelian gauge fields, since there appear non-zero terms of the form F F F , F FF , etc. This will be briefly discussed in Sec. 5.

A reasonable restriction
In Ref. [19], the authors have analyzed the stability and causality conditions for vector theories whose Lagrangian density is an arbitrary function of F 2 and FF , i.e. L vec (X, Y ). Although general conclusions could not be drawn, it appeared that nonlinear functions of only X, or only Y , are excluded. This motivates our fourth restrictive assumption: we only consider models which are at most linear in X, Y , and Z, i.e., at most quadratic in the vector field. Thus, in the remainder of this article, we consider a dark-sector Lagrangian density of the form and investigate under which conditions on the four functions f 0,1,2,3 a model is both stable-Hamiltonian bounded by below-and causal-hyperbolic equations of motion. Our analysis can also be considered a starting point for more ambitious ones, where some of our four restrictive assumptions would be dropped (see, e.g., appendix B for elements about Lagrangian densities L DS which are nonlinear in X, Y , Z).

Stability of the models
In this section, we turn to the study of the stability of a dark sector defined by the Lagrangian density (2.16). After having computed the associated Hamiltonian density (Subsec. 3.1), we investigate in details the conditions under which it is bounded by below (Subsec. 3.2), that is necessary for the stability of the quantum theory, and we summarize the results in Subsec. 3.3. In this last subsection, we also discuss why all the results, though derived in Minkowski spacetime, are also completely valid in the presence of gravity.

Canonical momenta
The canonical momentum conjugate to the scalar field φ is where an overdot stands for a time derivativeφ ≡ ∂ t φ; as usual, greek indices run from 0 to 3, while latin indices run from 1 to 3. The canonical momentum π φ can be expressed in terms of the electric and magnetic parts E, B of the field strength tensor, defined by (we use bold fonts for spatial vectors) as where ∇ ≡ (∂ i ) is the spatial gradient, an in-line dot and a cross respectively denote the Euclidean scalar product U · V ≡ δ ij U i V j and vector product ( The canonical momenta conjugate to the vector field A µ are When expressed in terms of the electric and magnetic fields, the latter reads (3.6)

Constraint on the nondynamical field component
Since the total Lagrangian density L does not involve anyȦ 0 term, A 0 is a nondynamical degree of freedom. The associated (Euler-Lagrange) equation of motion, is therefore a constraint. Equation (3.7) can be rewritten using ∂L/∂Ȧ 0 = 0 and that ∂ i A 0 only appears within terms of the form F i0 , thus it always comes with −Ȧ i . As a consequence 8) and the constraint reads Had we considered terms breaking the gauge invariance in the action, then this constraint would have been altered on its right hand side.

Hamiltonian density
Since the dark sector is decoupled from the other fields, its contribution to the Hamiltonian density is obtained by The canonical term π iȦ i can be rewritten in the following way: 11) and the spatial divergence ∂ i A 0 π i can be dropped, since it would disappear in a boundary term while integrating H DS to build the Hamiltonian. Using the constraint (3.9) then yields Finally, we inject the expression of the canonical momenta and of the Lagrangian density, and reorganize the various terms to obtain (3.14) In Eq. (3.13), the Hamiltonian density is not expressed in terms of its natural variables, which are π φ , ∇φ, π, and ∇A µ . Here, we actually choose to describe a physical state of the theory using the time derivatives of the fields instead of the canonical momenta. This is perfectly licit, since there exists a one-to-one and onto relation between both descriptions, and this choice will turn out to make the discussions of the following sections easier.

Hamiltonian stability of the theory
In this subsection, we study the necessary conditions on the functions f 0,1,2,3 for the Hamiltonian density (3.13) to be bounded by below. Our method relies on proofs by contradiction: given some properties of f 0,1,2,3 , we look for configurations of the fields φ, A µ , such that H DS can be made arbitrarily negative. If at least one such state can be exhibited, then the theory is unstable, thus forbidden.

Conditions on f 0
For this paragraph only, and without loss of generality, we consider states for which E = B = 0, so that the contributions H 1,2,3 of the Hamiltonian density do not enter into the discussion. There are two necessary conditions on f 0 for H DS to be bounded by below, namely: ∂f 0 /∂K ≥ 0. If there existed a state (φ, K) of the scalar field so that ∂f 0 /∂K < 0, then we could takeφ, |∇φ| → ∞ while keeping K constant, which would make the Hamiltonian density diverge towards −∞. Such a situation is therefore excluded.
f 0 (φ, K ≥ 0) must be bounded by below. If there existed a positive value of K so that φ → f 0 (φ, K) was not bounded by below, then we could set the derivatives of φ so thaṫ φ = 0, hence H DS = f 0 /2, which would not be bounded by below.
Note that the above reasoning does not apply for negative values of K, since the terṁ φ 2 ∂f 0 /∂K can possibly compensate the divergence of f 0 . As an example, f 0 (φ, K) = φ 2 K is clearly not bounded by below for K < 0, but its contribution in the Hamiltonian density is hence completely admissible.

Conditions on f 3
There are two conditions on f 3 for the Hamiltonian density to be bounded by below: If there existed a configuration (φ, K) so that this derivative was not zero, then we could always tune E, B so that Therefore, we consider f 3 (φ, K) = f 3 (φ) ≥ 0 from now on.

Conditions on f 1
The conditions on f 1 turn out to be the same as those on f 3 , although their proofs are slightly subtler due to the difficulty of controlling the possible compensations between terms.
If there existed a configuration (φ, K) so that ∂f 1 /∂K < 0, then we could choose a state where E, B, ∇φ are all orthogonal to each other, and so that Therefore, we consider f 1 (φ, K) = f 1 (φ) ≥ 0 from now on.

Condition on f 2
Just as f 1,3 , f 2 cannot depend on K for the theory to be stable. Indeed, if there existed a configuration (φ, K) so that ∂f 2 /∂K = 0, then we could consider a state for which E, B, ∇φ are aligned, with sgn(E · B) = sgn(∂ K f 2 ), and set for instance In this situation, the Hamiltonian density would become so thatφ 2 , B 2 → ∞, while keeping K constant, would imply H DS → −∞. Hence, we can consider f 2 (φ, K) = f 2 (φ) from now on. Note that, contrary to f 1,3 , there is no restriction on the sign of f 2 , since the function itself does not appear in the Hamiltonian.

Summary and discussion
We have proved that, among the various couplings between the scalar field and the vector fields, many bring uncompensated instabilities in the theory, by making the Hamiltonian unbounded by below. In the framework chosen in this article, the most general Lagrangian density for the dark sector leading to a stable theory is where f 1 , f 3 are positive functions, φ → f 0 (φ, K ≥ 0) is bounded by below, and ∂f 0 /∂K ≥ 0. So far, our analysis has been performed on a Minkowski spacetime. Nevertheless, our conclusions remain valid in the presence of gravity, thanks to the equivalence principle. Indeed, the divergences underlined in the previous paragraphs are local properties, namely, they regard the Hamiltonian density rather than the Hamiltonian itself. Suppose one wishes to perform the same study in an arbitrary spacetime. Then, in the vicinity of any event E, one is free to work in a free-falling frame, where spacetime is locally Minkowskian, and thus where the above calculations are valid (modulo negligible gravitational tidal effects). In other words, in the vicinity of E, one could construct a configuration of the fields so that the Hamiltonian density is arbitrarily negative. Note that this reasoning would not be true if the fields were non-minimally coupled to gravity, or more generally in any scenario where the equivalence principle is not respected.

Causality of the models
A field theory is considered causal if it admits an unambiguous notion of time evolution; any initial condition of the fields-i.e., their state on a spacelike hypersurface-must generate a unique final state through the equations of motion. In other words, time evolution must be a well-posed Cauchy problem. This is equivalent to the mathematical statement that the equations of motion must be hyperbolic, that is, whose second-order part involve a differential operator G µν ∂ µ ∂ ν with signature (−, +, +, +), where the (−)-direction is timelike with respect to spacetime's metric.

Equations of motion
Let us first determine the equations of motion induced by the action (3.24). On the one hand, stationarity with respect to variations of the scalar field implies where primes in f 1,2,3 denote derivatives with respect to φ. On the other hand, stationarity of the action with respect to variations of the vector field implies Equations (4.1) and (4.2) form a coupled system of second-order differential equations for φ and A µ , which can be formally written as where ≡ ∂ µ ∂ µ denotes the d'Alembertian.

Diagonalizing the system of equations of motion
As it appears clearly in the expression (4.5) of D, the presence of f 3 couples the equations of motion of the scalar and vector fields even in their second-order part. As a consequence, we cannot investigate their hyperbolicity independently from each other; instead, we must consider the whole system (4.3), diagonalize 5 it, and study the hyperbolicity of each "eigenequation". In practice, we proceed by diagonalizing the principal symbol σ D (p µ ) of the system, defined as the matrix-valued polynomial obtained from the principal differential operator D by replacing ∂ µ with an abstract variable p µ . In the expression of σ D , there naturally appear three vectors, namely p µ , ∂ µ φ, and B µ ≡ p αF αµ , from which we can construct an orthonormal tetrad (e a ) a=1...4 ; assuming that p µ is not null-like, we define indeed e µ 4 ≡ ε αβγµ e 1α e 2β e 3γ . (4.10) the orthogonality between e µ 1 and e µ 3 being ensured by the antisymmetry ofF µν . Let us use these notations to rewrite the various contractions involved in the symbol σ D , (4.14) Since the above expressions exhibit projections over the tetrad vectors (such as e 1σ e σ 1 ), we expect the symbol to be much simpler if it is written in the tetrad basis 6 (e 4 , e 1 , e 2 , e 3 ) instead of the coordinate basis (∂ µ ) µ=0...3 ; and indeed the result is The five eigenvalues of σ D are therefore λ 1 = 0, λ 2 = f 1 p 2 , λ 3 = f 1 p 2 + f 3 (∂ ⊥ φ) 2 p 2 , and the two solutions (λ 0 , λ 4 ) of the second-degree equation where we used that e 2 1 e 2 2 e 2 3 e 2 4 = −1, since an orthonormal tetrad has only one timelike vector. If p µ is null, the construction of the tetrad is slightly different. One can set, e.g., e µ 1 = p µ , and e µ 2 = ∂ µ φ/(p ν ∂ ν φ) − (∂φ) 2 p µ /2(p ν ∂ ν φ) 2 , so that both e µ 1 and e µ 2 are null vectors, and e µ 1 e 2µ = 1. The other two ones, e µ 3 and e µ 4 are defined similarly as before, except that B µ ⊥ must be B µ − (e ν 2 B ν )p µ . One can check that the expression of the symbol in the basis (e 4 , e 1 , e 2 , e 3 ) is then the same as in Eq. (4.15) with p 2 = 0.
In principle, the second-order differential operators involved in the eigenequations of motion are obtained from the eigenvalues of the principal symbol σ D by using the correspondance p µ ↔ ∂ µ . This can be directly achieved for λ 1,2,3 , yielding the eigenoperators (4.20) The fact that one operator is zero is not suprising, because it translates that one of the degrees of freedom of the vector field is non-dynamical. As solutions of Eq. (4.19), the last two eigenvalues of the principal symbol generally involve square roots, and it is unclear how one should interpret them in terms of differential operators. In the case f 3 = 0, however, σ D as written in Eq. (4.15) is already diagonal, and the remaining differential operators read

Hyperbolicity of the eigenequations
It turns out that the third eigenoperator D 3 is actually sufficient to rule out the f 3 -term. Indeed, consider for instance a state with a purely homogeneous scalar field, then where ∆ ≡ ∂ i ∂ i is the Laplacian. We know from the stability analysis that f 3 ≥ 0, but if there exists a value of φ so that f 3 (φ) > 0, then forφ large enough, D 3 becomes elliptical. Therefore, it is necessary to have f 3 = 0 for the theory to be both stable and causal. For f 3 = 0, the eigenoperators D 1...4 are all proportional to , which is hyperbolic with a timelike (−)-direction, thus the causality requirement does not impose further constraints on f 1 . Concerning f 0 , additionnally to the condition ∂f 0 /∂K ≥ 0 imposed by the Hamiltonian stability requirement, we must also have for the eigenoperator D 0 of Eq. (4.23) to be hyperbolic. We propose, in Appendix A, a simple proof of the above condition, which is well known 7 in the the context of k-essence [68][69][70][71][72].

Conclusion and further remarks
We have derived necessary conditions for the stability and causality of models built from one scalar field and one vector field coupled to each other. Under the restrictions stated in Sec. 2, we showed that the most general action describing a stable vector-scalar dark sector, whose dynamics is ruled by hyperbolic equations of motion, reads with K ≡ ∂ µ φ∂ µ φ, and where the coupling functions obey: • ∂f 0 /∂K ≥ 0, and f 0 (φ, K ≥ 0) bounded by below (stability); • f 1 (φ) ≥ 0 (stability); 7 Note, by the way, that the discussion about hyperbolicity in one of the first reference article [68] is partially wrong. Indeed, the authors claim that the hyperbolicity conditions are (a) ∂f0/∂K > 0, and (b) ∂ 2 f0/∂K 2 ≥ 0, which is not really the case: (a) is rather imposed by the stability condition, and (b) does not exist at all. They also mention Ineq. (4.26), but as a condition which "assures the stability of the Cauchy problem-that is, small changes in the Cauchy data cannot produce large changes in the solution arbitrarily close to the initial surface", whereas it is actually necessary for the Cauchy problem to exist.
There are no further restrictions over the coupling function f 2 . It is remarkable that the class of models satisfying the assumptions of Sec. 2 are so constrained by the basic principles of stability and causality. However, it is worth noting that the theories excluded by our analysis are really ruled out only if one considers them as fundamental. If, on the contrary, they represent the effective behavior of a more fundamental but healthy theory, then the only requirement is a reasonable domain of stability and causality. By essence, the present work cannot draw any definite conclusion within the world of such effective theories.
Gauge invariance was a central assumption in our analysis. We shall mention that an important issue with this property in vector-field models was pointed out in Ref. [73], where the authors considered a possible generalization of electromagnetism in Minkowski spacetime, inspired from scalar Galileon theories. Their conclusion came in the form of a "no-go theorem" for generalized vector field Galileons, which states that it is impossible to construct more general theories than standard electromagnetism, because all possible extensions following the Galileon construction procedure lead to topological or boundary terms, and are thus nondynamical. In order to escape this theorem, one can however build models with multicomponent gauge-invariant vector fields, or couple the vector field with another field, e.g., a scalar field as done in this article. The coupling of different types of fields with non-trivial dynamics was addressed earlier in Refs. [62,74], while Ref. [75] proposed a complete study of scalar Galileons with gauge symmetries.
One may then wonder what models can be built once the condition of gauge invariance is removed. References [64,65] have recently addressed the problem of gauge-invariance breaking for single-vector-field models, in the spirit of Galileon theories. These analyses conclude that, for some particular combinations of the non-gauge invariant terms, ghost-like instabilities disappear and it is possible to obtain a well-behaved Galileon-type generalization of the Proca theory with three physical propagating degrees of freedom. In general, dropping gauge invariance in a vector-scalar theory leads to a system with more physical degrees of freedom, the dynamics of which can be governed by a huge variety of terms in the action, corresponding to all the possible contractions formed out of A µ , ∂ µ φ, ∂ µ A ν , such as A µ A µ , Given that the structure of each of these terms is very different, there is a priori no general procedure to deal with all of them together, so that one should probably perform a dedicated analysis of the stability and causality for each model. Nevertheless, we can point out that some particularly interesting terms are often studied, and admit a simple analysis; for instance, vector potentials of the form V (A 2 ), or couplings of the form A µ ∂ µ φ, ∂ µ A µ . A very specific class of models of the form f (F 2 ) + V (A 2 ) was recently studied in Ref. [76], and where shown to have hyperbolic equations of motion for some special regions of phase space.
Finally, we should mention that our analysis admits a straightforward generalization to multiple vector fields which are gauge invariant under a non-Abelian gauge group. However, in the non-Abelian case we have an important difference with respect to the U(1) Abelian case. As shown in Subsec. 2.2, in the U(1) case any term of the action involving the vector can always be reduced to even powers of the Faraday tensor or its dual, possibly contracted with the derivatives of the scalar field, which can be further reduced to products of X, Y and Z. All the odd products of the Faraday tensor and its dual are identically zero. In presence of a non-Abelian gauge group, this is no longer the case, and there appear other combinations which add non-trivial dynamics to the system. Among the lowest order terms, there appear for instance cubic combinations of the form F F F and F FF : where C bc a are the structure constants of the group and a, b, c are Lie algebra indices. These terms are consistent with gauge symmetries and are dynamical. The term F F F appears generically in non-Abelian gauge theories and terms like F FF appear for instance in QCD when discussing CP-violations originated by gluonic operators of dimension six [77]. Recently, the dynamics of such terms was also studied in the context of leptogenesis in non-Abelian gauge fields populated models [78]. Despite of the numerous ways in which the the fields can interact when gauge invariance is broken or when non-Abelian gauge groups are considered, it is expected that the stability and the causality analysis would impose constraints over all those possible interactions and it would be interesting and valuable to extend the methods followed here to those cases.
3. K = 0. This case is slightly less trivial than the previous two ones. Up to a spatial rotation, we can write ∂ µ φ =φ(−δ µ t + δ µ 1 ). The differential operator G µν ∂ µ ∂ ν can then be diagonalized using the two vectors where we denoted f 0 ≡ ∂f 0 /∂K for short. Using the basis (∂ − , ∂ + , ∂ 2 , ∂ 3 ), the differential operator indeed becomes Thus, the differential operator is always hyperbolic. The question is now whether the (−)-direction is timelike or spacelike; it is immediate to check that so that ∂ − is timelike (and ∂ + is spacelike) if and only if f 0 = ∂f 0 /∂K ≥ 0, which is consistent with the condition found in the previous two cases, for K = 0.
B Beyond linearity in X, Y , Z In this appendix, we provide some results that can be a starting point for the analysis of more general models than the ones described by action (2.16). Namely, consider a dark-sector Lagrangian density L DS (φ, K, X, Y, Z) which is not necessarily linear in X = F 2 , Y = FF , and Z = (F ∂φ) 2 . In the following, we drop the 'DS' label to alleviate notation.

B.1 Hamiltonian density
The canonical momenta conjugate to the scalar and vector fields are respectively π φ = −2L ,Kφ + 2L ,Z îφ B 2 + det(∇φ, E, B) ó , (B.1) where a coma denotes a derivative. The constraint is unchanged compared to the case where L is linear in X, Y , Z, namely ∇ · π = 0. The Hamiltonian density then reads It would be tempting to conclude that L ,K L ,X ≤ 0, and L ,Z ≥ 0 are necessary conditions for H DS to be bounded by below, but unfortunatelyφ 2 , E 2 , (φB − E × ∇φ) 2 , K, X, Y , Z are not independent variables, so that one cannot take, e.g.,φ 2 → ∞ while keeping the others finite. The actual stability conditions could be much subtler, for instance they could involve combinations of the derivatives of L, and thus require a dedicated study.