The Maxwell-Proca theory: definition and construction

We present a systematic construction of the most general first order Lagrangian describing an arbitrary number of interacting Maxwell and Proca fields on Minkowski spacetime. To this aim, we first formalize the notion of a Proca field, in analogy to the well known Maxwell field. Our definition allows for a non-linear realization of the Proca mass, in the form of derivative self-interactions. Consequently, we consider so-called generalized Proca/vector Galileons. We explicitly demonstrate the ghost-freedom of this complete Maxwell-Proca theory by obtaining its constraint algebra. We find that, when multiple Proca fields are present, their interactions must fulfill non-trivial differential relations in order to ensure the propagation of the correct number of degrees of freedom. These relations had so far been overlooked, which means previous multi-Proca proposals generically contain ghosts. This is a companion paper to arXiv:1905.06968 [hep-th]. It puts on a solid footing the theory there introduced.

In the second half of the 19th century, the foundations of classical electromagnetism were laid. The Lorentz force law, together with Maxwell's equations, provided a unified framework to explain the origin and propagation of electric and magnetic fields. Light was then understood to be an electromagnetic phenomenon. This body of work has since become one of the cornerstones of theoretical physics. It is therefore not surprising that ongoing efforts are made so as to generalize Maxwell's theory. In the following, we discuss the subset of extensions that provides the relevant ingredients to our new proposal: the Maxwell-Proca theory. (If the reader arrives from [1], we now elaborate on the spin-one void noted in that introduction.) Early extensions. As a field theory, Maxwell electrodynamics describes an Abelian massless vector field and its linear interactions with sources. During the first half of the 20th century, extensions to this setup were introduced in order to formalize new observations. For instance, consider nuclear phenomena. As pointed out by Yukawa [2], the forces there implicated do not obey a Coulomb law, but show a faster decay over distance. This insight was in good agreement with the theory of a massive vector field introduced by Proca [3]. In his paper, Proca proposed a relativistic massive wave equation for a vector field (analogous to the Klein-Gordon equation [4]). Subsequent work [5] refined the implications of his proposal: it was to be understood as a model of a massive spinone force carrier, i.e. a massive photon. Of course, the phenomenologies of the Maxwell and Proca fields are different. For example, the Proca field explicitly breaks the * Electronic address: veroerdi@mppmu.mpg.de † Electronic address: brageg@mppmu.mpg.de ‡ Electronic address: julioamz@mpp.mpg.de § Electronic address: angnissm@mpp.mpg.de gauge symmetry and thus propagates one more degree of freedom than the Maxwell field. Alternatively, this helicity-zero mode can be excited without spoiling gauge invariance. Shortly after Proca's work, Bopp [6] and Podolsky [7] presented a model of a U (1)-invariant Abelian vector field with a modified kinetic structure. Such novel electrodynamics (and generalizations [8]) was introduced to amend some of the theoretical concerns at the time. Specifically, this was a consistent proposal describing short distance interactions and thus resolving the self-energy problem -a neat review of which can be found in [9]. Even though the proposed Lagrangians contain higher order self-derivative interactions, the equations of motion remain second order 1 . However, some other peculiar characteristics appear, such as non-locality.
Many other theories of relativistic vector fields emerged. A prominent representative is the theory in which the Maxwell kinetic term is promoted to an arbitrary smooth function of itself. Concretely, an early proposal of this kind was Born-Infeld electrodynamics [11]. In this model, the gauge symmetry is preserved and the equations of motion are again second order. Besides, Born-Infeld theory (and closely related constructions like the Heisenberg-Euler model [12]) manages to resolve the infinities in the self-energy of the electron. This is achieved at a price, though: the field equations become non-linear. In fact, all theories sharing such properties are commonly labeled non-linear electrodynamics [13]. We stress that mass, derivative self-interactions and nonlinearity are constituents of the Maxwell-Proca theory here constructed.
Modern proposals. In recent times, attempts have been made in pursuit of the most general theory of a massive vector field whose self-interactions extend the standard Proca mass term. Such theories are collectively known as generalized Proca or vector Galileon -a name that will be clarified later. This picture was first suggested in [14]. The approach is based on the inclusion of derivative selfinteractions, carefully chosen so as to not excite further modes than the three expected in four spacetime dimensions. Throughout the paper, we refer to any further degree of freedom carried by the action as a ghost. The motivation behind these models goes beyond extending previous electromagnetic-like theories: when considered over a cosmological background, generalized Proca fields can give rise to self-acceleration and fit data from latetime observations [15]. Evidently, a coupling to gravity is required for such applications. Even though this coupling is already used in the initial references [14], a rigorous counting of degrees of freedom is not carried out there. This issue has been further investigated [16,17], but still much is to be done in this direction 2 .
In closer detail, the above mentioned theories are defined by the Lagrangian density of an Abelian spin-one field that lives on Minkowski spacetime. Schematically, the Lagrangian can be split into two parts: The first piece L (0) is an arbitrary smooth function built from the field itself and its (dual) field strength. As such, it generalizes both non-linear electrodynamics and the original Proca theory. On the other hand, the second piece L admits a closed form and contains selfinteractions composed of contractions between the field and its derivatives. Contrary to L (0) , the derivatives within L do not come as field strengths exclusively, but they are such that the Lagrangian remains first order. In fact, it has been the work of a couple of years to track down the complete set of interactions comprising the soclaimed most general L [14,19,20]. The interested reader can find thorough reviews on these theories in [21]. The generalized Proca theories introduced so far share a characteristic that has been exploited very actively. When inspecting the longitudinal mode's dynamics in the decoupling limit, one finds that it matches the acclaimed scalar Galileon. Based on this feature, generalized Proca theories are interchangeably labeled vector Galileons. First introduced in [22], Galileon fields also arise in the context of modified gravity, for instance in models that give an effective mass to the graviton [23]. Both scalar and vector Galileons are close relatives of 2 There is a wider collection of papers related to vector Galileons in the presence of gravity. Their aim is not to prove ghost-freedom. Instead, submodels are taken as a source to tackle specific goals, such as finding black hole or star configurations [18]. While the models that deal with a single vector field and yield second order equations of motion are automatically ghost-free, we highlight both the pertinence and non-triviality of the pending proof of ghost-freedom in the remaining proposals.
Horndeski's construction. The latter gives the most general theory of a self-gravitating scalar [24] or vector (gauge) field [25]. Indeed, when Galileon-like theories are promoted to evolve in the presence of gravity, they must couple to the curvature according to Horndeski's prescription. Conversely, Galileons can be understood as the zero-curvature limit of a Horndeski theory. Finally, theories of multiple spin-one fields have also been investigated, albeit they have not received that much attention. Interacting Maxwell fields were alluded to in the search of one-form Galileons in [26]. In [27,28], theories that comprise several copies of generalized Proca fields were introduced. There, a global rotational symmetry was imposed on the field space, which drastically reduced the allowed interactions. The motivation behind this restriction was to provide an innate source capable of supporting isotropic cosmologies. To our knowledge, there are no models incorporating couplings between massless and massive vector fields.
Summary of results. In this work, we focus on ghost-free interactions of multiple real Abelian fields on Minkowski spacetime. Loosely speaking, these can be either massless (Maxwell) or massive (Proca). We formalize both notions by providing their canonical definitions in terms of Hamiltonian constraints. Then, we present a bottom-up construction of the most general first order Lagrangian involving such fields. Our result, the Maxwell-Proca theory, is axiomatically complete. It naturally contains Maxwell electromagnetism, non-linear electrodynamics and Proca's theory, together with its modern generalizations, as particular subcases. The general structure of the Maxwell-Proca interactions is similar to the generalized Proca action in (1), but it is subject to highly restrictive differential conditions on the form of L when multiple Proca fields are involved. The origin of these conditions lies at the secondary level of the Lagrangian constraint analysis, which we explain in detail.
Organization of the paper. In Sec. II, we introduce the mathematical formalism behind our setup. All information here presented is based on the Lagrangian-Hamiltonian constraint analysis in Appendix A. In Sec. III, we systematically generate the exhaustive set of ghost-free interactions. We continue in Sec. IV, discussing direct strategies to enhance the generality of the Maxwell-Proca theory. At last, in Sec. V, we summarize the key ideas behind the technical details of our construction.

