Finite-Volume QED Corrections to Decay Amplitudes in Lattice QCD

We demonstrate that the leading and next-to-leading finite-volume effects in the evaluation of leptonic decay widths of pseudoscalar mesons at $O(\alpha)$ are universal, i.e. they are independent of the structure of the meson. This is analogous to a similar result for the spectrum but with some fundamental differences, most notably the presence of infrared divergences in decay amplitudes. The leading non-universal, structure-dependent terms are of $O(1/L^2)$ (compared to the $O(1/L^3)$ leading non-universal corrections in the spectrum). We calculate the universal finite-volume effects, which requires an extension of previously developed techniques to include a dependence on an external three-momentum (in our case, the momentum of the final state lepton). The result can be included in the strategy proposed in Ref.\,\cite{Carrasco:2015xwa} for using lattice simulations to compute the decay widths at $O(\alpha)$, with the remaining finite-volume effects starting at order $O(1/L^2)$. The methods developed in this paper can be generalised to other decay processes, most notably to semileptonic decays, and hence open the possibility of a new era in precision flavour physics.


I. INTRODUCTION
For many physical quantities relevant for studies of flavour physics, recent improvements in lattice computations have led to such a precision that electromagnetic effects and isospin breaking contributions cannot be neglected anymore (see e.g. Ref. [2] and references therein).
For light-quark flavours, important examples include the calculations of the leptonic decay constants f K and f π and of the form factor f + (0) in semileptonic K 3 decays. These are used to determine the CKM matrix element |V us | and the ratio |V us |/|V ud | at high precision.
For such quantities, which have been computed with a precision at the sub-percent level, the uncertainty due to the explicit breaking of isospin symmetry (of the order of (m u − m d )/Λ QCD ∼ 0.01) and to electromagnetic corrections (of the order of α ∼ 0.007) is similar to, or even larger than, the quoted QCD errors [2].
The question of how to include electromagnetic effects in the hadron spectrum and in the determination of quark masses in ab-initio lattice calculations was addressed for the first time in [3]. Indeed, using a variety of different methods, several collaborations have recently obtained remarkably accurate results for the hadron spectrum, for example in the determination of the charged-neutral mass splittings of light pseudoscalar mesons and baryons [4][5][6][7][8][9][10][11][12][13][14] (see [15][16][17] for reviews on the subject).
In a recent paper, a new proposal to include electromagnetic and isospin-breaking effects in the non-perturbative calculation of hadronic decays was presented [1]. As an example of the new method, the procedure to compute O(α) corrections to leptonic decays of pseudoscalar mesons was described in detail. This can then be used to determine the corresponding CKM matrix elements.
There is an important point that needs to be stressed here. Whereas in the computation of the hadron spectrum there are no infrared divergences, in the calculation of the electromagnetic corrections to the hadronic amplitudes infrared divergences are present and only cancel for well defined, measurable physical quantities. This requires diagrams containing different numbers of real and virtual photons to be combined [18]. The presence of infrared divergences in intermediate steps of the calculation requires the development of a strategy which is different, and more complicated, than the usual approaches followed to compute the electromagnetic corrections to the spectrum. We proposed such a strategy in Ref. [1].
There we envisaged that at O(α) the physical observable is the inclusive decay rate of the pseudoscalar meson into a final state consisting of either −ν or −ν γ, with the energy of the emitted photon in the rest frame of the pion smaller than an imposed cut-off ∆E. Here − is a charged lepton and ν the corresponding neutrino. The cut-off ∆E on the energy of the final-state (real) photon should be sufficiently small that for photons with such an energy we can neglect the structure of the meson and treat it as an elementary point-like can be evaluated in perturbation theory directly in infinite volume and the result has been presented in Ref. [1]. Γ pt (∆E) is infrared finite and independent of the scheme used to regulate the divergences which are present separately in Γ pt 0 (L) and Γ pt 1 (∆E, L); its explicit expression is reproduced in Eq. (13) below.
Γ 0 (L) is infrared divergent and depends on the infrared regularisation. Since all momentum modes of the virtual photon contribute to Γ 0 , it depends on the structure of the meson and is necessarily non-perturbative. It should therefore be computed in a lattice simulation.
In Ref. [1] we stressed that the infrared divergence cancels in the difference Γ 0 (L) − Γ pt 0 (L). In this paper we show that the 1/L finite-volume (FV) corrections are also universal, that is they are independent of the structure of the pseudoscalar meson and hence cancel in the difference Γ 0 (L) − Γ pt 0 (L). We do this in Appendix A using the QED skeleton expansion, in which the meson propagator and the vertices to which the photon couples, are defined in terms of QCD correlation functions and then inserted into one-loop diagrams. Combining the skeleton expansion with the electromagnetic Ward identities of the full theory, we are able to demonstrate explicitly that the leading and next-to-leading FV effects are universal. This allows us to calculate Γ pt 0 (L) in perturbation theory with a point-like pseudoscalar meson up to and including the 1/L corrections and present the result expanded in inverse powers of L Γ pt 0 (L) = C 0 (r ) +C 0 (r ) log (m P L) + where r = m /m P and m P and m are the masses of the pseudoscalar meson and the lepton respectively. Throughout this paper we will refer to the first two terms on the right-hand side of Eq. (4) as the leading FV effects, and the third term, C 1 /m P L, as the next-toleading correction. The explicit expression for Γ pt 0 (L) is given in Eqs. (97) and (98) below. The coefficients C 0 (r ),C 0 (r ) and C 1 (r ) are universal, although C 0 (r ) and C 1 (r ) depend on the infrared regulator.C 0 (r ) is universal and does not depend on the regularisation.
C 0 (r ),C 0 (r ) and C 1 (r ) cancel the corresponding terms contained in Γ 0 (L). In this way Γ 0 (L) − Γ pt 0 (L) is infrared finite and independent of the infrared regularisation up to terms of O(1/L 2 ). Higher order FV terms are not universal and thus cannot be corrected with an analytic computation. We do not discuss them further, beyond showing in Sec.V that they are indeed non-universal.
The discussion in the previous paragraph has parallels in the calculation of the FV corrections to the spectrum [4,5,7,9,12,13,19]. In that case there are no infrared divergences and the O(1/L) and O(1/L 2 ) FV corrections are universal, but the O(1/L 3 ) corrections are structure-dependent. For matrix elements the leading dependence on the volume in Γ pt 0 is an infrared divergence of the form log(m P L) and the next-to-leading term is of O(1/L).
Both of these are universal.
Below we present the perturbative one-loop calculation of Γ pt 0 (L) on a finite volume using QED L [9] as the infrared regulator. We will describe in detail the method developed for the calculation of the perturbative corrections to decay amplitudes in a finite volume; these calculations are more difficult than the corresponding evaluation of the corrections to hadron masses. In addition to the presence of infrared divergences, even in the rest-frame of the meson there is a dependence on the three-momentum of the final-state lepton from the diagram in which the photon is emitted from the meson and absorbed by the lepton. When evaluating the FV corrections, the summand in the summation over the spacial momentum modes of the photon, k, depends not only on k but also on p · k, i.e. on the direction of the final state lepton's momentum, p , with respect to the axes of the cubic lattice.
This complicates the calculation significantly and leads to results which also depend on the direction of p . We believe that the techniques developed in this paper, which extend those of Ref. [20], have a wider applicability and will be useful for many other processes.
Although our explicit discussion is limited to the leptonic decay rates of pseudoscalar mesons, the method is general and can be extended to many other processes including, for example, to semileptonic decays. We should add however, that although the results presented in this paper are valid in principle for both heavy and light pseudoscalar mesons, there may be a practical limitation in the case of the heavy D and B mesons. In that case it is likely that in order to make experimental measurements feasible, ∆E may have to be sufficiently large that the structure dependence of the meson can no longer be neglected and therefore that the emission of real "hard" photons, with energies of E γ ≥ Λ QCD should be implemented in the lattice simulation. We do not discuss the prospects for this further in this paper.
The plan for the remainder of this paper is the following. In the following section (Sec. II) we present the decay rate without electromagnetic corrections and introduce some basic notation used in the subsequent sections. Section III contains a discussion of the regularisation of the ultraviolet divergences and the W -renormalisation scheme. Since the ultraviolet divergent terms are unaffected by FV effects, we simply sketch the renormalisation procedure referring to our previous paper for further details [1]. In Sec.IV we review the method proposed in Ref. [1] for the cancelation of infrared divergences and present an extended discussion of the different proposals for their regularisation; The universality (or non-universality) of the infrared divergences and FV corrections is explained in Sec. V. In particular we sketch the demonstration that the leading and next-to-leading FV effects are universal. The perturbative calculation of the electromagnetic corrections to the leptonic decay amplitude and meson mass on a FV, including the O(1/L) corrections, is presented in full detail in Sec. VI.
All the results are expressed in terms of a few master integrals. The evaluation of the oneloop master integrals is performed in Sec. VII. The calculations described in this section are of general use and can be applied to many other cases of phenomenological interest.
Finally, in Sec. VIII we present our final result, our conclusions and the outlook for the implementation of our method. There is a single appendix in which the universality of the leading and next-to-leading FV effects is proved using the skeleton expansion.
the circles are displaced for convenience and the index L represents left.

