The order $p^8$ mesonic chiral Lagrangian

We derive the chiral Lagrangian at next-to-next-to-next-to-leading order (NNNLO) for a general number $N_f$ of light quark flavours as well as for $N_f=2,3$. We enumerate the contact terms separately. We also discuss the cases where some of the external fields are not included. An example of a choice of Lagrangian is given in the supplementary material.


Introduction
Chiral perturbation theory (ChPT) [1,2,3] is an effective field theory (EFT) of QCD at low energies. An introduction and further references can be found in [4]. For N f flavours of light quarks the effective degrees of freedom are the N 2 f −1 lightest pseudoscalar mesons. These arise as the (pseudo-)Goldstone bosons of the spontaneously broken chiral symmetry SU(N f ) L × SU(N f ) R to SU(N f ) V . For N f = 2, i.e. when only u and d quarks are considered, these are the pions, and for N f = 3 when the s quark is included, the additional degrees of freedom are the kaons and the eta meson. The power counting is done in terms of a generic momentum, p, and, as in any EFT, a calculation is done to a fixed order in this variable. This requires a chiral Lagrangian of the corresponding order. The leading order (LO), or p 2 , next-to-leading order (NLO), or p 4 , and next-to-next-toleading order (NNLO), or p 6 , chiral Lagrangians are known. The LO is basically current algebra, the NLO was derived in [2,3]. The p 4 Lagrangian in the anomalous sector contains no new free parameters [5,6]. The NNLO or p 6 Lagrangian was studied in [7,8] and for the anomalous sector in [9,10]. In addition the divergence structure is known for all of these [2,3,9,11].
The basis of [8] for a general number of flavours and N f = 3 is generally accepted as minimal. For two flavours, an extra relation was found in [12] and another one for the case of vanishing scalar and pseudoscalar external fields [13]. In general, it is quite difficult to be sure one has a minimal Lagrangian. The only foolproof method we are aware of is to find Green functions of the external fields that allow to determine all free parameters. This program has been done to NNLO for the case with vanishing scalar and pseudoscalar external fields in [14]. A very systematic approach to NNLO was done in the Diplom thesis [15]. This latter reference was the inspiration for the present work.
In this paper we present a derivation of the next-to-next-to-next-to-leading order (NNNLO), or p 8 , chiral Lagrangian L 8 in the nonanomalous or positive intrinsic parity sector. We do this for two, three and a general number of light quark flavours. We treat the cases of only Goldstone bosons, including scalar and pseudoscalar external fields, including vector and axial vector external fields and including all four sets of external fields. As a byproduct we clarify the number of terms in the two-flavour case at NNLO.
There are several underlying reasons to do this work. It is of general interest to know how many terms there actually are at a given order in the EFT expansions. For the case of extensions of the standard model there is a large literature, we quote here only the Hilbert series approach [16] and a software to find the minimal Standard Model EFT Lagrangian [17]. The remaining literature can be traced via these two. Our approach can in principle be extended to these cases as well. The second reason is that at tree level there are a number of works studying effective fields theories using amplitude methods [18] and references therein. Knowing the solution to p 8 provides a cross-check for extensions of that work to higher derivatives. The third motivation to find the Lagrangian is that in a previous calculation [19], analytical expressions for the pion mass and decay constant at order p 8 were obtained for two flavour ChPT. However, as the NNNLO Lagrangian was not yet known, the tree level contributions were parametrised in terms of two unknown parameters. These parameters are in fact combinations of low energy constants (LECs) multiplying the monomials in the minimal set we derive here, and their numerical values can in principle be found from experiment. The present status of their values was reviewed in [20] and their determination from lattice QCD in [21].
The Lagrangians are of the generic form The c (2n) i are the free constants, usually referred to as low-energy-constants (LECs). The O i are monomials in terms of the (pseudo-)Goldstone bosons and external fields that are invariant under chiral symmetry, hermitian and invariant under charge conjugation and parity 1 . The set of O i should be complete and minimal.
We first write down all possible single term operators 2 in terms of a number of basic building blocks. At this level we impose parity (P ) and invariance under chiral symmetry. This is described in section 2.1. This basis is the simplest one to derive all the relations described below. At this stage we also take into account the different ways to write each operator due to cyclicity of the trace, relabelling of indices and those that are zero because of contraction of antisymmetric and symmetric pairs. The next step is then to combine these operators into monomials that are invariant under charge conjugation (C) and hermitian conjugation (H). The method we use here is discussed in section 2.2. An extension was needed to find all possible contact terms, i.e. the operators only depending on external fields. This is discussed in section 2.3.
In general this leads to a very large number of operators, but here the importance is to be complete. It is also easier to make sure we have obtained all relations when the basis at this level contains all possible single term operators, even those that can obviously be removed. The linear relations between operators that we use are all partial integrations, finding operators that can be removed using field redefinitions by using the equation of motion, the Bianchi identity and relations that follow from commuting derivatives. These are discussed in sections 3. 1-3.4 and are what we use to determine the general N f flavour Lagrangian. For two-and three flavours there are additional relations following from the Cayley-Hamilton relation. These are discussed in section 3.5. Again, there are a number of extra conditions needed when determining the relations between possible contact terms. These are discussed in section 3.6.
The analytical work is done using FORM [22] and Python to rewrite some of FORM's output back into FORM commands. Identical relations are removed within the FORM programs but we still end up with a very large number of operators and relations. The number of independent relations was determined by using Gaussian elimination. This we implemented using C++ with exact arithmetic via the GNU multiple precision library (GMP) [23].
The results are discussed in section 4. We give here the numbers of independent monomials for the various cases of number of flavours and inclusion of external fields. Explicit lists of the monomials involved are relegated to the appendix for the contact terms and to the supplementary material for the others. As an example of the use of our Lagrangians we discuss meson-meson scattering in section 5 and show that the number of LECs agrees with general considerations on the form of the amplitude.