II. AXIOMATIZATION
The Maxwell and Proca actions are the elemental theories of Abelian vector fields: a massless one in the former case and a massive one in the latter. In d-dimensional Minkowski spacetime and free of sources, the Lagrangian densities of these theories have the form respectively. We take both fields to be real. Here, Similarly, B µν is that of the Proca field and the constant m = 0 parameterizes its mass. We begin by convening the essentials of these theories, since they serve as building blocks for the more elaborate Maxwell-Proca theory we propose. Maxwell electrodynamics is manifestly invariant under the U (1) gauge transformation of the second kind where Φ is an arbitrary real scalar field. This means that A µ transforms as a connection, in the affine representation of the gauge group. We highlight that both the Lagrangian (off-shell) and the field equations (on-shell) are invariant under (3). The explicit mass term that tells apart both Lagrangians in (2) has the crucial effect of spoiling such gauge symmetry for the Proca field. One of the central points throughout the paper is the count of the degrees of freedom propagated by a given Lagrangian, as a means to ensure the absence of ghosts. It is well known that L M propagates d − 2 whereas L P propagates d − 1 degrees of freedom. In four spacetime dimensions, Maxwell's theory accounts for two transversal polarization states, while the Proca field incorporates an additional longitudinal mode.
Both L M and L P are singular. In other words, the determinant of their primary Hessians, vanishes. Here,Ȧ µ := ∂ 0 A µ stands for the velocity of the massless field whileḂ µ corresponds to that of the massive field. As a consequence, some of the Euler-Lagrange equations following from (2) are first order. Namely, they are not equations of motion -in the sense that they do not involve accelerations of the fields-but relations between the field variables and their velocities. These relations are the primary Lagrangian constraints 3 of the 3 This discussion applies to field theories with irreducible con-system (if tautologically true, they can be further distinguished as identities). They are responsible for the mismatch between the number of a priori degrees of freedom and those that are physically meaningful. Given any singular Lagrangian, the existence of primary constraints ϕ = 0 that are not identically zero starts a renowned iterative procedure to obtain the remaining constraints. This algorithm is suitably referenced and self-consistently reviewed in Appendix A 1. The gist of it is as follows. The demand that the primary constraints be stable in timeφ ! = 0 yields the so-called secondary Euler-Lagrange equations. Those which are independent of accelerations, if any, are the secondary Lagrangian constraints φ = 0. When extant, the stability of those secondary constraints that are not identically vanishing is not guaranteed. Rather, requiringφ ! = 0 leads to the third iteration. The procedure terminates when no Lagrangian constraints appear or when they are all identically zero, since then their stability need not be assured.
Going back to the Maxwell and Proca theories in (2), the just described constraint analysis tells us that Maxwell electrodynamics is endowed with a single primary (Bianchi) identity, while its massive counterpart contains two Lagrangian constraints: one primary and one secondary. Furthermore, the generator of the gauge symmetry inherent to Maxwell's theory can be regarded as an additional primary identity.
Knowledge of the identities and Lagrangian constraints in a singular theory can be used to determine the degrees of freedom it propagates. Most often, this is done by means of the Hamiltonian formula proposed by Dirac long ago: where M 0 is the number of field variables and (N 1 , N 2 ) denote the number of first and second class constraints, respectively. Recall that a first class constraint is one that yields a vanishing Poisson bracket with all of the constraints in the theory, whereas a second class constraint does not. It is worth noting that in [29] an alternative expression to (5) was developed for particle systems, which depends solely on Lagrangian quantities, including the off-shell gauge generators. The result was later on adapted to field theories in [30]. Indeed, Appendix A 2 is devoted to the exposition of a long established systematic procedure to calculate N 1 and N 2 from the Lagrangian constraints. As is widely known, the Bianchi and gauge identities in Maxwell's theory generate two first class constraints. On the Proca side, the Lagrangian constraints translate to two second straints. All L M , L P and the later introduced Maxwell-Proca theory pertain to this class. The case involving reducible constraints is briefly reviewed in Appendix A 1. class constraints. A quick inspection of (5) then yields the anticipated number of degrees of freedom.
Considering all the above exposed, we introduce the following two definitions: Definition 1. A Maxwell field is any real Abelian vector field whose Lagrangian is associated with two first class constraints.
Definition 2. A Proca field is any real Abelian vector field whose Lagrangian is associated with two second class constraints.
From this point of view, the distinction between a Proca field [3,5] and a generalized Proca field [14,19,20] is insubstantial. Subsequently, we do not discern between the two.
Notice that, in the above definitions, the Poisson bracket is the relevant operation: it distinguishes between first and second class constraints and thus defines the constraint algebra. Bilinearity of the Poisson bracket implies linearity of the constraint algebra. Consequently, a direct sum of whatever number of copies of the Maxwell and Proca Lagrangians will preserve the constraint patterns in Definitions 1 and 2. This is generically not true when there is an explicit interaction involving (some of) the Proca fields. Specifically, we prove in Appendix A 1 that coupling massive vectors alters the stability of their primary Lagrangian constraints. What happens is that, at the secondary level, the constraints φ = 0 are absent. Instead, we get independent Euler-Lagrange equations responsible for exciting unwanted ghost-like degrees of freedom.
One of our key results lies in showing that one can always find befitting interactions that prevent the aforementioned deviation from the Proca constraint scheme. This is uniquely achieved by setting the secondary Hessian to zero. In more detail, consider some Lagrangian density L ′ encoding the interaction of a set of M Proca fields {B µ }. The antisymmetric matrix conforms the non-vanishing part of the secondary Hessian. Enforcing induces the appearance of the secondary constraints that restrain the ghost. We stress that only when (7) is fulfilled, can the vectors {B µ } composing L ′ be identified as Proca fields in the sense of Definition 2. This prescription has been irreconcilably overlooked in the existing literature [27,28]. As a consequence, the models there submitted are ghost-full. The preceding observations motivate the proposition of our main result: ii The Lagrangian is at most first order.
By a first order Lagrangian density we mean up to boundary terms. Then, the variational principle will produce at most second order field equations. This is a sufficient condition to safeguard our theory from Ostrogradsky-like instabilities, as already noted in footnote 1. When a single Proca field is considered (i.e. M = 1), this is the necessary and sufficient condition.
The most concise form of the Maxwell-Proca theory is where the kinetic part can always be expressed as and the interaction part encompasses every possible Lorentz scalar built out of the fields and their first derivatives, while respecting the defining axioms. We identify three sectors: massless-massless, massless-massive and the massive-massive. Explicitly, We reassert that the Maxwell-Proca theory, in its final form of Sec. III C, is only ghost-free when regarded as a whole together with the condition that its secondary Hessian vanishes. Whilst the massless interactions L (AA) trivially fulfill this condition, L (AB) and L (BB) don't: they are severely restricted by it.