II. THE DECAY RATE WITHOUT ELECTROMAGNETIC CORRECTIONS
At lowest order in electromagnetic perturbation theory (i.e. at O(α 0 )), the processq 1 q 2 → −ν can be written in terms of the amplitude of an effective four-fermion local Hamiltonian where G F is the Fermi constant, the subscripts i = 1, 2 on q i denote the flavour of the quarks and V 12 is the corresponding element of the Cabibbo-Kobayashi-Maskawa (CKM) matrix.
We illustrate the Feynman diagram for the leptonic decay of a π − meson in pure QCD in Fig. 1. In the absence of electromagnetism the non-perturbative amplitude for the decay of a pseudoscalar meson P − is defined in terms a single number, the corresponding decay constant f P : where P − is composed of the valence quarksq 1 and q 2 , and the axial current in (6) is composed of the corresponding quark fields. From Eqs. (5) and (6) one readily derives the tree level decay rate Since we aim to determine the width up to and including O(α) contributions, m P in Eq. (7) is the physical mass of the meson.
The calculation of electromagnetic corrections leads to an immediate difficulty: Γ 0 contains infrared divergences and by itself is therefore unphysical. The well-known solution to this problem is to include the contributions from real photons. The physical, infrared safe, experimentally measurable observable is then the partial width given in Eq. (1) in the Introduction. Γ(∆E) is free from infrared divergences, and will be computed following the procedure proposed in Ref. [1] and briefly reviewed in section IV below. We start however, with a brief discussion of the renormalisation of the ultraviolet divergences which arise from virtual photon exchanges.