Building blocks for chiral Lagrangians
ChPT is constructed from the chiral symmetry G = SU(N f ) L × SU(N f ) R of the QCD Lagrangian with N f massless quarks. Due to a non-vanishing quark vacuum expectation value, this symmetry is spontaneously broken to H = SU(N f ) V and the corresponding N 2 f − 1 Goldstone bosons live in the coset space G/H and are the lightest pseudoscalar mesons.
Gasser and Leutwyler [2,3] added scalar s, pseudoscalar p, vector v µ and axial-vector a µ external fields to the QCD Lagrangian. They transform in such a way that the chiral group G can be made local and leave the QCD Lagrangian invariant. q is an N f column vector, s, p, v µ , a µ are N f × N f matrices in flavour space. We assume that v µ and a µ are traceless. The Goldstone boson fields φ can be organized in a unitary matrix u(φ) [24,25]. The transformation under a chiral symmetry transformation ( The top line in (3) is the definition of h(g L , g R , u(φ)) which is an element of the unbroken subgroup H. The remaining lines in (3) define the combinations χ, ℓ µ and r µ . In the remainder we drop the arguments of u and h. The constant B in the definition of χ is related to the chiral limit value of the quark vacuum expectation value and pion decay constant B = − qq /N f F 2 . We also use the notation X = tr (X) for any matrix in flavour space. The effective Lagrangian is expanded in powers of p 2 according to where the first four contributions are the LO, NLO, NNLO and NNNLO Lagrangians, respectively. Each of these is a sum of order p 2n monomials O where the c (2n) i are LECs and N 2n is the number of monomials. We use here the standard counting where the fields s, p are counted as p 2 , v µ , a µ as order p. The monomials O (2n) i depend on the field content, and each is invariant under chiral symmetry, Lorentz transformations and the discrete transformations parity, charge conjugation and hermitian conjugation.
The easiest way to construct invariants is to produce objects that transform linearly under a simple symmetry group. Here one can use different choices, all transforming purely left or purely right-handed or purely under the vector group SU(3) V . The former choices were made in [7,10,26], the latter in [8,9]. We will use the latter except for the determination of the contact terms. Table 1: Transformation properties of the first chiral building blocks. ε(0) = −ε(i = 1, 2, 3) = 1.

Chiral building blocks
The p 4 Lagrangian was constructed in [3] using We will refer to F µν L,R as field strengths. We can define covariant derivatives D µ acting on these building blocks that transform as the building blocks themselves via for any operator O transforming in one of the four ways. The respective momentum orders of the above structures are U ∼ p 0 , D µ U ∼ p and χ, F µν L,R ∼ p 2 . Their transformation properties under the discrete transformations can be found in table 1. This set of building blocks we will refer to as the first basis.
As mentioned above it is easier to use building blocks that all transform the same way. We choose to use the transformation under h ∈ H as defined in (3) and all building blocks X transforming as X −→ hXh † . It can be checked that the following blocks transform in the desired way and form a complete set: Table 2: Transformation properties of the main chiral building blocks, where ε(0) = −ε(i = 1, 2, 3) = 1.
Note that when no external fields are included, only u µ with ℓ µ = r µ = 0 needs to be considered. The covariant derivative ∇ µ is defined as where the chiral connection is Equations (7) and (9) highlight one of the reasons this set of building blocks is more convenient than the first set. We will refer to this as the main set in the remainder.
Since ∇ µ X transforms like X under chiral symmetry, we can always take covariant derivatives of the structures in (8) and obtain another building block of one order higher in p-counting. The respective orders of the above structures are u µ , ∇ µ ∼ p and χ ± , f ±µν ∼ p 2 . Their transformation properties under P , C and hermitian conjugation are given in table 2. Further note that u µ and f ±µν are traceless since r µ and ℓ µ are.
All possible operators of a given order in p that are invariant under the chiral group can be generated by taking traces of the building blocks and their covariant derivatives. At this stage we make sure that terms that are the same from cyclicity of traces or relabelling indices are identified. We also remove terms that are obviously zero because a pair of antisymmetric indices is contracted with a pair of symmetric ones.
We look here at the nonanomalous sector so we do not need the Levi-Civita tensor ε µναβ . Lorentz-indices are always contracted pairwise so parity can be easily implemented as well by keeping only terms with an even number of u µ , χ − and f −µν .

A basis of odd and even monomials
In addition to chiral invariance and parity as implemented in the single term operators built directly from the building blocks we need to impose charge conjugation and have terms that are hermitian. One advantage of the building blocks is that charge conjugation and hermitian conjugation always relate the same operators. This can be seen from where X i is any of the building blocks in (8) or their covariant derivatives.
We can now construct the combinations that are hermitian and odd or even under C. Any operator transforms according to where O j can be any operator in the basis and λ C,h.c. ± = ±1 are the eigenvalues under the respective discrete transformations. There are thus four possibilities of the signs, which we denote by (λ C ± , λ h.c. ± ), in addition we can have i = j or i = j. We define C-even and C-odd monomials O ± i for the case j = i: For j = i we instead have The final Lagrangian should only contain the monomials O + i .

Contact terms
The monomials constructed using the main building blocks in (8) are sufficient to construct a complete basis and determine the total number of terms in the Lagrangian. However, there are in general terms possible that only depend on the external fields, i.e. the contact terms. These cannot be simply measured in physical processes but depend on the precise definitions of the external fields used in QCD. For determining the number of physically relevant parameters at a given order it is therefore desirable to also know the number of contact terms.
The contact terms only depend on external fields and are therefore best handled in the first basis using O = χ, F µν L , F µν R and their covariant derivatives. The relations between the building blocks in this and the main basis are This shows that the contact terms are included when using the main basis.
The contact terms are found by writing down all possible chiral invariant trace structures of the building blocks O and their covariant derivatives, taking into account that the transformations under the left and right group should match. One also needs to be more careful with defining the P, C and hermitian conjugation even operators generalizing the method of the previous subsection.
The third complication is that there is for a given number of flavours another combination with a well-defined chiral transformation [27]: This allows to take care of the invariant operator det χ + det χ † of chiral order p 2N f and the extra building block for the two-flavour case in [2], calledχ there. χ is perfectly regular when s, p go to zero. This can be seen from the explicit expressions using χ ij = x ij : Terms involvingχ can be rewritten in the other building blocks so they are included in our main basis. However, structures involvingχ are not fully independent of structures with χ, so a careful check that they are really independent is needed. An example is χχ † = χχ † for N f = 2. The easiest way to check is simply to use the explicit forms given in (17).

Constructing the chiral Lagrangian
Equipped with the chiral building blocks of the main set basis it is now possible to construct the Lagrangians. By demanding that the operators be invariant under chiral symmetry, Lorentz transformations and the discrete symmetries of QCD, the only possible monomials at LO are This gives the LO Lagrangian in the standard form F is a LEC corresponding to the LO pion decay constant and the second LEC B is hidden inside χ as given in (3). At higher orders many more structures must be included, and it is important to note that the operators O i can be linearly related. The relations relevant for the NNNLO Lagrangian will be discussed in more detail below, but in short they arise from • the vanishing of total derivatives (i.e. from integration by parts), • terms that can be removed using field redefinitions, this is equivalent to removing combinations of operators that vanish due to the LO equation of motion (EoM), for a proof see [28] and appendix A of [8], • the Bianchi identity of the field-strength tensor • the Cayley-Hamilton theorem (for N f = 2, 3).
The Schouten identity [29] resulting from the fact that any completely anti-symmetric tensor of rank higher than the number of dimensions must vanish, yields no extra relations here since Lorentz invariance in the nonanomalous sector implies that there can be at most four independent indices. We first construct all operators using the building blocks of the main set defined in section 2.1. Then we determine all relations on this set of operators since it is easier to handle the single term operators in the relations. We then rewrite the operators in terms of the C-even and odd operators O ± i defined in section 2.2 and keep the C-even part of all relations which involve only the N 2n monomials O + i . We remove identical relations and write the remainder in matrix form A ij O + j = 0. The number of linearly independent equations is thus equal to the rank of A ij , and the minimal set of monomials at order p 2n therefore has N 2n − rank(A) elements.
Note that also the contact terms are linearly related, so that a similar procedure must be done in order to find the minimal number of such terms. The differences are discussed in section 3.6.
We used FORM [22] to produce the original operator basis and to produce all linear relations. The output of this was then used as input for a C++ program calculating the rank of A ij by Gaussian elimination. In order to avoid numerical imprecision, we used the library GMP [23] which allows for exact arithmetic. Given the size of the final matrices, we implemented a method using sparse matrices allowing the Gaussian elimination to work without swapping. The general N f case was handled by choosing a number of values for N f and checking that we got the same result, we used N f = 17, 49, 73, 199 and obtained the same basis for all of them.
As a check, we reproduce the known numbers for NLO and NNLO in [2,3,8]. We have also two independent implementations of the procedures as a check. A minor check is that the relations can be rewritten separately in C-even and C-odd operators.

Partial integration or addition of total derivatives
The order p 2n action is invariant under the addition of a total derivative to the Lagrangian. This leads to partial-integration relations. It is very difficult to be sure one has all such relations using partial integrations. However all total derivatives can be constructed in a fashion similar to the construction of all basis operators. We construct all operators of the form O µ that are chiral and parity invariant. We can then add a relation ∂ µ O µ = 0 since we can add a term ∂ µ O µ to the action.
One therefore has to write down all possible structures of order p 2n−1 and then take the derivative of each these. Generating all possible structures O µ can be implemented in the same way as generating all possible operators O i thus guaranteeing that we have all possible partial-integration relations.
As an example at NNNLO we get the relation so that one of the operators can be written in terms of the other two.

Field redefinitions or equations of motion
Field redefintions also allow to remove terms from the Lagrangian. This is at the level of the terms that can be removed equivalent to deriving relations between the higher order terms using the LO EoM [28,8].
The LO Lagrangian in (19) has EoM The way we implement this is to take the set of all possible operators generated using the building blocks. We then look for all occurrences of ∇ µ u µ or covariant derivatives thereof and then use (21) to produce a relation. An example is the operator χ + u ρ χ + ∇ ρ ∇ µ u µ which using the EoM leads to the relation We replace all possible occurrences of ∇ µ u µ but only one at a time if it occurs twice in an operator. Using it twice in the same term does not yield new relations. Since we have constructed all possible operators this method catches all possible relations from field redefinitions.

The Bianchi identity
The field strength tensor Γ µν is defined through the equation and can be expressed in our building blocks via Γ µν satisfies the Bianchi identity The Bianchi identity is of order p 3 . Using the definition of Γ µν this can be written We will use it in this form which makes it clear that in the absence of external fields the Bianchi identity does not give extra relations. B µνρ is cyclic in µνρ and transforms as our main set of building blocks. It is also automatically zero under a contraction of any two indices. Based on chiral dimension, the only structure the Bianchi identity can be traced with at NLO is u µ , but as each term in the Bianchi identity is of even parity (this is easily seen from table 2) such an operator would be odd under parity. Therefore, the Bianchi identity becomes non-trivial only beyond NLO. At NNLO, it can be traced with ∇ µ f +νρ , f −µν u ρ or ∇ µ u ν u ρ [8]. At NNNLO there are additional complications, since it is possible to have several traces and that it can be traced with tensors of rank three or five. However, we can treat B µνρ as a separate building block and use the same methods as before to construct all possible operators involving B µνρ . Each of these gives after inserting (26) a relation between the operators O i .
It is sufficient to have one insertion of B µνρ since a relation involving it twice can be constructed from those only involving it once. In addition since we have all possible partial-integration relations there is also no need to consider covariant derivatives of B µνρ .

Other relations for general N f
There are also other identities that can be used to obtain relations. These essentially follow from the fact that partial derivatives commute or that the commutator of partial derivatives is related to field strengths.
The first is This implies that ∇ µ u ν is symmetric in its two indices when neglecting external fields. This is implemented in a very similar way as done for the EoM. We look for all occurrences of ∇ µ u ν and covariant derivatives of it in our complete list of operators and replace it with the relation (27). The second such equation is the defining relation of the field-strength tensor in (23) and its expression in terms of the building blocks (24). We implement all these relations by looking inside the list of operators for quantities that have two covariant derivatives acting on them. We then replace ∇ µ ∇ ν X by (28) and replace Γ µν by its expression (24). This leads to a number of relations between the original operators. We also have cases where more than two covariant derivatives act on an object X. One could in principle also use the relation then on two "inside" covariant derivatives. This does not lead to new relations since we have included all possible operators and already have all partial-integration relations. The relation produced by "inside" covariant derivatives will be produced by the relation from the term where the "outside" derivatives have been partially integrated first.

The Cayley-Hamilton theorem
The characteristic polynomial of any N f × N f matrix A can be written where λ is an eigenvalue of A and I is the unit matrix. This can be written and when expanded in 1/λ this yields the respective polynomial coefficients in p(λ). The Cayley-Hamilton theorem states that any N f × N f matrix satisfies its own characteristic equation, i.e. p(A) = 0 as a matrix equation. For N f = 2 and N f = 3 one finds For a given N f the matrix A can be split into a sum of N f terms, and the Cayley-Hamilton theorem thus gives a relation between them. Letting A = B + C for N f = 2 and using that both C and D also satisfy the Cayley-Hamilton theorem yields For N f = 3 and A = B + C + D the relation is Tracing each of the above relations with any other matrix of the corresponding size thus produces a relation between the monomials. In practice how we use these relations is that we look for all possible operators constructed using our building blocks and for each possible trace we need to look for all possibilities of having blocks BC. B and C can be a single building block but also products of them. Once we have a choice, we replace BC with the relation (32) giving a relation between operators. The three-flavour case is dealt with in an analoguous manner but looking for a possible block BCD inside traces and replacing it with (33).

Contact terms
As noted in section 2.3 here we have to work in the first basis. We classify the possible operators into three classes: Terms 1. with only χ, χ † ,χ andχ † and covariant derivatives, 2. containing both χ,χ and field strengths.
3. with only F µν L and F µν R and covariant derivatives, For the first two classes there were sufficiently few terms that it could be handled by hand. The last class can be dealt with in the same fashion as before, but with some minor modifications. Partial integrations are immediately generalizable. Field redefinitions are not relevant for the contact terms since the external fields are given. The Bianchi identity is valid separately for the left and right handed field strengths 3 The relations from commuting covariant derivatives split up according to the four cases in (7) depending on the transformation property of the building block: The Cayley-Hamilton relations remain the same, however one must make sure that the (products of) building blocks that form B, C, D in (32) or (33) transform either all under SU(N f = 2, 3) L or all under SU(N f = 2, 3) R only. In practice this means the Cayley-Hamilton relation is only useful for classes 2 and 3.
We wrote a separate set of programs to determine the number of independent contact terms of class 3 and a final set to remove terms that can be removed when the minimal set of contact terms is included.

The NNNLO Lagrangian
In this section we present the NNNLO Lagrangian obtained from the considerations in the previous sections. The number of independent monomials and the number of contact terms is unique.
The original papers constructing the different Langrangians made a number of choices of which monomials to keep and which to consider superfluous. However, the actual choice of the independent set of monomials is very much a matter of choice.
We follow here a number of guiding principles to keep terms: 1. We always keep the maximal number of independent contact terms.
2. We remove terms that vanish when external fields vanish as much as possible.
3. We remove terms involving covariant derivatives in favour of those involving external fields.
4. We remove terms that contribute to processes with a low number of mesons as much as possible. This is done by counting occurrences of u µ , χ − and f −µν .
5. Terms involving scalar-pseudoscalar external fields are placed before those with only vector-axial-vector external fields.
6. After that we preferentially keep terms with lower number of flavour traces. This is to make the large N c counting [30] of the monomials explicit, only leading in N c is equivalent to keeping only single trace monomials.
In addition in [8] covariant derivatives were always symmetrized and χ ±µ was used rather than ∇ µ χ ± . In the lists of terms given we can always choose instead of the monomial listed those with symmetrized covariant derivatives (including the index µ of u µ ) and replace ∇ µ χ ± by χ ±µ . These are always equivalent choices of monomials. The cases given below always include the (pseudo-)Goldstone bosons but have different sets of external fields included. The cases with only external fields correspond to the columns giving the number of contact terms.
Let us first look at the case with all external fields included. The number of independent monomials we find at each order and the number of those that are contact terms are given   Table 4: Number of monomials in the minimal basis for the case with no scalar or pseudoscalar external fields included. Also listed is how many of them are contact terms.
for general N f flavours and for N f = 2, 3 in table 3. We agree with the known results for order p 2 , p 4 , p 6 with following caveat: in [15] there is one more relation listed for the p 6 N f = 2 case, however if correctly rewritten in terms of C-even operators this relation is equivalent to an earlier one. So after correcting that we agree there as well. The monomials our program produces at order p 6 are slightly different from those in [8]. We find a basis that has fewer terms with only the Goldstone bosons. We have however checked that the basis given in [8] is a minimal basis for us as well, i.e. the 115, 94, 55 [12] monomials translated into ours form a complete set.
As an example of the size of the calculations involved we give the intermediate numbers for the general N f case. The main basis consists of 9740 C-even monomials and we find 12444 different relations. The Gaussian elimination gave 7878 linearly independent relations, and thus a final minimal basis of 1862 monomials.
The second case that we study is when we have no scalar or pseudoscalar external fields and is given in table 4. This case was studied at order p 6 in [14] and in an AdS/QCD context in [13]. The number of operators at orders p 2 , p 4 can be easily checked against the known full Lagrangians. The p 6 case for N f = 2 agrees with [14].
The third case is when we only include scalar and pseudoscalar external fields. This is given in table 5. The results for p 2 and p 4 can be easily determined from the known full   Lagrangians. The final case is when we do not include any external fields. Here there are obviously no contact terms. The number of monomials at each order is given in table 6. The number of monomials at p 2 , p 4 can again be easily seen from the known Lagrangians. In [8] the monomials were chosen to follow roughly the same criteria as given above, however this was done by hand and no check whether the basis without external fields was minimal was done. From table 2 in [8] one sees 21 (N f ), 12 (N f = 3) and 6 (N f = 2) monomials that do not vanish when all external fields are put to zero at order p 6 . The extra relation found in [12] reduces the N f = 2-number to 5 in agreement with what we obtain. One can go back to the unpublished notes underlying [8] or use the list of monomials and relations given explicitly in [15] removing all external fields to derive that for general N f at p 6 we have 19 independent monomials, in agreement with our result.
The number of terms is in general rather large and finding relations between processes in general will be quite some work. To order p 6 this was done in [31]. However we can look at the types of terms that appear for the case with no external fields. Our ordering scheme removes terms preferentially with a lower minimal number of mesons involved. We can thus check how many monomials require at least 4, 6 or 8 mesons in the vertices to contribute. This is done by explicitly looking at the minimal Lagrangians our programs produce. The result is given in table 7. Note that we have a slightly lower number of monomials contributing to terms with only 4 mesons than the monomials quoted in [8]  of LECs contributing to 4-meson processes is compared with a general amplitude analysis in section 5.

Meson-meson scattering
The general amplitude M(s, t, of meson-meson scattering for a general N f can be written as [32,33] M(s, t, u) = X a X b X c X d + X a X d X c X b B(s, t, u) + δ ab δ cd C(s, t, u) + δ ac δ bd C(t, u, s) + δ ad δ bc C(u, s, t) , where s, t and u are the Mandelstam variables satisfying s + t + u = 0 in the massles case considered here. The X a are SU(N f ) generators normalized to 1. The functions satisfy C(s, t, u) = C(s, u, t) and B(s, t, u) = B(u, t, s). This follows from crossing and the SU(N f ) V symmetry. For pion-pion scattering and N f = 2 the amplitude can be written M(s, t, u) = δ ab δ cd A(s, t, u) + δ ac δ bd A(t, u, s) + δ ad δ bc A(u, s, t) , where A(s, t, u) = A(s, u, t). The traces with four generators can here be written in terms of the Kronecker delta terms.
For N f = 3 the Cayley-Hamilton relation (33) shows that the sum of all 6 four generator traces can be rewritten in terms of the Kronecker delta terms. The part of B(s, t, u) that is fully symmetric can be transferred to the C(s, t, u) using that relation.
The tree level contributions from our Lagrangian become polynomials in s, t, u, so we can expand the two functions as Chiral invariance requires α 0 = β 0 = 0. For N f = 2 we only have A(s, t, u) which can be expanded as C(s, t, u). The number of free parameters at order p 2n matches exactly that in the corresponding column of table 7. For general N f the number of parameters doubles w.r.t. N f = 2 since we now have both functions C(s, t, u) and B(s, t, u) present. This again matches the numbers in table 7.
The N f = 3 case is a bit more subtle, we can transfer the fully symmetric part in s, t, u away from B(s, t, u) so we have to determine how many of these there are at each power taking into account s + t + u = 0. At orders p 4 , p 6 , p 8 there is only one such combination at each order. This agrees with the numbers in the N f = 3 column in table 7. If one only considers meson-meson scattering no results follow from chiral symmetry beyond those of SU(N f ) V in the massless limit except that α 0 = β 0 = 0. This is not very surprising, derivative couplings in the massless limit are always soft for all legs for a four meson vertex, see [18].

Conclusions
In this manuscript we have determined the order p 8 or NNNLO chiral perturbation theory Lagrangian for the purely mesonic case of positive intrinsic parity. We have reproduced the known results at lower orders and discussed the cases of general number of light flavours and with two and three light flavours. We also discussed the cases where the types of external fields are restricted to either scalar-pseudoscalar or vector-axial-vector as well as the case with no external fields. We separately determined the number of contact terms.
Our main results are the number of terms given in tables 3-7. The list of monomials is included as supplementary material but the contact terms have been listed in the appendix in table 8.
As an example of how our results can be used we looked at meson-meson scattering and the general form of the amplitude and compared it with the freedom our Lagrangian allows.   Table 8: The list of contact terms and the number in the minimal basis for a general number of flavours N f and for N f = 2, 3.