III. CONSTRUCTION OF THE GHOST-FREE INTERACTIONS
In this section, we systematically construct the interactions (11) building on the definitions given in the previous section. We first review the necessary and sufficient conditions to ensure the existence of the suitable constraints at the primary level. Albeit slightly more involved, this construction is very similar to that for a single Proca field [14]. It relies on the antisymmetric structures first recognized in [26] as a necessary and exhaustive tool for developing healthy Lorentz invariant theories. After implementing consistency at the primary level, we present the highly non-trivial conditions that arise at the secondary level of the constraint algebra. We discuss how these conditions notably restrict multi-Proca interactions and apply them to generate ghost-free interactions in several examples.

A. Primary constraints
For our present purposes, it is useful to introduce a collective notation for all Abelian vector fields X (A) µ , where A ranges from 1 to M + N . This means that, for the moment, we set aside the issue of gauge invariance.
As already stated, healthiness of the Maxwell-Proca theory necessarily relies on a suitable rank reduction of its primary Hessian. This is uniquely 4 achieved by Notice that (12) must hold for all µ. This is true also for the single Proca limit (N = 0 and M = 1), as described, for instance, at the very beginning of page three in the review paper [20]. We refer to the above as the primary constraint enforcing relations, which are not to be confused with the primary constraints themselves. Together with the demand for Lorentz invariance, they will serve as a guiding principle to forge the interactions. They are explicitly derived in the appendix. Specifically, they are equivalent to (A7). Next, recall that the Maxwell-Proca theory is of the form (9). The kinetic part L kin given in (10) automatically satisfies the relations in (12). We thus focus on the interactions L int . In order to fulfill (12), they can have at most a linear dependence on all the velocities associated with temporal componentsẊ (A) 0 . It is then natural to split L int into terms L (0) which do not depend onẊ (A) 0 and terms L int which carry a linear dependence on these velocities: By definition, L (0) is an arbitrary real smooth function of the fields and their field strengths: Due to these loose restrictions, it contains infinitely many terms.
We now turn to the composition of L int . For a first order Lagrangian and under the assumptions of smoothness and reality, L int can be formally expanded as with Here, n counts the number of derivatives and (A 1 , A 2 , . . . ) label all Abelian fields, as (A, B, . . . ) before. The expansion coefficients T ν1...ρn A1...An have an arbitrary dependence on the fields but contain no derivatives. Note that, if all derivatives appear only in field strengths, the corresponding term belongs in L (0) by definition. Observe that the exchange of any triple (A i , ν i , ρ i ) ←→ (A j , ν j , ρ j ) leaves (16) invariant. Hence, without loss of generality, we take the expansion coefficients to be symmetric under this replacement.
In order for L int to satisfy (12), it is necessary and sufficient that does not contain any time derivatives. For L (1) this is automatically true, while for with n > 1, it must be enforced. The only Lorentz covariant way to do so is to entertain three distinct antisymmetrizations in (16): i Antisymmetrization of all n derivatives with each other.
ii Antisymmetrization of all but one derivatives with each other and with all vector field indices. Then, the respective n − 1 fields appear only inside field strengths, while the derivative on the nth field can be contracted in an arbitrary manner.
iii Antisymmetrization of k < n derivatives with l < n fields, together with the antisymmetrization of the remaining n − k derivatives with the left-over n − l fields.
All these options restrict but do not fully determine the coefficients T ν1...ρn A1...An . Hence, the interactions contain free tensorial functions.
On four-dimensional Minkowski spacetime, the antisymmetrizations are uniquely accomplished using the totally antisymmetric Levi-Civita tensor ε µνρσ . The above considerations show that all L (n>1) must contain at least one Levi-Civita tensor. There is another rather drastic consequence: since we always need to antisymmetrize at least n indices with each other, the in principle infinite series expansion (15) truncates at n = 4 in four spacetime dimensions. Thus, the full interaction Lagrangian reads For concreteness, we explicitly write out all independent Lorentz contractions that respect the primary constraint enforcing relation.
As stated earlier, for L (1) , there is no condition on the expansion coefficients: where henceforth the spacetime tensors {F } are arbitrary smooth real functions 5 . Generically, they do not possess any symmetry properties and they are built from the fields X (A) µ and the Lorentz invariant objects η µν and ε µνρσ . The allowed coefficients in L (2) are of the form where each line corresponds to the cases i-iii, respectively.
σ1σ2 . (22) Here, the first term implements case i. Case ii cannot be realized for n > 2 in four dimensions. The remaining terms all belong to and exhaust case iii. Finally, in L (4) , where the first term is associated to case i while the remaining lines fully cover case iii.
To summarize, (19)-(23) define the complete set of interactions that satisfy the primary constraint enforcing relations.