SCHEME
The standard method used in weak leptonic and semileptonic decays to renormalise the theory is to work in the so called W -renormalisation scheme [21]. We refer the reader to Ref. [1] for more details of the applications to the present case.
When including the O(α) corrections, the ultraviolet divergences are removed by defining the Fermi constant G F in the W -renormalisation scheme. Its value is then given in terms of the physical muon decay rate Γ µ = 1/τ µ , where τ µ is the lifetime of the muon: In the same scheme the effective Hamiltonian of eq. (5) gets a finite correction where the four fermion operator (q 1 γ µ (1 − γ 5 )q 2 ) (¯ γ µ (1 − γ 5 )ν ) is also renormalised in the same scheme.
The renormalisation of the weak Hamiltonian requires the evaluation of the diagrams in Fig. 2 (shown for the case P − = π − ) with the photon propagator defined in the Wregularisation scheme Since we are not able to implement the W -regularisation directly in present day lattice simulations in which the inverse lattice spacing is much smaller than M W , the relation between the operator in eq. (9) in the lattice and W regularisations can be computed in perturbation theory. For lattice actions in which chiral symmetry is broken by the Wilson term (or for related actions), this relation takes the form: where α/4π × γ em = α/4π × 2 is the one-loop, regularisation independent electromagnetic anomalous dimension of the four fermion operator and The numerical values of the coefficients C 1 . . . C 5 corresponding to the Wilson action for both the quarks and gluons can be found in Ref. [1]. This concludes the discussion of the treatment of the ultraviolet divergences and of the definition of a finite four-fermion operator in the W -renormalisation scheme.

