A tmQCD mixed-action approach to flavour physics

We discuss a mixed-action approach in which sea quarks are regularised using non-perturbatively ${\rm O}(a)$ improved Wilson fermions, while a fully-twisted tmQCD action is used for valence quarks. In this setup, automatic ${\rm O}(a)$ improvement is preserved for valence observables, apart from small residual ${\rm O}(a)$ effects from the sea. A strategy for matching sea and valence is set up, and carried out for $N_\mathrm{f}=2+1$ CLS ensembles with open boundary conditions at several simulation points. The scaling of basic light-quark observables such as the pseudoscalar meson decay constant is studied, as well as the isospin splitting of pseudoscalar meson masses.


Introduction
Flavour physics is one of the most promising sectors for the search of physics beyond the Standard Model (SM). In particular, processes involving heavy quarks play a crucial role in the ongoing searches of New Physics phenomena. Current (LHCb) and upcoming (Belle II) experiments will considerably improve the amount and precision of available results. At this stage, the verification of the unitarity of the CKM matrix [1,2] as well as the study of rare decays, appear as meaningful ways to examine the compatibility between the SM predictions and the experimental measurements.
Lattice QCD (LQCD) computations in the heavy quark sector are essential to achieve reliable determinations of SM observables. Our main target is to increase the precision of LQCD calculations of hadronic matrix elements necessary for the determination of CKM matrix elements and of other quantities of phenomenological interest in the charm sector. While the accuracy of LQCD computations of decay constants relevant for D (s) -meson decays are in the same ballpark as the current experimental precision of leptonic channels, this is still not the case for other charmed mesons, nor in semileptonic decays. In particular, an accurate computation of the form factors involved in D → π ν and D → K ν decays is essential to making the most of the new experimental results.
We have developed a setup that aims at optimising the control of systematic uncertainties in computations involving charm quarks. In the sea sector, we employ CLS N f = 2 + 1 configurations [3], with open boundary conditions in the time direction [4][5][6], and fine values of the lattice spacing.
Concurrently, the use of a fully-twisted [7][8][9] valence sector for the computation of observablesimplying automatic O(a) improvement up to sea quark mass effects -along with the non-perturbative O(a) improvement of the sea action, is expected to result in excellent scaling properties without the need to tune operator-dependent improvement coefficients.
Here we will concentrate in the development of the mixed-action setup, by identifying an efficient strategy for the matching of the sea and valence sectors and by studying the scaling properties of basic light quark observables. First results will be provided for CLS symmetric point ensembles with degenerate quark masses, m u = m d = m s . We refer to Ref. [10] for a study of charm quark observables using CLS ensembles with the O(a) improved Wilson formulation.

Sea Sector
The set of gauge ensembles used in our study were produced within the CLS initiative [3]. The lattice action involves the Lüscher-Weisz tree-level improved gauge action [11], where p and r run over the plaquettes and rectangles of the lattice and their respective coefficients are c 0 = 5/3 and c 1 = −1/12. The fermionic action comprises a Wilson Dirac operator [12] for N f = 2 + 1 flavours, including the Sheikholeslami-Wohlert term [13], where ∇ µ and ∇ * µ are the forward and backward covariant derivatives and a non-perturbative determination of the coefficient c SW is used [14].
The present analysis is based on four values of the lattice spacing ranging from 0.087 fm to 0.050 fm, as shown in Table 1. The ensembles lie along lines of constant trace of the bare quark mass matrix M q , trM q = 2m q,u + m q,s = const , where m q,f = m 0,f − m cr . In this way, cutoff effects of O(atrM q ), appearing, in particular, in the Symanzik expansion of the bare coupling, can be kept constant when changing the sea quark masses. In the approach towards the physical point, it is in practice beneficial to consider a renormalised chiral trajectory in terms pion and Kaon masses that is related at lowest order in chiral perturbation theory to the condition in eq. (3) in terms of bare quark masses. Specifically, the following dimensionless quantities, depending on the gradient flow observable t 0 defined in section 3.2, are used to establish the chiral trajectory [3,15]. Further details about the choice of the chiral trajectory will be given in section 3.1. In addition to the chiral trajectories targeting the physical point, further ensembles with mass degenerate u, d and s quarks -i.e. along a symmetric line that approaches the three flavour massless theoryare also available for β = 3.46, as indicated in Table 1.  [3] used in the present study. The values of the inverse bare coupling, β = 6/g 2 0 , correspond to the following approximate values of the lattice spacing, a = 0.087 fm, 0.077 fm, 0.065 fm and 0.050 fm [15]. In the third and fourth columns, N s and N t , refer to the spatial and temporal extent of the lattice. Approximate values of the pion and Kaon masses are provided.
The simulations use periodic boundary conditions in the spatial directions and open boundary conditions (OBC) in time. The efficiency of the HMC algorithm to sample the topological sectors rapidly deteriorates when the lattice spacing is reduced below 0.05 fm. The use of OBC allows to connect the field space in the continuum theory by letting the topological charge to flow smoothly in and out of the bulk of the lattice through the boundaries [4][5][6]. As a result, the scaling law of the autocorrelation times is expected to improve significantly. In this way, the results from fine values of the lattice spacing -which are instrumental to control lattice artefacts proportional to the heavy quark mass -should also benefit from a reliable estimate of the statistical uncertainties.