B. Secondary constraints
Our most crucial result, arising from the constraint analysis outlined in Appendix A 1, is the unveiling of a set of secondary constraint enforcing relations. Namely, the 5 They need not be tensors with respect to the field labels. For example, F A νρ = η νρ for all A is a valid choice for (20).
axiomatically consistent closure of the constraint algebra requires the vanishing of the secondary Hessian with L int as in (13). Together with (12), the above gives the necessary and sufficient conditions to avoid the propagation of ghosts. Then, the algebra closes automatically at the tertiary level. The differential relations (24) further restrict the internal structure of the interactions constructed in the previous section, which will be shortly exemplified.

C. Gauge invariance
We now identify each Abelian vector X (A) with either a Maxwell or a Proca field. Then, the interactions (19)-(23) can be classified as in (11): L (AA) describing purely massless interactions, L (BB) for purely massive ones and L (AB) containing the interactions between massless and massive fields. In doing so, we need to respect gauge invariance.
All N Maxwell fields A (ᾱ) µ transform as connections under their own U (1) gauge group, as in (3). It has long been established that consequently the Maxwell vectors can only show up inside their field strengths. As a result, the purely massless interactions are a subset of the general zeroth interactions: µ , on the other hand, are gauge invariant by definition. Therefore, the purely massive sector L (BB) is straightforwardly obtained by replacing the fields X (A) with Proca fields: . (25) The coupling L (AB) between massless and massive fields is obtained in a similar way: replacing X (A) with either A (ᾱ) or B (α) , with the Maxwell fields appearing only inside field strengths. In particular, this means that the functions F depend exclusively on the Proca fields.
In terms of Maxwell and Proca fields, the secondary Hessian takes the form since Maxwell fields, which appear only inside field strengths, automatically satisfy (24). For the same reason, the interactions L (0) are not restricted at the secondary level. In other words, (24) solely constrains the n > 0 terms in L (AB) and L (BB) by Notice that only when two or more Proca fields are coupled, does (27) become non-trivial. The derivation of the above relations can be found in the appendix. Concretely, (27) is a rewriting of (A20). The nature and relevance of these relations can be grasped as follows. Equation (12) enforces the existence of a set of suitable (primary) constraints. Then, equation (27) guarantees their stability. The (secondary) constraints that follow from equation (27) hold true at all times, i.e. they are automatically stable. Therefore, no further constraint enforcing relations apply.
The internal structure of the expansion coefficients in (20)-(23) is so rich that implementing (27) in full generality seems intractable, because these conditions form a set of coupled non-linear partial differential equations. Hence, we do not incorporate these relations into a building principle like we did at the primary level. Instead, we provide a few examples of how to use them to extract ghost-free interactions.

D. Explicit examples
First, consider L (BB) (1) , which remained unconstrained at the primary level. It can be written in the form where the coefficient is an arbitrary function of the Proca fields obeying F α µν = F α νµ , so as to exclude terms that belong in L There are two straightforward solutions to this equation. The functions F α µν can be constants or they can depend exclusively on the field with label α: F α µν = F α µν B (α) for all α. In the first case, the interaction becomes a boundary term. In the second case, we are left with M decoupled Proca theories. While it is difficult to find the most general form of F α µν , one can always construct particular solutions. A non-trivial solution fulfilling (29) is given in [1].
Next, we take an example from [27]. In our notation, their term reads and belongs in L (BB) (2) . For this interaction (which passes the primary level by completely antisymmetrizing the derivatives), the secondary level (27) requires Clearly, this condition is not satisfied for generic field configurations. The interaction term spoils the closure of the constraint algebra and propagates additional degrees of freedom. We conclude that (30) and various other interactions proposed in [27,28] carry ghosts. We finish by providing a few examples of ghost-free interactions. Any Lagrangian compatible with the primary constraint enforcing relations (12) that contains a single Proca field (M = 1) immediately passes the secondary level as well. This is true for any number N of Maxwell fields. Accordingly, a large number of healthy terms can directly be obtained from (19)- (23). A purely multi-Proca example (M > 1, N = 0) is which belongs in L (BB) (2) and has An L (AB) (2) example with M, N > 1 is given by where the expansion coefficients are set to Here, the K are arbitrary constants and It is easy to check that the interactions (32) and (34) satisfy (27) and thus do not propagate any ghost-like degrees of freedom. We emphasize that this holds true because (33) and (36) have been tailored to be in agreement with the secondary constraint enforcing relation.
In closing this section we recall that, for interactions involving more than one Proca field, one must always ensure that (27) is satisfied. This imposes highly nontrivial restrictions on the interactions. ii All vector fields are taken to be real.
iii The Lagrangian is at most first order.
iv We work on Minkowski spacetime.
To conclude, we comment on the scenarios that unfold upon concrete relaxations of the assumed postulates. In all cases, the materialization of our subsequent proposals into an exhaustive Lagrangian is non-trivial and lies beyond the scope of this paper.
Generalization to non-Abelian fields. As is widely known, classical Yang-Mills theory [31] can be regarded as an extension of Maxwell electrodynamics to non-Abelian vector fields. Conversely, electromagnetism can be understood as a Yang-Mills theory with gauge group U (1). This relationship has already been exploited to construct generalized Yang-Mills theories starting from specific non-linear electrodynamics [32]. On the other hand, non-Abelian extensions of the Proca field are not possible. Indeed and as already noted in [21,27], attempts at constructing a massive version of the Yang-Mills theory boil down to the consideration of multiple interacting Abelian Proca fields in the presence of a global rotational symmetry in the field space. Putting together both observations, the non-Abelian analogue to the Maxwell-Proca theory could be dubbed Yang-Mills-Proca. Notice that the constraint enforcing relations (12) and (27), as well as the purely massive interactions L (BB) in (11), are the same as in the Maxwell-Proca theory.
Generalization to complex fields. While the extension of the massless sector to the complex field case is straightforward, the study of a complexified massive sector is not. For example, allowing for Proca fields to be complex makes it possible to charge them under the fundamental representation of the U (1) gauge group of any of the Maxwell fields: with e a coupling constant and Φ the same real scalar field that appeared in (3) before. Then, one can generate new interaction terms by replacing partial derivatives ∂ µ acting on the Proca fields by covariant derivatives D µ := ∂ µ + ieA µ , as long as the resulting term is invariant under the simultaneous transformations (3) and (37). Observe that such interactions imply the Maxwell fields are no longer restricted to appear as part of a field strength. This generalization does not affect the form of the constraint enforcing relations (12) and (27). Generalization to higher order field theory. If the Maxwell-Proca theory is promoted to depend on arbitrarily high order time derivatives of the fields, i.e.
the Lagrangian and Hamiltonian formalisms underlying Appendix A no longer constitute the suitable framework to study the constraint algebra. Indeed, the development of an appropriate formalism has been investigated for a long time and continues to be pursued, see e.g. [33]. The foundations of a consistent approach amenable to geometrization were laid in [34], based on the notions put forward in [35] for second order theories. The main idea consists on introducing so-called intermediate spaces that connect the tangent bundles T 2k−1 C, natural to the Lagrangian analysis, to the cotangent bundles T * (T k−1 C), which appear on the Hamiltonian side. Here, k = 1, . . . , n for a nth order Lagrangian and C is the configuration space. The intermediate spaces sustain a systematic inspection of the higher velocities and a step-by-step inference of the corresponding momenta. Moreover, iterative algorithms that ensure the stability of the constraints in the theory are naturally defined in these spaces, under mild requirements on the form of (38). Obviously, the constraint enforcing relations (12) and (27) do not apply to higher order cases. Generalization to curved backgrounds. Relaxing our last axiom iv means allowing for general background configurations, such as (anti-)de-Sitter or Friedmann-Robertson-Walker spacetimes. The complete set of permitted metrics can be obtained by coupling the Maxwell-Proca theory to gravity. As discussed elsewhere [1], this leads to diverse physical applications.
Regarding free vector fields with a standard quadratic mass, a consistent coupling to General Relativity is obtained by direct covariantization. For instance, the Proca mass term in Minkowski space m 2 η µν B µ B ν is promoted to m 2 g µν B µ B ν , where g µν is (the inverse of) the dynamical spacetime metric. For terms containing two or more derivatives, this procedure becomes ambiguous due to the presence of the derivatives ∂ µ , which do not commute anymore when replaced by covariant derivatives ∇ µ .
When covariantizing the Maxwell-Proca action one needs to be particularly careful to not destroy the underlying constraint algebra. This algebra should contain the constraints for the massless spin-two together with those of the Maxwell-Proca theory. In other words, we expect the constraint enforcing relations (12) and (27) to be supplemented by additional conditions coming from the gravitational sector.
In fact, already for the single Proca theory in [14], nontrivial curvature couplings appear in the action. The precise form of these terms was motivated by focusing on the longitudinal mode of the vector field and requiring that its couplings to gravity be the same as in Horndeski's theory, which is known to be ghost-free. This is a necessary but not sufficient condition for a (first order) Lagrangian to propagate the correct number of degrees of freedom.
The authors of [16] took a more rigorous approach and initiated an ADM analysis [36] of the interactions proposed in [14]. Their results show that the non-minimal curvature couplings are precisely such that both the B 0 component of the Proca field and the lapse function in the metric appear without time derivatives in the action. This demonstrates that the Lagrangian in [14] satisfies another necessary condition for consistency. A full constraint analysis of the theory that provides a solid proof of ghost freedom, however, is still pending. Based on these partial results for the single Proca field, it is natural to expect that the ghost-free covariant version of our Maxwell-Proca theory will also contain non-minimal curvature couplings.
We finish by pointing out two recent developments pertinent to the beyond Maxwell-Proca scenarios above described. First, it is interesting to note that a particular example of a consistent theory that simultaneously relaxes axioms iii and iv arises in [37]. One of the decoupling limits of non-linear massive gravity [38] there considered distributes the theory's five propagating modes into two plus three degrees of freedom, corresponding to a massless graviton and a Proca field, respectively. The resulting vector field Lagrangian contains higher order terms and is defined on four-dimensional anti-de Sitter spacetime, meanwhile it is ghost-free by construction, as proven in [39]. Second, when relaxing axiom iv and coupling the Maxwell-Proca theory to gravity, one should be careful not to miss the kind of terms unveiled shortly after this work, in [40]. Such terms reduce to total derivatives in flat spacetime.