IV. INFRARED DIVERGENCES AND FINITE VOLUME CORRECTIONS
Having performed the renormalisation, all the expressions for widths and amplitudes are explicitly free from ultraviolet divergences. Moreover both terms appearing on the right-hand side of Eq. (2) are separately infrared finite and also independent of the infrared regularisation. We will demonstrate below that in the first term, Γ 0 (L) − Γ pt 0 (L) is independent of the regularisation up to terms of O(1/L 2 ). The independence of the last term from the infrared regularisation can be readily demonstrated in perturbation theory; Γ pt (∆E) ≡ Γ pt 0 + Γ 1 (∆E) is a well-defined physical quantity, corresponding to the decay rate of a point-like particle calculated at O(α) in infinite volume. It is given by the following expression [1]: where r E = 2∆E/m P and 0 ≤ r E ≤ 1−r 2 . Note that the terms in square brackets in eq. (13) vanish when r E goes to zero; in this limit Γ pt (∆E) is given by its eikonal approximation.
Since Γ pt (∆E) is itself independent of the infrared regulator, it can be computed with a different regularisation from the one used in computing the difference Γ 0 − Γ pt 0 . This implies that, provided that we use the same infrared regulator for Γ 0 and Γ pt 0 , the infrared divergences, which are universal, cancel in the difference Γ 0 − Γ pt 0 , leaving an O(α) finite term which is independent of the regulator. We stress that the regulator-dependent finite terms also cancel in the difference Γ 0 − Γ pt 0 . We now discuss FV effects. Let V = L 3 be the spacial volume and for simplicity we take the length L in each direction to be the same. Whereas in QCD the FV corrections are exponentially suppressed, for electromagnetic corrections, because of the presence of a massless photon, FV effects are particularly important since they are only suppressed by inverse powers of the volume. The 1/k 2 term in the photon propagator in Eq. (10) implies the presence of a zero-mode in the finite-volume summation over momenta. Several suggestion have been proposed in the literature for the treatment of the zero mode of FV QED: 1. In the first proposal the four momentum zero mode of the photon field is removed, i.e.
A µ (k = 0) = 0. This is denoted as QED TL [3]; 2. The second proposal, denoted by QED L , is to remove the three-momentum zero modes of the photon field, i.e. to set A µ (k 0 , k = 0) = 0 for all k 0 [9]; 3. A traditional way to regulate infrared divergences in QED is to give the photon a small mass. This is denoted by QED γ [10]; 4. Finally, the fourth proposal is to enforce C * boundary conditions for all fields along the spatial directions, i.e. to require that the fields are periodic up to charge conjugation.
In this theory, which we refer to as QED C , the zero-modes of the gauge field are absent by construction because A µ (x) is anti-periodic in space [19].
Although, at first sight, it may appear that regularising the theory by giving a mass to the photon is the safer option, with presently available lattice volumes this approach has several major drawbacks [17] and we prefer to use the finite volume itself as the infrared regulator.
For the hadron spectrum, both the O(1/L) and the O(1/L 2 ) corrections are universal, that is they are the same for point-like and composite hadrons [4,12,19] and can be analytically computed (they do however depend on the regulator). We have chosen to work in QED L .
In this regularisation one finds for the electromagnetic mass shift for a pseudoscalar meson of charge q in a large, but finite, volume [4]: where κ = 1.41865 is a universal constant. Eq. (14)  that the FV corrections are universal, the authors of Ref. [12] used an approach based on non-relativistic effective field theories, whereas in Ref. [4] an independent demonstration of the universality of the corrections, based on the electromagnetic Ward identities of the theory, was used. We will discuss this in more detail in the following.
Using similar approaches we can demonstrate, see Sec. V and Appendix A below, that for the amplitude with the virtual photon, the coefficients C 0 (r ) ,C 0 (r ) and C 1 (r ) in Eq. (4) are universal so that the FV corrections to the difference . This should be compared to the O(1/L 3 ) at which the structure dependent FV corrections begin to contribute in the spectrum; the difference, as explained in detail below, is due to the different behaviour of the integrands as the momentum of the photon goes to 0.
We note also that, since the sum of all the terms in Eq. (2) is gauge invariant, as is the perturbative rate Γ pt (∆E), the combination ∆Γ 0 (L) is also gauge invariant, although each of the two terms on the right-hand side of Eq. (15) in general is not.
In the remainder of this paper, all the calculations are performed using the QED L prescription for handling the zero mode. If the photon were massive so that then the FV corrections would be exponentially suppressed in the mass of the photon. The power-law FV corrections arise with a massless photon from the singularities of the summand at k = 0 in the sum over the momentum k of the photon. The power of the corrections depends on the degree of the singularity in k of the summand. Thus in order to study the FV corrections we consider the soft region in which k ≈ 0.
We start by introducing our notation and some basic definitions. Let ξ be the difference between the finite and infinite volume result of some one-loop expression where k is the momentum of the photon, p represents the external momentum or momenta and the prime on ξ indicates that the contribution from the spacial zero mode ( k = 0) is absent from the sum. Although Eq. (17) includes both summations and integrations, in the following we will refer to I(p, k 0 , k ) in expressions for ξ as the integrand. A practical rule summarising the relation between the power of the finite-volume corrections and the leading singularity of the integrand at k = 0 is as follows: Thus for example, if the integrand I in Eq. (17) has a term which behaves as O(1/k 3 ) as k → 0, then the corresponding FV corrections is of O(1/L). In Sec. VII we will demonstrate the scaling rules in Eqs. (18) and (19) explicitly for the integrals which appear in the oneloop graphs for the leptonic decays of pseudoscalar mesons. This includes the extension of previous applications to contributions which depend on external three-momenta and in particular on their directions (in our case this is the momentum of the final-state lepton in the rest frame of the meson). Up to the next-to-leading-order to which we work, these scaling rules hold with k being the momentum of the photon. At higher orders, other regions of phase-space may contribute to power corrections in the volume; these are still given by Eq. (18), but k is no longer the momentum of the photon. An example is the contribution to O(1/L 3 ) corrections to the spectrum found by the BMW collaboration [13], which arises from the region of small k where the photon's momentum is written as (m P + k 0 , k ).
It is convenient to discuss the different Feynman diagrams within the framework of the skeleton expansion and in terms of hadronic vertices and propagators and these are sketched in Fig. 3. The full propagators and vertices are defined explicitly in terms of correlation functions in Appendix A; for the purposes of this section we can view them as one-particle irreducible subgraphs. The amplitude itself is obtained from the lattice computation of correlation functions in QCD+QED as described above. We now discuss the strategy for the evaluation of the leading and next-to-leading FV corrections; the detailed evaluation is performed in Secs. VI and VII. The key observation, already mentioned in the previous paragraph, is that (up to this order) these effects are determined by the behaviour of the integrands as k → 0, where k is the momentum of the photon. In order to calculate these FV effects, we need the vertices in Fig. 3 for small photon-momenta k.
In their studies of the FV corrections to the hadronic spectrum, the authors of Refs. [7,12,13] used the following non-relativistic effective theory of a charged (pseudo)scalar particle interacting with soft photons: In Eq. (20) φ represents the field of the pseudoscalar meson P using the non-relativistic normalisation, D µ = ∂ µ − ieA µ is the covariant derivative, r 2 is the mean squared charged radius of P andα E andβ B are related to the electric and magnetic polarizabilities of P .
The terms proportional to r 2 ,α E andβ B depend on the structure of P and so additional information about these parameters is required if the corresponding effects are to be included.
These terms all include two derivatives on the photon field, and hence two powers of k in momentum space, and therefore only enter at Next-to-Next-to-Leading Order at small k and are therefore suppressed by two powers of 1/L. In this paper we do not attempt to evaluate these corrections analytically, but can envisage fitting the behaviour numerically if necessary or appropriate. The key point to notice is that the first two terms on the right-hand side of Eq. (20), which include zero or one derivative on the photon field, are universal, i.e. they do not depend on the structure of P and can be evaluated in perturbation theory for a point-like charged particle. A more formal discussion of this universality, based on the QED Ward Identities rather than relying on the effective theory is presented in Appendix A.
At this point we should stress the limitations of the effective theory at O(1/L 3 ). Consider a generic contribution which is non-singular in the infrared. From the Poisson summation formula one would expect that the corresponding FV effects to be suppressed exponentially in the volume. In QED L however, the contribution from the zero mode is removed and this leads to a 1/L 3 effect which is not captured by simple power counting. Since in this study we limit our interest in FV corrections to O(1/L 2 ) for the spectrum and O(1/L) for matrix elements we are not effected by this limitation.
We now illustrate the implications of the discussion in this section, starting with the self-energy diagrams in Figs. 3(b) and 3(c). In the discussion of the universality of the FV corrections which follows, and in the explicit calculations with a point-like pseudoscalar  14).
In Appendix A we show explicitly that the Ward identities of scalar QED are sufficient to demonstrate the universality of the leading and next-to-leading FV corrections. Once it is established that the two leading terms as k → 0 are universal and independent of the structure of the meson, they can be computed in one-loop perturbation theory with a point-like meson using relativistic scalar QED or with the non-relativistic theory. For the remainder of this paper we perform the calculations in scalar QED, restricting the results, of course, to the universal terms.
When evaluating the contribution to the decay amplitude from the wave-function renormalisation of the pseudoscalar meson P there is an additional consideration. At O(α) this is obtained from the coefficient of p 2 of the one-loop diagram in Fig. 3(b) evaluated on-shell, i.e. at p 2 = m 2 P . This requires differentiating the integrand where p is the momentum of the external pseudoscalar meson, and subsequently going on shell. In this way we recover the well-known infrared divergent behaviour, i.e. the leading behaviour as k → 0 corresponds to n = 4 in Eq. We now have all the ingredients necessary to calculate Γ pt 0 (L) and hence to determine the coefficients C 0 (r ),C 0 (r ), and C 1 (r ) of Eq. (4). We have shown that we need to evaluate the diagrams of Fig. 3 using scalar QED for the vertices Γ and the propagators and the effective weak Hamiltonian given in Eq. (21) below for the weak vertices W . In the next section we calculate these diagrams, presenting the results in terms of master integrals which are subsequently evaluated in Sec. VII.
In this section we describe the calculation of the second term on the right-hand side of Eq. (15), Γ pt 0 (L), at O(α). We start however, by briefly recalling the calculation of Γ pt 0 at O(α 0 ), i.e. without electromagnetism. In the following p, p and p ν are the momenta of the meson P , the charged lepton and the neutrinoν respectively.
Here Φ P represents the field of the meson and V ij is the corresponding element of the CKM matrix. For compactness of notation we have dropped the labels i, j from both Φ ij P and the decay constant f P ij .
Without electromagnetic corrections we need to compute the Feynman diagram of Fig. 3 (a), which is a standard calculation. Since the leptonic terms are factorised from the hadronic ones, the amplitude is simply given by obtained by using the equations of motion and neglecting the mass of the neutrino, so that The decay width is then readily obtained by squaring the amplitude A 0 and integrating over the phase space leading to the result in Eq. (7).
For later convenience we write We now consider the one-photon exchange contributions to the decay P − → −ν starting from the crossed diagram in Fig. 3(d). In infinite volume this takes the form Using the equations of motion, the numerator In a finite volume the momentum integration is replaced by a summation over the momenta which are allowed by the boundary conditions. We introduce an infrared cutoff λ which will be useful in intermediate steps of the calculation. We envisage that λ 1/L, but otherwise the specific choice of λ is immaterial since the final result in a finite volume is independent of λ in the limit λ → 0. We then write The contribution from N 11 , which we write as X FV 11 , can readily be written as a multiplicative correction to the lowest-order amplitude, where r = m /m P and we have introduced the master integrals In Eq. (29) and in the following the master integrals are all dimensionless.
The contribution of N 12 is more difficult to evaluate, but we do not need its explicit form when evaluating the rate. Convoluting this term with the (complex conjugate) of the lowest-order contribution we obtain a contribution to the decay rate which is proportional where By combining together the results in Eqs. (28) and (30) we obtain the effective correction to the lowest order amplitude The calculation and the explicit expressions of the integrals S i can be found in sec. VII.
We now discuss the remaining Feynman diagrams.