Valence Sector
The mixed action uses in the valence sector a Wilson twisted mass Dirac operator [7-9, 16, 17] where, in addition to a standard mass term and to the chirally rotated mass term, the Sheikholeslami-Wohlert term is also included, This choice of the lattice action enforces that the renormalisation factors of sea and valence sectors coincide.
Maximal twist is achieved once the bare standard mass is set to the critical mass, m 0 = m cr 1, while the twisted mass, µ 0 = diag µ 0,u , µ 0,d , µ 0,s , µ 0,c is associated to the physical bare quark masses. We note that this target physical quark mass can be achieved by using both signs of the twisted mass parameter µ 0,q . At maximal twist, the symmetries of the valence action guarantee the absence of O(a) lattice artefacts proportional to a hadronic scale Λ and to the valence quark masses µ 0,q [8]. While residual lattice artefacts of O(ag 4 0 trM q ) -arising from the sea quark mass matrix M q -can appear, this mixed action setup should profit in the heavy quark sector from the absence of the leading lattice artefacts proportional to the heavy quark masses.

Matching of Sea and Valence Quark Masses
As mentioned in the previous section, at maximal twist, the standard quark mass has to be set to its critical value in order to enforce that its renormalised value vanishes in continuum limit. This can for instance be achieved by imposing that a given observable preserves a continuum symmetry at finite values of the lattice spacing. More specifically, the twisted vector symmetry can be restored by tuning to zero the valence PCAC quark mass, m R 12 v , in the light-quark sector. 1 In practice, this tuning process turned out to be straightforward since, for all the considered ensembles, a short linear interpolation in am 0,u v = 1/ (2κ u | v ) − 4 around the vanishing PCAC quark mass, m R 12 v = 0, was sufficient to fulfil the maximal twist condition. In order to recover the unitarity of the theory in the continuum limit, it is essential to match the valence and sea quark masses on all the ensembles. This can be achieved by enforcing that the renormalised twisted mass is equal to the renormalised PCAC quark mass of the sea, µ R 1 v ≡ m R 12 s . Including O(a) counterterms in the evaluation of m R 12 s [18], this matching condition reads, where the Z P renormalisation factors appearing on both sides of the equation are equal and can therefore be dropped, while the Z A factor has been determined non-perturbatively in Ref. [19]. The bare sea quark mass m 12 includes the O(a) improvement term proportional to c A , which has been computed non-perturbatively in Ref. [20]. The mass-dependent b-type improvement coefficients appearing in eq. (6) have been computed to one-loop in perturbation theory [21] while recent work to determine them non-perturbatively has been reported in Refs. [22,23]. Alternatively, the quark mass matching can be done by tuning the sea and valence pion masses to be equal, m π | v ≡ m π | s . In this approach, both side sides of the equation are automatically O(a) improved.

Observables
The presence of open boundary conditions in the Euclidean time direction modifies the theory in the neighbourhood of the boundaries. The ground state in the boundary |Ω is different from the ground state in the bulk |0 , which asymptotically corresponds to the vacuum of the theory. In the transfer matrix formalism, a typical two-point function is given by, where the partition function is, Z = Ω|e −TĤ |Ω , andĤ stands for the Hamiltonian of the system. In order to address non-trivial boundary effects induced by the open boundary conditions, fermionic observables are extracted through appropriate ratios that cancel these effects. Examples of such ratios will be provided in section 3.2.