V. CONCLUDING REMARKS
The key ideas behind the results here obtained are as follows. We begin by formalizing the notion of a Proca field in Definition 2, putting it on the same footing as that of the renowned Maxwell field (see Definition 1). Based on these two concepts, we construct the Maxwell-Proca Lagrangian: the most general first-order theory in Minkowski spacetime for an arbitrary number of interacting Maxwell and Proca fields propagating the correct number of degrees of freedom.
In more detail, the Maxwell-Proca Lagrangian is most conveniently expressed as the sum of kinetic plus interaction terms (9). The kinetic piece is (10). The interactions are obtained by systematically converting the uniquely defined primary constraint enforcing relations (12) -derived in the appendix (A7)-into as definite a Lagrangian as possible. We distinguish two types of interactions in (19): those that trivially satisfy (12) and those that don't. The first set can succinctly be written as (14), while the coefficients of the second set are made explicit in (20)- (23).
Additionally, when multiple interacting Proca fields are considered, this Lagrangian must be forced to satisfy the uniquely defined secondary constraint enforcing relations in (27) -derived in the appendix (A20)-to avoid ghosts. No further conditions on the Lagrangian apply because, once equation (27) is fulfilled, the constraint algebra automatically closes at the tertiary level. Notice that, in full generality, (27) is a set of coupled non-linear partial differential equations whose solution is not known. This explains our inability to turn (27) into a Lagrangian building principle, unlike (12) before.
The vast and diverse applicability of the Maxwell-Proca theory is discussed and exemplified elsewhere [1].

Acknowledgments
The authors thank Axel Kleinschmidt for a fruitful discussion, which gave us the idea for this work. We are grateful to Chrysoula Markou for enlightening conversations in the early stages of the project. VED thanks Nicolás Coca López for lucid explanations regarding superconductors. JAMZ is thankful to E. A. Ayón-Beato and D. F. Higuita-Borja for the computing power granted at ZymboLab (CINVESTAV). The work of JAMZ is partially supported by the "Convocatoria Max-Planck-CONACyT 2017 para estancias postdoctorales" fellowship. This work is supported by a grant from the Max Planck Society. VED lovingly dedicates this work to her mother.