C. The rainbow diagrams
In this subsection we evaluate the diagrams in Fig. 3(e) and 3(f). The expression for the diagram in Fig. 3(f), in which the photon is emitted by the pion and absorbed by the weak Hamiltonian, is given by where In this way one obtains We now consider the diagram of Fig. 3(e) in which the photon is emitted from the weak Hamiltonian and absorbed by the charged lepton : where One then readily obtains where S 5 is a new integral and sum defined by

D. Wave function renormalisation of the pseudoscalar meson and charged lepton
The wave-function renormalisation constant of the pseudoscalar meson, Z P , is determined from the diagram in Fig. 3(b) which is given by the expression: The renormalisation constant Z P is given by: After performing the integration over k 0 , the result is particularly simple: where for convenience we have introduced the infrared cutoff λ as in the previous cases.
Here, as in the evaluation of the master integrals in general (see Eqs. (51) and (52)) we organise the calculation as follows: The infinite volume result δZ IV P is ultraviolet divergent and must be regularised in the W -scheme The result of the integration over k is ξ P is defined as Using the techniques of Ref. [20] we find where the numerical value of K P is K P = 4.90754. The explicit integral expression for ξ P is given in Sec. VII E. Note that the spacial zero mode is included in the sum in ξ P whereas it is not included in the definition of δZ P in Eq. (42). We therefore subtract this contribution explicitly; this is the term −∆Z P = −e 2 /(2L 3 λ 3 ) on the right-hand side of Eq. (43).
Collecting all the terms together we obtain In spite of the title of this subsection we do not need to evaluate the wave function renormalisation constant of the charged lepton . As explained in Ref. [1], its contribution simply cancels in the difference between the perturbative and non-perturbative contributions to the decay rate in Eq. (15).

VII. MASTER INTEGRALS IN A FINITE VOLUME
In this section we discuss the calculation of the master integrals introduced in sec. VI, which is valid if the integrand f i k 0 , k does not have singularities in the limit L → ∞, we decompose the master integrals as follows The three terms in Eq. (52) correspond to the three terms on the right-hand side of Eq. (51).
S IV i is the infinite-volume result which may have logarithmic ultraviolet or infrared divergences. The ultraviolet divergences are eliminated by the W -regularisation, as explained above, whereas the infrared divergences are regulated by the introduction of λ. The difference ξ for a general integrand I is defined as a simple modification of ξ defined in Eq. (17) to include the contribution from the spacial zero mode: In this section we evaluate the ξ i (i = 1-5) corresponding to the S i defined in Sec.VI. The term ∆S i accounts for the fact that the spacial zero mode k = 0 is included in the sum in the Poisson summation formula (50) but not in Eq.(51). The differences ξ i contain power divergences of the form 1/λ n , where n is a positive integer, which cancel those in ∆S i , whereas the logarithmic divergences of the form log λ cancel between ξ i and S IV i . Thus the final result does not depend on λ in the limit λ → 0. The logarithmic infrared divergence appears instead as a logarithm of the volume, namely as log(m P L). There remain of course finite terms which depend on the infrared regularisation, i.e. on the definition of the photon propagator in finite volume (e.g. whether we use QED TL , QED L or some other definition).
As was shown in Sec. V however, these terms are universal, in that they are the same in the point-like perturbative calculation and in the non-perturbative computation of the hadronic amplitude in a numerical simulation of QCD. Thus the dependence on the regularisation cancels in the finite difference of Eq. (15) and only remains in the power-suppressed, nonuniversal FV terms which go to zero as 1/L 2 or faster, as L → ∞.
We now discuss the evaluation of the master integrals in QED L , starting with the most complicated one, S 3 , which has three different denominators.