Mass Corrections
The renormalised chiral trajectory is chosen to correspond to a fixed value of φ 4 , in eq. (4), that passes through the physical point. This typically requires the application of small mass corrections to the simulated bare quark masses, that were set via eq. (3) in the dynamical simulations. This small variation of the bare quark masses can be applied to any observable through a low order Taylor expansion [15]. For a generic derived observable defined as a function f Ā 1 , . . . ,Ā n of the mean values of primary observablesĀ i = A i , the expansion with respect to a shift of the mass quarks reads, where the derivative terms read, and S stands for the action.
Notice that for Wilson twisted mass (Wtm) valence fermions only terms proportional to derivative of the action can contribute to the expansion.

Computation of Observables
We use the previously discussed mixed action approach to compute light quark observables involved in a continuum limit scaling study of our lattice formulation. The computation of sea sector observables is relevant for the matching procedure and to compare relative cutoff effects between sea and valence sectors.
We employ the gradient flow quantity t 0 as a relative reference scale. It is defined via the dimensionless observable t 2 E(t, x 0 ) [4,24], where the energy density E(t, x 0 ) is evaluated at a fixed value of the flow time t, The superscript "av" stands for the Euclidean time average over the plateau range in x 0 . Pseudoscalar meson masses and decay constants as well as PCAC quark masses are extracted from the two-point correlation functions, where P and A 0 are the pseudoscalar density and the time component of the improved axial current, A 0 =ψγ 0 γ 5 ψ + ac A ∂ 0ψ γ 5 ψ, respectively. We set the source of the two-point correlation function close to each boundary and symmetrise with respect to x 0 by taking advantage of the time reversal symmetry. Pseudoscalar masses are extracted by averaging over a plateau range the effective mass, am eff (x 0 ) = log f P (x 0 , a) f P (x 0 + a, a) .
The determination of the light PCAC quark masses m 12 is required for the matching procedure. We determine the PCAC masses through the following expression, where∂ 0 stands for the symmetric definition of the numerical derivative with respect to the Euclidean time.
Decay constants can be obtained through the Euclidean time average of a ratio, where X is a fermionic operator. The ratio R X (x 0 ) is constructed in order to cancel boundary effects and to isolate the desired matrix elements from a fit to the plateau region. Specifically, the pseudoscalar meson decay constant of the valence sector can be computed in the following way, By exploiting the PCVC relation, the valence decay constant can be obtained without an explicit use of renormalisation constants. Isospin breaking effects induced by the twisted mass regularisation can be monitored through the difference between charged pion m π ± and neutral pion m π 0 masses. In our mixed action approach, these effects are confined to the valence sector and, therefore, the neutral connected pion mass m π (0,c) already provides a good probe to analyse the size of these O(a 2 ) lattice artefacts. At leading order in Wilson chiral perturbation theory, this mass splitting can be related to the low energy constant w 8 [25], in the following way, m 2 π ± − m 2 π (0,c) = a 2 w 8 .
For various observables considered in our study, a criterion for choosing the Euclidean time plateaux is required in order to remove boundary effects and excited state contributions. Following the proposal of Refs. [15,26], we define the plateau regions by imposing that statistical errors are at least four times larger than the systematic errors from the first excited state at the starting point of the plateau, δ stat (x 0,min ) 4δ syst (x 0,min ). For the ensembles considered in this study, it turned out that for pseudoscalar mesons the excited state contamination can be neglected for x 0 4 √ 8t 0 in the Wilson case, whereas x 0 5 √ 8t 0 is needed with Wilson twisted mass fermions.