On-shell
Off-shell Primary Secondary Tertiary Euler-Lagrange eqs. In this appendix, we provide the explicit count of the degrees of freedom propagated by the Maxwell-Proca the-ory. We thus demonstrate the absence of ghosts claimed around equation (11) in Sec. II. Namely, given the La-grangian constructed in Sec. III, this count ensures that Definitions 1-3 in Sec. II hold true. Indeed, in sections A 1 and A 2 we obtain the results summarized in tables I and II, respectively.
For our present purposes, it suffices to consider this action in the schematic form Here, the {Q}'s are the M 0 = 4(N + M ) field variables that label the points in the configuration space C of the system, with I = 1, . . . , M 0 , and the dot denotes derivation with respect to timeQ I := ∂ 0 Q I . The velocities {Q} span the tangent space T C of the configuration space C. In all generality, we take C to be a differentiable Banach manifold. The Lagrangian L depends on the spatial derivatives of the field variables as well, L = L(∂ i Q I ). However, such dependencies are not relevant to the upcoming degree of freedom count and hence we omit them. Note that our L does not depend on time explicitly. This simplifies a bit the subsequent algebra, but it doesn't introduce nor remove any conceptual feature compared to the generic scenario L = L(t, Q I ,Q I ).
Our end point will be the well-known classical Hamiltonian formula that enumerates propagating degrees of freedom: where N 1 and N 2 stand for the number of first and second class constraints, respectively. Recall that, for physically sensible theories, N 2 is always an even integer. We shall find N 1 = 2N by construction and restrict the form of the Lagrangian through the demand N 2 ! = 2M , in order to obtain the desired count of degrees of freedom: two per Maxwell field and three per Proca field.
Needless to say, the formalism underlying this appendix is not new. It has been long established that any field theory, such as the Maxwell-Proca theory here studied, is always described by a singular Lagrangianfor example, see [41]. The equivalence between the Lagrangian and Hamiltonian pictures [42] implies that any field theory is always a constrained Hamiltonian system. Then, the famous Dirac-Bergmann (Hamiltonian) theory of constraints [43,44] -which was geometrized in [45]-, provides the canonical formalism to study such theories. In particular, it gives rise to our end point formula (A3).
Our appendix describes the fastest and easiest way to obtain the results of interest. We predominantly make use of the Lagrangian formalism for two fundamental reasons. First and foremost, this allows for the construction of the Maxwell-Proca theory. We anticipate that the uniquely defined primary constraint enforcing relations in (A7) are rewritten as equation (12), which is then systematically converted into as definite a Lagrangian as possible in section III A. Then, the also uniquely defined secondary constraint enforcing relations in (A20) are rewritten as equation (27). This latter set of equations is left as an additional requirement on the Lagrangian that must be satisfied to ensure ghost-freedom. Second, the applications of the Maxwell-Proca theory we propose in [1], as well as the cosmological applications that motivated the very notion of a generalized Proca field [14] use the Lagrangian (and not the Hamiltonian) of the theory as their starting point. We choose to make use of the Hamiltonian formalism only to explicitly count the degrees of freedom propagated by the Maxwell-Proca theory, thereby showing it is free of ghosts. The Hamiltonian itself is not worked out because this is not required for the counting process.
The fact that the Maxwell-Proca theory is a gauge theory introduces subtleties on both the Lagrangian and Hamiltonian side of the story. We will discuss and reference this important point where appropriate. For the benefit of the reader, the appendix is written as an upto-date and self-contained partial review of the theory of singular Lagrangians and constrained Hamiltonians in the presence of gauge symmetries.