A. Calculation of S 3
We begin the discussion of the master integrals with S 3 , defined in Eq. (29) of sec. VI because it contains all the main features and difficulties of the FV calculations. Since the infinite-volume integrals are straightforward to evaluate, we use the Poisson summation formula to evaluate the difference between S 3 and the corresponding integral, i.e. we evaluate where α, β are Feynman parameters. Note that since we use the Poisson formula, it is ξ 3 and not ξ 3 which we evaluate here; later we will subtract the spatial zero mode separately.
Changing integration variables, k → k + αp + βp , and performing the k 0 integration we obtain where we have taken the meson P to be at rest, p P = 0. The summation in Eq. (56) is over vectors of integers n, with n = 0. We now observe explicitly the technical complication mentioned above, namely the dependence of ξ 3 on p . We proceed by generalising the techniques developed in ref [20] to such a case, noting that for positive X and writing The integrand has been manipulated so that the integration over k is Gaussian and can readily be performed to obtain The next step is to recognise that the sum over n can be written in terms of the Jacobi elliptic theta function θ 3 : although we find it more convenient to present the calculation and results in terms of In terms of this θ function, ξ 3 is given by where p i is the i-th component of p (i = 1, 2, 3) and in the last line we have changed the integration variable to τ = 4πt/L 2 . We aim to extract the L and λ dependences and to do this we consider three separate contributions to the right-hand side of Eq. (63).
(i) We start by considering the region τ < 1, which gives an infrared convergent contri- up to exponentially small corrections in the volume. In Eq. (64) we have introduced the variable ρ = α + β and then rescaled it by L to absorb the factors of L 2 in Eq. (63) and finally set L → ∞. K 31 is a finite number, which can readily be evaluated numerically.
Note that it depends not only on the masses m P and m but also on the orientation of the momentum p with respect to the axes of the lattice.
(ii) The region τ > 1 is split further into two contributions, starting with where we have taken the −1 term from the square brackets in Eq. (63). It is possible to evaluate ξ 32 analytically and we obtain (up to terms which vanish exponentially with the volume) where γ 0.577216 is Euler's constant and r = m /m P .
(iii) Finally we have to evaluate ξ 33 , where Exploiting the Poisson summation formula for the θ-function and changing variables to t = 1/τ , ξ 33 can be rewritten as In writing Eq. (69) we have subtracted and added 1 on the second line. We will see that it is the first term which contains the expected 1/λ 3 contribution from the spacial zero mode, whereas the second term is finite and can be evaluated numerically. The first term is where we have changed integration variables from t to 1/t and from α, β to ρ = α + β and β = β/ρ. The ρ integration can now be performed to give where we have dropped the tilde on the integration variable β. The second term in the numerator of Eq. (71) gives a contribution which is exponentially suppressed in the volume, so that The infrared cut-off λ is needed in Eq. (72) to regulate the integration over t which can be performed to give where we have noted that (p − p ) 2 = p 2 ν = 0. The β integration can also be performed, where E is the energy of the charged lepton , so that finally We now return to Eq. (69) and evaluate the second term in the braces where The integral over t is infrared finite and so we have set λ = 0 in Eq. (77) and have also taken L → ∞ in Eq. (78). K 32 can be evaluated numerically and again depends on the direction of the lepton's momentum. This completes our calculation of ξ 33 = ξ 331 + ξ 332 .
This also completes our calculation of ξ 3 . Collecting the results from Eqs. (64), (66), (75) and (78) the result is where the expressions for the L and λ independent constants K 31 and K 32 can be found in Eqs. (64) and (79) respectively.
In order to obtain S 3 , in addition to ξ 3 we need to determine S IV 3 and ∆S 3 (see Eq. (51)). S IV 3 is the corresponding IV integral which can be evaluated using standard perturbative techniques: Since S 3 has no ultraviolet divergences, the terms obtained by subtracting the 1/(k 2 − M 2 W ) term in the W -regularised photon propagator are suppressed by a factor of 1/M 2 W and can be neglected.
The contribution from the spatial zero mode only requires an integration over k 0 and is The sum S 3 is then given by where we have replaced E by m P (1 + r 2 )/2. The dependence on λ has disappeared as anticipated. The 1/λ 3 term in ξ 3 (see Eq. (80)) is simply the term from the spatial zero mode and is indeed cancelled by ∆S 3 , whereas the log λ 2 terms cancel between ξ 3 and S IV 3 .

B. Calculation of S 1
In the calculation of S 1 defined in Eq. (29) we follow the same steps as for S 3 in Sec VII A with the simplification that in this case we only have two propagators instead of three. On the other hand, this integral is ultraviolet divergent and must be regularised: We find where K 11 0.0765331 and K 12 0.0861695 are mass-independent dimensionless constants which are defined in Eq. (92). Collecting these terms together we obtain