Tuning to Maximal Twist
The tuning to maximal twist can be performed by a linear interpolation to a vanishing valence PCAC quark mass, m 12 | v , as a function of κ u | −1 v . This interpolation is illustrated in Figure 1(a). The tuning can thus be achieved through a few simulation points in the neighbourhood of m R 12 v = 0. At lowest order in chiral perturbation theory, the pion mass squared is given by, π v in the neighbourhood of maximal twist. The horizontal yellow band shows the corresponding result from the Wilson formulation. In both panels, the data points correspond to the ensemble H105 with coarsest lattice spacing (see Table 1). The sea and valence quark masses were matched by imposing, µ R 1 v ≡ m R 12 s , as explained in section 2.3.
where m R and µ R are the standard and twisted renormalised quark masses, respectively. For fixed values of µ R , m 2 π is expected to follow a parabolic behaviour in m R in the neighbourhood of its minimum, m R = 0. Figure 1(b) shows an example of such a parabolic behaviour around maximal twist. The horizontal band provides a comparison to the measurement of the pion mass through the Wilson regularisation when using the same statistics as in the Wtm case. Although relative cutoff effects of O(a 2 ) could be non negligible at this coarsest value of the lattice spacing, we observe that the sea and valence quark masses agree within the errors. The fact that the minimum of the parabola coincides with the vanishing of the PCAC mass in Figure 1(a) provides a cross-check of the identification of the maximal twist point.
As indicated in Section 3 and in Table 1, the three ensembles at β = 3.46 lie along the symmetric line which crosses the three flavour massless theory point. By extrapolating our determination of the bare standard mass at maximal twist along the symmetric line, we obtain an estimate of the critical mass, κ cr , which is in the same ballpark as an independent determination coming from the Schrödinger Functional scheme [27] (see also Ref. [28] for a study using large volume CLS ensembles).  [15] where the ensembles were tuned to φ 4 = 1.11, as shown by the grey discontinuous line. Red points are preliminary determinations from our fully twisted tmQCD mixed action setup. The renormalised light quark mass in the valence and in the sea were matched such that, µ R 1 v ≡ m R 12 s .

Continuum Limit
residual O(a) cutoff effects proportional to sea quark masses. The continuum limit scaling is carried out along a line of constant physics composed of symmetric point ensembles, m u = m s , for which φ 4 | s = 1.11 [15]. As mentioned in section 2.3, the matching can been performed by imposing that renormalised quark masses of sea and valence quark masses are equal. As a result, relative cutoff effects of O(a 2 ) are expected between the pseudoscalar meson masses of the sea and valence sectors. From Fig. 2, we observe that these relative cutoff effects are negligible for the lattice spacings, a = 0.077 fm and a = 0.050 fm, while they approximately amount to a 1.5% and a 3% effect at a = 0.065 fm and a = 0.087 fm, respectively. An ongoing complementary study, where the matching of sea and valence quark masses is done through the pion masses, will allow to monitor the continuum scaling of the light quark masses.
In Fig. 3, we present preliminary results for the continuum limit scaling of the pion decay decay constant f π in the the sea [15] and valence sectors. The scaling behaviour is consistent with O(a) improvement for both sea and valence regularisations. Furthermore, the expected agreement of the continuum limit results is verified.

Pion Mass Splitting
In Fig. 4, we compare our preliminary determination of the mass splitting between the charged pion and the connected contribution to the neutral pion with results from different lattice actions [29,30]. The mass difference, (m 2 π ± − m 2 π (0,c) ) r 4 0 /a 2 , in terms of the Sommer scale r 0 , quantifies the scale of the  [15]. Red points are preliminary determinations from our fully twisted tmQCD mixed action setup. The matching between sea and valence has been performed by imposing µ R 1 v ≡ m R 12 s .
O(a 2 ) isospin breaking effects. We observe that our estimates of this mass splitting fall on the same ballpark as those from other lattice formulations [29,30].

Conclusions and Outlook
We define a regularisation based on a mixed action approach, aimed to improve the control of systematic uncertainties in the charm quark sector. The setup combines a fully-twisted valence action matched through the quark masses with the O(a) improved Wilson sea fermionic action. The matching procedure is simplified since the sea and valence actions share the same renormalisation factors. The tuning to maximal twist is performed by linearly interpolating the valence PCAC quark mass to the vanishing point. A first study of the continuum limit scaling of the pion decay constant has shown the agreement of the continuum results with respect to the Wilson formulation. This provides an important validation of this mixed action formulation. We plan to extend the present study to include an alternative matching condition of sea and valence quark masses through the use of pseudoscalar meson masses in order to profit from their automatic O(a) improvement. The next phase of the project will concern the extension of this mixed action setup to the strange and charm sectors.  [29,30] in units of the Sommer scale r 0 . The quantity in the horizontal axis is related to the low energy constant w 8 of Wilson chiral perturbation theory, parametrising the scale of O(a 2 ) lattice artefacts in this pion mass splitting. discussions and for providing valuable input for our analyses. We thankfully acknowledge the computer resources at FinisTerrae II, MareNostrum and Calendula and the technical support provided by CESGA, BSC and FCSCL (FI-2016-3-0022, FI-2017-1-0041, FI-2017-2-0029). We acknowledge support through the Spanish MINECO project FPA2015-68541-P, the Centro de Excelencia Severo Ochoa Programme SEV-2016-0597 and the Ramón y Cajal Programme RYC-2012-10819.