Lagrangian constraint analysis
As already anticipated, we begin the study of the constraint structure of the Maxwell-Proca theory in the Lagrangian formalism. In this section, we obtain the results conveyed in table I. Namely, we count the number of independent Euler-Lagrange equations of motion, Bianchi identities and Lagrangian constraints present in the theory and classify them into primary, secondary and tertiary. We show that no tertiary Lagrangian constraints appear and hence no quaternary level arises. We also find the gauge identities in the theory.
What follows is not just a counting exercise, though. In the process, we will find the necessary and sufficient conditions to make our count match that of the auxiliary theory This is the theory of N Maxwell and M Proca free fields, where for simplicity we have taken the mass m = 0 of all the Proca fields to be the same. The aim is to show that the (self-)interactions between all N + M fields that we proposed in (19) do not introduce further unphysical degrees of freedom, compared to L aux . We will see this does not happen automatically: for M > 1, the constraint enforcing relation (7) must be fulfilled. Because we adopt a coordinate-dependent approach, the subsequent analysis depends on the form of the La-grangian: the number of identities and Lagrangian constraints varies with the field basis (unlike the number of degrees of freedom, which is basis-independent). The results here reported correspond to the mass eigenbasis formulation of the Maxwell-Proca theory, which was introduced in Sec III C.
Our analysis is largely based on the iterative algorithm presented in [46]. However, there exists a vast literature on the topic that could be consulted equivalently: [47], to mention but a few of the standard references.
Iterative analysis. Given the classical field theory Lagrangian density (A1) of interest, the principle of stationary action leads to the following Euler-Lagrange equations of motion: where we have defined In (A5), the right-hand side is an on-shell statement, while W IJ is a M 0 × M 0 matrix known as the Hessian. By construction 6 , our L is such that for all J and for I = 4z − 3, with z = 1, . . . , N + M . The above is the necessary and sufficient condition to ensure ghost-freedom of the theory at the primary level.
In the language of the main text, this condition reads (12). It was until now believed [21,27] that (A7) is the only condition needed to avoid ghosts. We will show that this expectation is too naive and that, for M > 1, (A7) must be supplemented by an additional condition (arising at the secondary level) in order to truly remove ghosts. Note that (A7) makes the rank of the Hessian be smaller than its dimension and so its determinant vanishes. This implies that the Lagrangian density under consideration is singular, as we said it should for a field theory. What is more, (A7) yields a Hessian of rank 3 (N + M ). This is what is needed to reproduce the primary level of (A4). The desired rank reduction of the Hessian can only happen by setting entire columns (or, equivalently, rows) of this matrix to zero for any first order Lagrangian, i.e. as in (A7). Rank reductions in the form of N + M columns being linear combinations of the remaining set of linearly independent 3(N + M ) columns can be absorbed into a field redefinition.
Given (A7), it is clear that not all the equations in (A5) are independent second order differential equations; only M 0 − M 1 number of them are. The remaining M 1 equations do not involve the accelerations {Q}. Instead, they are relations between the generalized field variables {Q} and their time derivatives {Q}. Here, Among these M 1 relations, we will further distinguish between Bianchi identities and primary Lagrangian constraints, whose number we denote by L 1 and K 1 , respectively.
In order to split (A5) into the above mentioned three subsets of equations, one must first obtain the M 1 linearly independent zero vectors {γ I } associated to the Hessian: These null vectors must be chosen so that they satisfy Additionally, one must also compute an auxiliary M 0 × M 0 matrix M IJ . This matrix always exists and is uniquely determined from the relations: Then, the M 0 − M 1 independent second order (primary) Euler-Lagrange equations are the non-zero subset in while the M 1 relations between the Q's and their velocities are ϕ a := γ a I α I = 0 ∀a.
Notice that the right-hand side is, in principle, an onshell equality. For the Maxwell-Proca theory, N number of the ϕ's are identically zero: These are the L 1 = N Bianchi identities associated to the Maxwell fields, which will be later related to first class constraints. It is important to realize that (A14) holds true off-shell: it does not require the right-hand side of (A13). In general, any identity at the primary level is offshell by definition. On the other hand, the imposition that the remaining ϕ's vanish, is an on-shell statement that gives rise to K 1 = M primary Lagrangian constraints on the Proca fields 7 . These will be later connected with second class constraints. It is worth pointing out that, in agreement with the standard notation, "≈" stands for a weak equality. In more detail, the requirement that those ϕ's which are not identities vanish defines a subspace C 0 ⊆ C in the configuration space, generally referred to as the constraint surface. A weak equality is one which only holds true in C 0 , not necessarily in the entire of C. We here finish the study of the primary level of the Maxwell-Proca theory. A significant remark on nomenclature follows. Equation (A14) is referred to as Bianchi identities because these are the natural generalization of the familiar Bianchi identity in classical electromagnetism to the full Maxwell-Proca theory. In the language of differential forms, one can succinctly write these full Maxwell-Proca Bianchi identities as dF −J = 0, where d denotes exterior derivative, F is the field strength of the usual Maxwell field and J is the most general possible source term within the defining axioms of the theory. That is, J consists of a set of sources coming from the non-linearities of the usual Maxwell field, from additional Maxwell fields and from all the Proca fields. In the well-known Maxwell electrodynamics limit (N = 1 and M = 0 with no interaction terms), J = 0 and the usual Bianchi identity dF = 0 is recovered. Before proceeding to the secondary level, we will briefly discuss the increased algorithm needed for reducible constrained systems. These are theories for which the onshell relations (A13) are not functionally independent. In our setup, we encounter this feature in two instances: i If working with a mixed eigenstate basis, obtained from the mass eigenbasis {A, B} by simple linear field redefinitions.
ii If studying the charged complex Proca proposal outlined in Sec IV.
For any reducible theory, (A13) satisfy L 1 relations of the form Z α a ϕ a ≡ 0 with Z α ∈ T C ∀α.
The above are the analogues of the Bianchi identities (A14) of the irreducible case. The case i is trivial, in the sense that (A16) amount to linear relations between the reducible constraints {ϕ}. On the other hand, for the case ii, the Z's include differential operators and can be determined iteratively order by order 8 , as explained in [30]. Once that (A16) has been established, only the M 1 −L 1 irreducible Lagrangian constraints {ϕ} analogue to (A15) are to be used to continue with the iterative algorithm.
Going back to the Maxwell-Proca theory, the demand that (A15) holds true at all times (i.e.φ α ≈ 0 for all α) initiates the study of the secondary level 9 , which is an innately on-shell level. Specifically, it gives rise to K 1 = M secondary equations of motion: where we have defined in close analogy to the primary quantities in (A6). As happened with (A5) before, generally not all the above secondary equations of motion are independent second order differential equations. Instead, they can again be divided into three subsets: K 1 − M 2 secondary Euler-Lagrange equations, L 2 identities and K 2 secondary Lagrangian constraints, with In order for the Maxwell-Proca theory to have the desired Lagrangian conformation (i.e. that of the auxiliary theory in (A4)), M 2 = K 1 is needed, implying that there are no independent secondary Euler-Lagrange equations. It is easy to see that this is uniquely achieved by setting which is the sufficient and necessary condition that ensures ghost-freedom at the secondary level. Fulfillment of (A20) implies the constraint enforcing relations (27) on the Maxwell-Proca interaction terms. Written in this manner, it is easy to see that the above identical vanishing happens by construction for any number of Maxwell fields and when a single Proca field is considered. In other words, (A20) trivially holds true if M ≤ 1. However, when two or more Proca fields are present in the theory, (A20) -equivalently, (27)-non-trivially restricts the proposed interaction terms. We highlight that this is precisely the anticipated condition that was missed out in all [21,27] and that supplements the known condition (A7). To avoid the propagation of ghosts in any multi-Proca setup, both (A7) and (A20) must hold true. As noted in the main text, we regard this realization as our main result. Given (A20), repetition of the logic exposed for the primary level requires us to obtain the M 2 null vectors { γ α } 10 associated to W αβ , in such a way that they satisfy the orthonormality relation γ r α γ α s = δ r s . Further, it is easy to see that the always extant M 1 × M 1 auxiliary matrix M αβ , which is analogue to (A11) at the secondary level and which is uniquely determined from identically vanishes in this case: M αβ ≡ 0. Subsequently, the M 2 secondary relations defined in C 0 among the generalized coordinates and their time derivatives are given by 11 We see that none of the φ's vanish identically, implying L 2 = 0. To obtain the K 2 = M 2 secondary Lagrangian constraints in the theory, one must indeed demand that they are zero: We will soon see that these relate to second class constraints for the Proca fields. For clarity, we have here introduced the symbol " 1 ≈" to indicate that the requirement that the φ's vanish defines a subspace C 1 ⊆ C 0 in the constraint surface. As a result, a weak equality that only holds true in C 1 (and not necessarily in C 0 or the entire of C) shall be denoted by 1 ≈ in the following. One more iteration is required. Again, all statements at this tertiary level hold on-shell only. Demanding that (A24) holds true at all times (namelyφ r ≈ 0 for all r), M 2 = M number of tertiary equations of motion are obtained: with W rs and α r , the tertiary counterparts of (A18) before, defined as W rs := γ r α γ α I ∂φ s ∂Q I (A26) and respectively. The M 2 equations (A25) can generically be divided into M 2 − M 3 independent tertiary Euler-Lagrange equations, L 3 (≤ M 3 ) identities and M 3 − L 3 tertiary Lagrangian constraints, where However, M 3 = 0 for the Maxwell-Proca theory. This implies that all equations in (A25) are independent tertiary Euler-Lagrange equations. It also means L 3 = 0 = K 3 , which signals the termination of the iterative procedure: an on-shell quaternary level would unfold from guaranteeing the stability of the K 3 tertiary Lagrangian constraints, but there are none. Gauge identity analysis. We turn our attention to gauge identities, which are known to escape the above algorithm and require separate consideration. We shall use g to denote their number. The key point is to notice that the subset of N number of α's defined in (A6) stemming from the Maxwell fields vanish identically: This is obvious from their very definition, once one takes into account that our theory depends on the Maxwell fields only through their field strengths, as pointed out in Sec. III C. This is clearly an off-shell identity, since it is independent of (A5). Matter of fact, (A29) are precisely the g = N gauge identities we were looking for. We will soon link them with first class constraints. Notice that the Bianchies in (A14) are projections of the gauges in (A29) to the kernel of the Hessian W IJ defined in (A6). In this way, the Bianchi identities can be viewed as relics of the gauge invariance of the theory. Accordingly, we will shortly relate them to first class constraints as well 12 .
We have now fulfilled our initial promise and completed the derivation of table I.