C. Calculation of S 2
The calculation of S 2 , defined in Eq. (29), is similar to that of S 1 . The integral is ultraviolet divergent and must be regularised We find where the constants K 21 and K 22 , which are dimensionless, are given in Eq. (92). Note that K 21 and K 22 depend on the direction of the lepton's momentum p with respect to the axes of the lattice. Collecting together the terms in Eq. (88) we obtain We note from Eqs. (32), (35) and (38) that S 4 and S 5 enter in the expressions for the diagrams with a factor of 1/(1 − r 2 ), cancelling the corresponding factors in Eqs. (90) and (91).

E. The auxiliary constants K ij
We present here the explicit expressions for the real constants K ij appearing in the expressions of the master integrals S 1 -S 5 . The θ-function is defined in Eqs. (60) and (61).
The separate appearance of K i1 and K i2 (i = 1, 2, 3) is a consequence of how we chose to organise the calculation; for example in the evaluation of ξ 3 we split the integration over τ in Eq. (63) into contributions from τ < 1 and τ > 1 with K 31 coming from the first region and K 32 from the second.
. For completeness we also give the expression for ξ P from Eq. (48): where we are using the notation of Ref. [20], 1 Writing the width in the point-like theory up to O(α) and including FV corrections as inserting the expressions for the master integrals from Sec VII into the result for X FV in Eq. (96) and adding the contribution from the wave-function renormalisation of P in Eq. (49) we find Our work has close parallels to the studies of electromagnetic corrections to the spectrum [4,5,7,9,12,13,19] where the leading (in that case O(1/L)) and next-to-leading (O(1/L 2 )) FV corrections are universal. In our case it is the coefficients C 0 (r ),C 0 (r ), and C 1 (r ) of Eq. (4) which are universal and which are obtained from Eq. (98).
In addition to the electromagnetic corrections studied in this paper, one also needs to account for comparable isospin-breaking effects due to the difference in the up-and down-quark masses. This is a technical complication, rather than a conceptual issue and we have not discussed it in this paper.
Although, the explicit expression presented in Eq. (98) corresponds to the leptonic decay of pseudoscalar mesons, the methods developed in Ref. [1] for the handling of infrared divergences and extended in this paper to evaluate the leading and next-to-leading FV corrections can be generalised to other decay processes, most notably to semileptonic decays.
We envisage that this will lead to a significant improvement in the precision of flavour physics.
coupled to the vertices.

a. The meson propagator
We define the two-point correlation function where φ P is an interpolating operator for the meson P and T represents time-ordering. The propagator is then defined by where |P ( p ) is the state of the meson P with three-momentum p. It is assumed that P is the lightest state which can be created by φ † P . The denominator on the right-hand side of Eq. (A2) is obtained in the standard way from the two-point correlation function in Eq. (A1) without integrating over the time and at sufficiently large times so that only the ground-state P contributes.
Inserting a complete set of states |n( p ) between the two operators in Eq. (A1) and performing the integration over x we obtain the following expression for the propagator where for the excited state we have used E 2 n ( p ) = E 2 n ( 0 ) + p 2 . In Eq. (A3) m P is the mass of the meson in QCD; its physical mass will be modified at O(α) as explained in Sec. A 3 a below. The second term in the braces in Eq. (A3) contains the effects of the excited states and vanishes on-shell, i.e. as p 2 → −m 2 P .

b. The meson-photon vertex
The coupling of the charged meson to a single photon, denoted by Γ in diagrams (b), (d) and (f) of Fig. 3 is defined in terms of the three-point correlation function In Eq. (A4) j µ is the electromagnetic current.
Of particular importance in the following will be the electromagnetic Ward Identities.
Under the infinitesimal gauge transformation on the quark fields of flavour f , the operators in C µ in Eq. (A4) transform as follows: and the QCD action transforms as S → S − i d 4 x λ(x) (∂ µ j µ (x)). From the generic nonanomalous Ward identity for a multilocal operator O: we obtain The result in Eq. (A8) is readily rewritten in terms of the meson-photon vertex and propagators as: c. The meson-two-photon vertex At O(α) we also need to consider the meson-two-photon vertex in diagram (c) in Fig. 3. This is defined from the four-point correlation function as follows: The Ward identities for this vertex can be derived as in Sec. A 1 b and give

The weak vertex
For the calculation of the decay amplitude at O(α) we also need to consider the proper vertices of the weak quark current with zero, one or two photons, denoted by W in the diagrams of Fig. 3. We start here with the vertex with no photons, which is obtained from the correlation function where the weak current J ρ W =q 1 γ ρ (1 − γ 5 )q 2 and q 1,2 are the fields of the valence quarks of the meson P . The weak vertex W ρ (p) is then defined by It will be useful to define F (p 2 ) by W ρ (p) = −p ρ F (p 2 ), and the leptonic decay constant f P is defined in the standard way by F (−m 2 ) = f P . We now consider the vertices with one or two photons and derive the corresponding Ward identities.
e. The weak vertex with a single photon The weak vertex with a single photon W µρ (p, k), which is an element in diagrams (e) and (f) of Fig. 3, is defined from the three-point correlation function as follows It satisfies the Ward identity (A17) f. The weak vertex with two photons The final element which we require is the weak vertex with two photons W µνρ (p, k, q) (see diagram (g) in Fig. 3) which is obtained from the four-point function The vertex W µνρ (p, k, q) is defined by where 2Λ µν (p, k, q) = Γ µν (p, k, q)+Γ µ (p, k)∆(p+k)Γ ν (p+k, q)+Γ ν (p, q)∆(p+q)Γ µ (p+q, k) . (A20) The corresponding Ward identities are now This completes the discussion of the elements which enter into the skeleton expansion and the corresponding Ward identities.

The Ward identities at small photon momenta
In this subsection we investigate the consequences of the Ward identities in Eqs. (A9), (A12), (A17) and (A21) on the structure of the vertices at low photon momenta. We start however, with a discussion of the meson propagator ∆(p+k) . In order to determine the wave function renormalisation we need to perform a double expansion of ∆(p + k) in 2 ≡ p 2 + m 2 P and k leading to and Z(p 2 ) is given in Eq. (A3), or equivalently The terms which are not exhibited explicitly in Eq. (A24) are not needed for the eventual evaluation of on-shell matrix elements or for the calculation of the leading and next-toleading FV effects as explained in Sec. V.
The Ward identity in eq. (A9) constrains the vertex Γ µ (p, k) to take the following form at low photon momenta k: Similarly in diagram (c) of Fig. 3 we need Γ µν (p, k, −k) which is constrained by the Ward identity in Eq. (A12) to take the form The weak vertices with one or two photons are similarly constrained by the Ward identities in Eqs.(A17) and Eqs.(A21) to take the form where the f n are the derivatives of F (p 2 ) (defined in Sec. A 1 d) (A28) The terms in the vertices which are not exhibited explicitly do not contribute to the leading or next-to-leading FV corrections.

Universality of leading and next-to-leading FV effects
We now use the vertices above to demonstrate that the leading and next-to-leading FV effects are universal and can be obtained from the calculation of one-loop saclar QED diagrams with point-like charged mesons. Eventually we wish to evaluate the ξ of Eq. (17) corresponding to the diagrams of Fig. 3. Even though the evaluation of the ξ involve both integrals and sums, for conciseness of the terminology for the remainder of this section we will refer simply to integrals and integrands.

a. FV effects in the meson mass
We start the discussion with the O(α) corrections to the meson mass, i.e. the calculation of the diagrams (b) and (c) in Fig. 3 at p 2 = −m 2 P . The leading behaviour as k → 0 of the integrand in diagram (b) is O(1/k 3 ) corresponding to an O( 1 m P L ) FV correction. The integrand is Γ µ (p, k)∆(p + k)Γ µ (p + k, −k)∆ γ (k 2 ) = ∆ γ (k 2 ) 4(−m 2 P + p · k) − 8z 1 m 2 P p · k 2p · k + k 2 where ∆ γ (k 2 ) = 1/k 2 is the photon propagator in the Feynman gauge. The first term in the braces in Eq. (A29) is the one we would obtain in the point-like theory.
We now determine the integrand in diagram (c). Taking  we see that the terms proportional to z 1 cancel and the total is precisely that of the point-like theory so that the electromagnetic shift in the mass is given by the integral over
When the meson is on-shell, i.e. when p 2 = −m 2 P , Σ(p, k) = −δm 2 P . and recalling that we take the external meson to be at rest, we obtain 1 2p 0 ∂Σ(p, k) ∂p 0 The first and second terms on the right hand side of Eq. For a point-like particle, the O(α) correction to square of the wave function renormalisation constant, i.e. ( √ Z P ) 2 , is simply given by Eq. (A33) with z 1 = 0. In the presence of QCD on the other hand, in addition to the second term on the right-hand side of Eq. (A33), there is an effective contribution to Z P from the term proportional to p 2 + m 2 P in the factor of Z present in the two meson propagators which are external to the loop in Fig. 3 (b). Recalling that on-shell Σ(p, k) = −δm 2 P , this contribution to the amplitude is the integral over k of which cancels the z 1 dependent term in Eq. (A33). Thus evaluating the wave-function renormalisation constant in the point-like theory reproduces the leading and next-to-leading FV effects of full QCD.
There is an analogous contribution to that in Eq. (A34) which arises from the expansion of F (p 2 ) at the weak vertex. We postpone discussing this contribution until we study diagram

c. Finite Volume Effects in Diagrams 3 (d) and (e)
We now turn to the diagrams 3 (d) and (e) in which the photon couples to the charged lepton. The coupling to the lepton and the lepton propagator are common in the two diagrams and so we focus on the remainder of the integrand which is not common. The lepton propagator behaves as O(1/k) at small photon momentum k and the photon propagator behaves as O(1/k 2 ). Thus we require the remaining terms up to and including O(k 0 ) in order to obtain the leading and next-to-leading FV effects. In Fig. 3 (d) this is The corresponding factor in diagram 3 (e) is W µρ (p, k) = f P δ µρ + 2f P f 1 p µ p ρ + O(k) .
Summing the results in Eqs. (A35) and (A36) we see that the f 1 -dependent terms cancel and we obtain precisely the expression of the point-like theory.

d. Finite Volume Effects in Diagrams 3 (f )
Consider now the diagram Fig. 3 (f). The product of the meson and photon propagators behave as 1/k 3 as k → 0, corresponding to O(1/L) FV corrections and so we can neglect terms proportional to k in the numerator of this diagram since we neglect corrections of O(1/L 2 ). This diagram therefore gives the integrand ∆ γ (k 2 ) Γ(p, k) ∆(p + k) W µρ (p + k, −k) = ∆ γ (k 2 ) f P p ρ 2 + 4f 1 m 2 P 2p · k + k 2 .
It is natural to combine this with the partial contribution from the diagram in Fig. 3 (b) arising from the expansion of the weak vertex to O(p 2 + m 2 P ). This contribution was mentioned at the end of Sec. A 3 b and the corresponding integrand is The result in Eq. (A38) cancels the f 1 -dependent term in Eq. (A37), leaving precisely the integrand one would obtain in the point-like theory.

Summary
In this appendix we have studied the implications of the electromagnetic Ward Identities on the contributions to the integrands of the diagrams in Fig