Hamiltonian constraint analysis
Next, we focus on the study of the Hamiltonian constraints in the Maxwell-Proca theory. In this section, we derive the results summarized in table II. In more detail, we first relate the identities and Lagrangian constraints we found in the previous section to the primary and secondary Hamiltonian constraints in the Dirac-Bergmann formalism [43,44]. We loosely refer to these as Dirac primaries and secondaries, respectively. Our analysis is primarily based on [46] again. Then, we reclassify these Dirac constraints as first or second class. To this aim, we stick to [48]. We restate, however, that an extensive bibliography exists elaborating on this topic -for example, see [49] and references therein. In this way, we obtain the number of propagating degrees of freedom in the theory, according to (A3). As already mentioned, the Hamiltonian formalism is an equivalent reformulation of the Lagrangian formalism in section A 1. Its central object of study is the Hamiltonian density H, whose explicit form we do not work out here (as it is not required to count the propagating degrees of freedom). The interested reader can consult in e.g. [50] the procedure to determine the H following from (A1). Formally, H relates to the Lagrangian density L as H :=Q I P I − L(Q I ,Q I ).

First class Second class
(A30) Generically, H is a function of three types of variables: the generalized coordinates Q I introduced in (A2), their time derivativesQ I and the momenta P I . The {P }'s span the cotangent space T * C of the configuration space C. For clarity, we note that P I encodes both the momenta P We denote as R the (3M 0 )-dimensional space comprising C, T C and T * C. The Hamiltonian density is independent of theQ's when restricted to the subspace R 0 ⊆ R defined by 13 In the case of the Maxwell-Proca theory of our interest, only M 1 = N +M of the above relations are independent. These are the Dirac primary constraints Ψ a . We shall only be concerned with their determination modulo the soon to be introduced Poisson bracket operation, so we define them as brackets of the constraints in the theory: This is a 2M 1 × M 1 matrix and, for the Maxwell-Proca theory, the rank equals the number of Proca fields under consideration: R = M . Note that R < M 1 , which reflects the singular nature of the theory in the Dirac-Bergmann approach. Then, it is possible to define two sets of linear combinations of the Dirac primaries { Ψ, Ψ}, the first set containing R elements and the second set containing S := M 1 − R = N elements: where we have chosen to reuse the indices introduced in (A14) and (A15) for simplicity 15 . Here, the parameters v α a , v α a ∈ R can and should be chosen so that the newly defined Ψ's and Ψ's are linearly independent and so that the Ψ's are first class (primary) constraints: {Ψ a , Ψ α } ≈ 0 ≈ {χ r , Ψ α } ∀a, α, r.
Clearly, the Ψ's are second class (primary) constraints. The procedure to distinguish between first and second class constraints among the Dirac secondaries is similar, albeit a little more cumbersome algebraically. First, we define the following set of linear combinations of all the Dirac secondaries and the second class Dirac primary constraints 16 : where the parameters J a b can and should be chosen so that J is a non-singular matrix: det(J) = 0. Then, in analogy to (A40), it is possible to form two sets of linear combinations of the χ ′ 's: Once again, the reuse of the (α, α) is due to the criterion explained in footnote 7. We denote by m 2 (≤ M 2 ) the 15 Again, their correspondence to the field labels in the main text is not coincidental: ( Ψα, Ψ α ) are constraints on the (Proca, Maxwell) fields. 16 It is possible to trivially extend the construction to include all Dirac primaries, but this is not necessary for our purposes.
number of χ's. Here, w α a , w α a ∈ R − {0} should be chosen so that the χ's are all first class {Ψ a , χ α } ≈ 0 ≈ {χ a , χ α } ∀a, α, for a maximal number m 2 . Then, the M 2 − m 2 Dirac secondaries χ's are second class constraints. For the Maxwell-Proca theory of our interest, m 2 = N . This completes the derivation of table II. As a double-check, we verify that (A46), the matrix ∆ of the Poisson brackets of all the constraints that we identified as second class, is non-singular: det∆ = 0. Note that ∆ is always an antisymmetric square matrix of dimension R + M 2 − m 2 (= 2M , for the Maxwell-Proca theory).
For completeness, we indicate the relation between the first and second class constraints and the Lagrangian picture next. As shown in the original works [43,44], linear combinations of the second class constraints ( Ψ, χ) generate the Lagrangian primary (A15) and secondary (A24) constraints. This happens through their time derivatives and within the suitable Hamiltonian subspaces, in agreement with (A34)-(A37). Observe that these constraints are all associated to the Proca fields. Similarly, linear combinations of the first class constraints ( Ψ, χ) generate the Bianchi identities (A14), all of which stem from the Maxwell fields. Surprisingly enough, the connection between first class constraints and gauge identities has been a subject of great controversy for decades and continues to receive attention. Although the correct relation was long ago suggested [44], this was for a long time believed to be a matter of interpretation with no physical implication. The subject was again pondered over some thirty years later, for instance see [51], but remained as a majoritarily misconceived issue. A lucid review of the origin and resolution of the main puzzlements to this respect appeared in [52] (see also references within). At last, the topic was shown to be of physical transcendence for the case of classical electrodynamics in [53]. The interested reader can consult this same reference for a meticulous historical review. In conclusion, only a suitably tuned linear combination between a primary and secondary first class constraints is able to generate a gauge symmetry consistently. For the Maxwell-Proca theory, this means particular compositions of ( Ψ, χ) yield (A29).
Finally, adding the number of first (second) class constraints coming from the Dirac primaries and secondaries, we see that N 1 = 2N (N 2 = 2M ), as anticipated. Using (A3), the number of propagating degrees of freedom for the proposed Maxwell-Proca theory is that of the auxiliary theory (A4): Namely, each Maxwell field propagates two degrees of freedom and each Proca field propagates three degrees of freedom. In other words, the interaction terms between the fields in Eq. (19) supplemented by the condition (27) are well-behaved, in the sense that they do not give rise to the propagation of ghost-like degrees of freedom.
As a final remark, note that we have worked in four spacetime dimensions for definiteness. However, all the results in this appendix can be straightforwardly generalized to arbitrary dimension d. The only change would be in the dimension of the configuration space C: M 0 = d(N + M ). This then would yield a degree of freedom count n dof = (d − 2)N + (d − 1)M .