Scalar and vector form factors of $D \to \pi(K) \ell \nu$ decays with $N_f=2+1+1$ twisted fermions

TL;DR

This paper presents a lattice QCD calculation of the vector and scalar form factors for $D o \pi(K) \\ell u$ decays, providing precise theoretical inputs for CKM matrix element extraction and testing CKM unitarity.

Contribution

It offers the first comprehensive lattice determination of these form factors with $N_f=2+1+1$ twisted fermions, including hypercubic effects and continuum extrapolation.

Findings

Calculated form factors across the full kinematic range.

Extracted CKM matrix elements $|V_{cd}|$ and $|V_{cs}|$ with uncertainties.

Confirmed CKM unitarity within current errors.

Abstract

We present a lattice determination of the vector and scalar form factors of the $D \to \pi(K) \ell \nu$ semileptonic decays, which are relevant for the extraction of the CKM matrix elements $|V_{cd}|$ and $|V_{cs}|$ from experimental data. Our analysis is based on the gauge configurations produced by the European Twisted Mass Collaboration with $N_f = 2+1+1$ flavors of dynamical quarks, at three different values of the lattice spacing and with pion masses as small as 210 MeV. The matrix elements of both vector and scalar currents are determined for a plenty of kinematical conditions in which parent and child mesons are either moving or at rest. Lorentz symmetry breaking due to hypercubic effects is clearly observed in the data and included in the decomposition of the current matrix elements in terms of additional form factors. After the extrapolations to the physical pion mass and to the continuum limit we determine the vector and scalar form factors in the whole kinematical region from $q^2 = 0$ up to $q^2_{max} = (M_D - M_{\pi(K)})^2$ accessible in the experiments, obtaining a good overall agreement with experiments, except in the region at high values of $q^2$ where some deviations are visible. A set of synthetic data points, representing our results for $f_+^{D \pi(K)}(q^2)$ and $f_0^{D \pi(K)}(q^2)$ for several selected values of $q^2$, is provided and also the corresponding covariance matrix is available. At zero 4-momentum transfer we get: $f_+^{D \to \pi}(0) = 0.612 ~ (35)$ and $f_+^{D \to K}(0) = 0.765 ~ (31)$. Using the experimental averages for $|V_{cd}| f_+^{D \to \pi}(0)$ and $|V_{cs}| f_+^{D \to K}(0)$, we extract $|V_{cd}| = 0.2330 ~ (137)$ and $|V_{cs}| = 0.945 ~ (38)$, respectively. The second-row of the CKM matrix is found to be in agreement with unitarity within the current uncertainties: $|V_{cd}|^2 + |V_{cs}|^2 + |V_{cb}|^2 = 0.949 ~ (78)$.

Tables

Table 1. Table 1: Summary of the simulated sea and valence quark bare masses, of the π 𝜋 \pi , K 𝐾 K and D 𝐷 D meson masses, of the lattice size L 𝐿 L and of the product M π ​ L subscript 𝑀 𝜋 𝐿 M_{\pi}L for the various gauge ensembles used in this work. The values of M K subscript 𝑀 𝐾 M_{K} and M D subscript 𝑀 𝐷 M_{D} do not correspond to the simulated strange and charm bare quark masses shown in the 5 th superscript 5 th 5^{\rm th} and 6 th superscript 6 th 6^{\rm th} columns, but to the renormalized strange and charm masses interpolated at the physical values m s p ​ h ​ y ​ s ​ ( M ​ S ¯ , 2 ​ GeV ) = 99.6 ​ ( 4.3 ) ​ MeV superscript subscript 𝑚 𝑠 𝑝 ℎ 𝑦 𝑠 ¯ 𝑀 𝑆 2 GeV 99.6 4.3 MeV m_{s}^{phys}(\overline{MS},2~{}\rm{GeV})=99.6(4.3)\,\mathrm{MeV} and m c p ​ h ​ y ​ s ​ ( M ​ S ¯ , 2 ​ GeV ) = 1.176 ​ ( 39 ) ​ GeV superscript subscript 𝑚 𝑐 𝑝 ℎ 𝑦 𝑠 ¯ 𝑀 𝑆 2 GeV 1.176 39 GeV m_{c}^{phys}(\overline{MS},2~{}\rm{GeV})=1.176(39)\,\mathrm{GeV} determined in Ref.  [ 18 ] .

ensemble	$\beta$	$V/a^4$	$a{\mu }_{sea}=a{\mu }_{\ell }$	$a{\mu }_s$	$a{\mu }_c$	$M_{\pi }\left(\mathrm{MeV}\right)$	$M_K\left(\mathrm{MeV}\right)$	$M_D\left(\mathrm{MeV}\right)$	$L\left(\mathrm{fm}\right)$	$M_{\pi }L$
A30.32	$1.90$	${32}^3\times 64$	$0.0030$	$\{0.0180,$	$\{0.21256,$	$275$	$569$	$2015$	$2.84$	$3.96$
A40.32			$0.0040$	$0.0220,$	$0.25000,$	$315$	$578$	$2018$		$4.53$
A50.32			$0.0050$	$0.0260\}$	$0.29404\}$	$351$	$578$	$2018$		$5.04$
A40.24		${24}^3\times 48$	$0.0040$			$324$	$584$	$2024$	$2.13$	$3.49$
A60.24			$0.0060$			$386$	$599$	$2022$		$4.17$
A80.24			$0.0080$			$444$	$619$	$2037$		$4.79$
A100.24			$0.0100$			$495$	$639$	$2042$		$5.34$
B25.32	$1.95$	${32}^3\times 64$	$0.0025$	$\{0.0155,$	$\{0.18705,$	$258$	$545$	$1950$	$2.61$	$3.42$
B35.32			$0.0035$	$0.0190,$	$0.22000,$	$302$	$556$	$1944$		$3.99$
B55.32			$0.0055$	$0.0225\}$	$0.25875\}$	$375$	$578$	$1959$		$4.96$
B75.32			$0.0075$			$436$	$600$	$1965$		$5.77$
B85.24		${24}^3\times 48$	$0.0085$			$467$	$611$	$1974$	$1.96$	$4.63$
D15.48	$2.10$	${48}^3\times 96$	$0.0015$	$\{0.0123,$	$\{0.14454,$	$220$	$526$	$1928$	$2.97$	$3.31$
D20.48			$0.0020$	$0.0150,$	$0.17000,$	$254$	$533$	$1933$		$3.83$
D30.48			$0.0030$	$0.0177\}$	$0.19995\}$	$308$	$547$	$1939$		$4.65$

Table 2. Table 2: Time intervals adopted for the extraction of the PS meson energies E D ​ ( P ) subscript 𝐸 𝐷 𝑃 E_{D(P)} and the matrix elements Z D ​ ( P ) subscript 𝑍 𝐷 𝑃 Z_{D(P)} from the 2-point correlators in the light ( ℓ ℓ \ell ), strange ( s 𝑠 s ) and charm ( c 𝑐 c ) sectors. The last column contains the values of the time distance t ′ superscript 𝑡 ′ t^{\prime} between the source and the sink adopted for the 3-point correlators ( 11 ).

$\beta$	$V/a^4$	${[t_{\min },t_{\max }]}_{(\ell \ell ,\ell s)}/a$	${[t_{\min },t_{\max }]}_{\left(\ell c\right)}/a$	$t^{\prime }/a$
$1.90$	${32}^3\times 64$	$[12,\mathrm{\hspace{0.17em}31}]$	$[8,\mathrm{\hspace{0.17em}16}]$	$18$
	${24}^3\times 48$	$[12,\mathrm{\hspace{0.17em}23}]$	$[8,\mathrm{\hspace{0.17em}17}]$	$18$
$1.95$	${32}^3\times 64$	$[13,\mathrm{\hspace{0.17em}31}]$	$[9,\mathrm{\hspace{0.17em}18}]$	$20$
	${24}^3\times 48$	$[13,\mathrm{\hspace{0.17em}23}]$	$[9,\mathrm{\hspace{0.17em}18}]$	$20$
$2.10$	${48}^3\times 96$	$[18,\mathrm{\hspace{0.17em}40}]$	$[12,\mathrm{\hspace{0.17em}24}]$	$26$

Table 3. Table 3: Values of the parameter θ 𝜃 \theta , appearing in Eq. ( 21 ), for the various ETMC gauge ensembles of Table  1 .

$\beta$	$V/a^4$	$\theta$
$1.90$	${32}^3\times 64$	$0.0,\pm 0.200,\pm 0.467,\pm 0.867$
	${24}^3\times 48$	$0.0,\pm 0.150,\pm 0.350,\pm 0.650$
$1.95$	${32}^3\times 64$	$0.0,\pm 0.183,\pm 0.427,\pm 0.794$
	${24}^3\times 48$	$0.0,\pm 0.138,\pm 0.321,\pm 0.596$
$2.10$	${48}^3\times 96$	$0.0,\pm 0.212,\pm 0.493,\pm 0.916$

Table 4. Table 4: Synthetic data points for the transition D → π → 𝐷 𝜋 D\to\pi representing our results for the vector and scalar form factors extrapolated to the physical pion point and to the continuum and infinite volume limits for eight selected values of q 2 superscript 𝑞 2 q^{2} in the range between q 2 = 0 superscript 𝑞 2 0 q^{2}=0 and q 2 = q m ​ a ​ x 2 = ( M D − M π ) 2 ≃ 3.0 ​ GeV 2 superscript 𝑞 2 superscript subscript 𝑞 𝑚 𝑎 𝑥 2 superscript subscript 𝑀 𝐷 subscript 𝑀 𝜋 2 similar-to-or-equals 3.0 superscript GeV 2 q^{2}=q_{max}^{2}=(M_{D}-M_{\pi})^{2}\simeq 3.0~{}\mathrm{Ge\kern-1.00006ptV}^{2} . The errors correspond to the uncertainties related to (statistical + fitting procedure + input parameters), chiral extrapolation, FSEs and discretization effects, respectively (see text). The errors in squared brackets correspond to the combination in quadrature of the statistical and all systematic errors.

$q^2\left({\text{GeV}}^2\right)$	$f_{+}\left(q^2\right)$	$f_0\left(q^2\right)$
$0.0$	$0.612\left(35\right)\left(4\right)\left(7\right)\left(1\right)\left[35\right]$	$0.612\left(35\right)\left(4\right)\left(7\right)\left(1\right)\left[35\right]$
$0.4286$	$0.715\left(31\right)\left(4\right)\left(6\right)\left(1\right)\left[32\right]$	$0.659\left(29\right)\left(4\right)\left(5\right)\left(1\right)\left[30\right]$
$0.8571$	$0.840\left(29\right)\left(3\right)\left(6\right)\left(1\right)\left[30\right]$	$0.713\left(24\right)\left(3\right)\left(3\right)\left(1\right)\left[24\right]$
$1.2857$	$0.991\left(29\right)\left(4\right)\left(6\right)\left(1\right)\left[30\right]$	$0.773\left(18\right)\left(4\right)\left(2\right)\left(1\right)\left[19\right]$
$1.7143$	$1.179\left(34\right)\left(10\right)\left(3\right)\left(1\right)\left[35\right]$	$0.842\left(15\right)\left(4\right)\left(5\right)\left(1\right)\left[17\right]$
$2.1429$	$1.415\left(43\right)\left(15\right)\left(8\right)\left(1\right)\left[47\right]$	$0.922\left(19\right)\left(5\right)\left(8\right)\left(1\right)\left[21\right]$
$2.5714$	$1.721\left(60\right)\left(21\right)\left(16\right)\left(1\right)\left[66\right]$	$1.017\left(29\right)\left(7\right)\left(13\right)\left(1\right)\left[32\right]$
$3.0000$	$2.130\left(86\right)\left(31\right)\left(27\right)\left(3\right)\left[96\right]$	$1.134\left(45\right)\left(10\right)\left(18\right)\left(1\right)\left[49\right]$

Table 5. Table 5: Synthetic data points for the D → K → 𝐷 𝐾 D\to K transition representing our results for the vector and scalar form factors extrapolated to the physical pion point and in the continuum and infinite volume limits for eight selected values of q 2 superscript 𝑞 2 q^{2} in the range between q 2 = 0 superscript 𝑞 2 0 q^{2}=0 and q 2 = q m ​ a ​ x 2 = ( M D − M K ) 2 ≃ 1.88 ​ GeV 2 superscript 𝑞 2 superscript subscript 𝑞 𝑚 𝑎 𝑥 2 superscript subscript 𝑀 𝐷 subscript 𝑀 𝐾 2 similar-to-or-equals 1.88 superscript GeV 2 q^{2}=q_{max}^{2}=(M_{D}-M_{K})^{2}\simeq 1.88~{}\mathrm{Ge\kern-1.00006ptV}^{2} . The errors correspond to the uncertainties related to (statistical + fitting procedure + input parameters), chiral extrapolation and discretization effects, respectively (see text). The errors in squared brackets correspond to the combination in quadrature of the statistical and all systematic errors.

$q^2\left({\text{GeV}}^2\right)$	$f_{+}\left(q^2\right)$	$f_0\left(q^2\right)$
$0.0$	$0.765\left(29\right)\left(11\right)\left(1\right)\left[31\right]$	$0.765\left(29\right)\left(11\right)\left(1\right)\left[31\right]$
$0.2692$	$0.815\left(28\right)\left(12\right)\left(1\right)\left[31\right]$	$0.792\left(26\right)\left(10\right)\left(1\right)\left[28\right]$
$0.5385$	$0.872\left(28\right)\left(13\right)\left(1\right)\left[31\right]$	$0.820\left(23\right)\left(10\right)\left(1\right)\left[25\right]$
$0.8077$	$0.937\left(28\right)\left(15\right)\left(1\right)\left[32\right]$	$0.849\left(21\right)\left(9\right)\left(1\right)\left[23\right]$
$1.0769$	$1.013\left(29\right)\left(17\right)\left(1\right)\left[34\right]$	$0.879\left(19\right)\left(9\right)\left(1\right)\left[21\right]$
$1.3461$	$1.102\left(32\right)\left(21\right)\left(1\right)\left[38\right]$	$0.911\left(17\right)\left(9\right)\left(1\right)\left[19\right]$
$1.6154$	$1.208\left(36\right)\left(26\right)\left(1\right)\left[44\right]$	$0.944\left(17\right)\left(8\right)\left(1\right)\left[19\right]$
$1.8846$	$1.336\left(43\right)\left(32\right)\left(1\right)\left[54\right]$	$0.979\left(17\right)\left(8\right)\left(1\right)\left[19\right]$

Table 6. Table 6: Values of the parameters appearing in the z-expansions of the vector and scalar form factors ( 68 - 69 ) in the case of the D → π → 𝐷 𝜋 D\to\pi transition.

$f^{D\to \pi }\left(0\right)$	$c_{+}^{D\to \pi }$	$P_V\left({\text{GeV}}^{-2}\right)$	$c_0^{D\to \pi }$	$P_S\left({\text{GeV}}^{-2}\right)$
$0.6117\left(354\right)$	$-1.985\left(347\right)$	$0.1314\left(127\right)$	$-1.188\left(256\right)$	$0.0342\left(122\right)$

Table 7. Table 7: Covariance matrix corresponding to the z-expansions of the vector and scalar form factors (A1)-(A2) in the case of the D → π → 𝐷 𝜋 D\to\pi transition.

	$f^{D\to \pi }\left(0\right)$	$c_0^{D\to \pi }$	$c_{+}^{D\to \pi }$	$P_S$	$P_V$
$f^{D\to \pi }\left(0\right)$	$1.25642\cdot {10}^{-3}$	$7.18296\cdot {10}^{-3}$	$6.77051\cdot {10}^{-3}$	$3.66997\cdot {10}^{-5}$	$2.87257\cdot {10}^{-5}$
$c_0^{D\to \pi }$	$7.18296\cdot {10}^{-3}$	$6.56690\cdot {10}^{-2}$	$6.30124\cdot {10}^{-2}$	$1.73569\cdot {10}^{-3}$	$8.43689\cdot {10}^{-4}$
$c_{+}^{D\to \pi }$	$6.77051\cdot {10}^{-3}$	$6.30124\cdot {10}^{-2}$	$1.20371\cdot {10}^{-1}$	$2.24220\cdot {10}^{-3}$	$3.25631\cdot {10}^{-3}$
$P_S$	$3.66997\cdot {10}^{-5}$	$1.73569\cdot {10}^{-3}$	$2.24220\cdot {10}^{-3}$	$1.48010\cdot {10}^{-4}$	$9.60595\cdot {10}^{-5}$
$P_V$	$2.87257\cdot {10}^{-5}$	$8.43689\cdot {10}^{-4}$	$3.25631\cdot {10}^{-3}$	$9.60595\cdot {10}^{-5}$	$1.60179\cdot {10}^{-4}$

Table 8. Table 8: Values of the parameters appearing in the z-expansions of the vector and scalar form factors ( 70 - 71 ) in the case of the D → K → 𝐷 𝐾 D\to K transition.

$f^{D\to K}\left(0\right)$	$c_{+}^{D\to K}$	$c_0^{D\to K}$
$0.7647\left(308\right)$	$-0.066\left(333\right)$	$-2.084\left(283\right)$

Table 9. Table 9: Covariance matrix corresponding to the z-expansions of the vector and scalar form factors ( 70 - 71 ) in the case of the D → K → 𝐷 𝐾 D\to K transition.

	$f^{D\to K}\left(0\right)$	$c_0^{D\to K}$	$c_{+}^{D\to K}$
$f^{D\to K}\left(0\right)$	$9.50493\cdot {10}^{-4}$	$6.92027\cdot {10}^{-3}$	$5.66397\cdot {10}^{-3}$
$c_0^{D\to K}$	$6.92027\cdot {10}^{-3}$	$7.99358\cdot {10}^{-2}$	$7.55735\cdot {10}^{-2}$
$c_{+}^{D\to K}$	$5.66397\cdot {10}^{-3}$	$7.55735\cdot {10}^{-2}$	$1.10925\cdot {10}^{-1}$

Equations

$$f_{+}^{D\to \pi }(0)=0.612(35),f_{+}^{D\to K}(0)=0.765(31).$$

$$|V_{cd}|$$

$$|V_{cs}|$$

$$|V_{cd}{|}^2+|V_{cs}{|}^2+|V_{cb}{|}^2=0.949(78),$$

$$⟨{\hat{V}}_{\mu }⟩\equiv ⟨P(p_P)⟩{\hat{V}}_{\mu }D(p_D)=(p_D+p_P{)}_{\mu }{f}_{+}(q^2)+(p_D-p_P{)}_{\mu }{f}_{-}(q^2),$$

$$f_0(q^2)=f_{+}(q^2)+\frac{q^2}{M_D^2-M_P^2}{f}_{-}(q^2),$$

$$⟨\hat{S}⟩\equiv ⟨P(p_P)⟩\hat{S}D(p_D)=\frac{M_D^2-M_P^2}{m_c-m_x}{f}_0(q^2),$$

$$C_2^{D(P)}\left(t^{\prime },{\vec{p}}_{D(P)}\right)$$

$$C_{\hat{\Gamma }}^{DP}\left(t,t^{\prime },{\vec{p}}_D,{\vec{p}}_P\right)$$

$$C_2^{D(P)}\left(t^{\prime },{\vec{p}}_{D(P)}\right)$$

$$C_{\hat{\Gamma }}^{DP}\left(t,t^{\prime },{\vec{p}}_D,{\vec{p}}_P\right)$$

$${\hat{V}}_{\mu }$$

$$\hat{S}$$

$$⟨S⟩\equiv ⟨P(p_P)⟩SD(p_D)=\frac{M_D^2-M_P^2}{{\mu }_c-{\mu }_q}{f}_0(q^2),$$

$$R_{\mu }(t,{\vec{p}}_D,{\vec{p}}_P)$$

$$R_S(t,{\vec{p}}_D,{\vec{p}}_P)$$

$${\tilde{C}}_2^{D(P)}\left(t,{\vec{p}}_{D(P)}\right)\equiv \frac{1}{2}\left[C_2^{D(P)}\left(t,{\vec{p}}_{D(P)}\right)+\sqrt{C_2^{D(P)}{\left(t,{\vec{p}}_{D(P)}\right)}^2-C_2^{D(P)}{\left(\frac{T}{2},{\vec{p}}_{D(P)}\right)}^2}\right],$$

$${\tilde{C}}_2^{D(P)}\left(t,{\vec{p}}_{D(P)}\right){}_{\overrightarrow{t\gg a}}\frac{Z_{D(P)}}{2E_{D(P)}}{e}^{-E_{D(P)}t},$$

$$R_{\mu }(t,{\vec{p}}_D,{\vec{p}}_P)$$

$$R_S(t,{\vec{p}}_D,{\vec{p}}_P)$$

$$p=\frac{2\pi }{L}\left(\theta +n_x,\theta +n_y,\theta +n_z\right),$$

$${⟨{\hat{V}}_0⟩}_{\mathrm{imp}}$$

$${⟨{\hat{V}}_i⟩}_{\mathrm{imp}}$$

$${⟨S⟩}_{\mathrm{imp}}$$

$${⟨{\hat{V}}_{\mathrm{sp}}⟩}_{\mathrm{imp}}\equiv \frac{1}{3}\left[{⟨{\hat{V}}_1⟩}_{\mathrm{imp}}+{⟨{\hat{V}}_2⟩}_{\mathrm{imp}}+{⟨{\hat{V}}_3⟩}_{\mathrm{imp}}\right].$$

$$⟨P(p_P)⟩{\hat{V}}_{\mu }D(p_D)$$

$$⟨P(p_P)⟩SD(p_D)$$

$${⟨{\hat{V}}_0⟩}_{\mathrm{imp}}$$

$${⟨{\hat{V}}_{\mathrm{sp}}⟩}_{\mathrm{imp}}$$

$${⟨S⟩}_{\mathrm{imp}}$$

$$q_{\mu }^E=\left(\vec{q},q_4\right)=\left(\vec{q},i{q}_0\right)$$

$$q^{[n]}{P}^{[m]}\equiv \sum _{\mu }{\left(q_{\mu }^E\right)}^n{\left(P_{\mu }^E\right)}^m.$$

$$⟨P(p_P)⟩{\hat{V}}_{\mu }^{E}D(p_D)={⟨{\hat{V}}_{\mu }^E⟩}_{\mathrm{Lor}}+{⟨{\hat{V}}_{\mu }^E⟩}_{\mathrm{hyp}},$$

$${⟨{\hat{V}}_{\mu }^E⟩}_{\mathrm{Lor}}=P_{\mu }^E{f}_{+}(q^2,a^2)+q_{\mu }^E\frac{M_D^2-M_P^2}{q^2}\left[f_0(q^2,a^2)-f_{+}(q^2,a^2)\right],$$

$${⟨{\hat{V}}_{\mu }^E⟩}_{\mathrm{hyp}}=a^2\left[{\left(q_{\mu }^E\right)}^3{H}_1+{\left(q_{\mu }^E\right)}^2{P}_{\mu }^E{H}_2+q_{\mu }^E{\left(P_{\mu }^E\right)}^2{H}_3+{\left(P_{\mu }^E\right)}^3{H}_4\right]$$

$$H_i(z)=d_0^i+d_1^{i}z+d_2^i{z}^2,$$

Full text

1 Introduction

The Cabibbo-Kobayashi-Maskawa (CKM) matrix [1, 2] describes the quark flavor mixing in the electroweak sector of the Standard Model (SM). It represents the only source of CP violation in the quark sector and the origin of its hierarchical structure, which describes so precisely the relative strength of the flavor-changing weak quark currents, is still unexplained. In the SM the CKM matrix is unitary, so we have that $V_{\mathrm{CKM}}\times V^{\dagger}_{\mathrm{CKM}}=\mathbbm{1}$. This gives rise to unitarity conditions between the elements of its rows and columns, that are represented by diagonal constraints and “unitarity triangles”, whose investigation has been the focus of much of the experimental and theoretical efforts in Flavor Physics during the recent years. On the one hand side inconsistencies in the CKM-picture would indicate the presence of new physics beyond the SM. On the other hand side, if all the precision tests of the SM performed so far are in agreement with the CKM paradigm, the absence of deviations provides stringent constraints on nonstandard phenomena and their energy scale. It is therefore important to determine all CKM matrix elements as precisely as possible by studying flavour-changing processes both experimentally and theoretically.

The golden modes for testing the unitarity of the second row of the CKM matrix, namely $\lvert V_{cd}\rvert^{2}+\lvert V_{cs}\rvert^{2}+\lvert V_{cb}\rvert^{2}=1$, are represented by the leptonic and semileptonic decays of charmed $D$ and $D_{s}$ mesons, which probe the $c\to d$ and $c\to s$ quark transitions, respectively. Combining the experimental measurements of the branching fractions of these processes with the theoretical calculations of the relevant hadronic matrix elements, i.e. the leptonic decay constants $f_{D}$ and $f_{D_{s}}$ and the semileptonic vector form factors at zero four-momentum transfer $f_{+}^{D\to\pi}(0)$ and $f_{+}^{D\to K}(0)$, the CKM entries $\lvert V_{cd}\rvert$ and $\lvert V_{cs}\rvert$ can be determined. The CKM matrix element $\lvert V_{cb}\rvert$, being of the order ${\cal{O}}(10^{-2})$, is irrelevant for the second-row unitarity test at the current level of precision. According to the $V-A$ structure of the SM weak currents, leptonic and semileptonic $D$ and $D_{s}$ decays should provide consistent results for the CKM elements $\lvert V_{cd}\rvert$ and $\lvert V_{cs}\rvert$.

As is well known, the theoretical calculations of hadronic matrix elements based entirely on first principles can be properly carried out by simulating the fundamental theory of the strong interaction, QCD, on a lattice. Thanks to the remarkable progress in algorithms and computing machines Lattice QCD (LQCD) has entered the precision era and the accuracy of numerical computations is becoming comparable to that of experiments. For some relevant hadronic quantities in Flavour Physics the goal of the per cent level of precision has been already achieved (see, e.g., the FLAG review [3]).

In this work we present the first $N_{f}=2+1+1$ LQCD calculation of the vector and scalar form factors $f_{+}^{D\pi(K)}(q^{2})$ and $f_{0}^{D\pi(K)}(q^{2})$ governing the semileptonic $D\to\pi(K)\ell\nu$ decays, using the gauge configurations generated by the European Twisted Mass Collaboration (ETMC) with $N_{f}=2+1+1$ dynamical quarks, which include in the sea, besides two light mass-degenerate quarks, also the strange and charm quarks with masses close to their physical values [4, 5]. The same ensembles have already been used to determine the leptonic decay constants $f_{D}$ and $f_{D_{s}}$, which have provided our leptonic determinations of $\lvert V_{cd}\rvert$ and $\lvert V_{cs}\rvert$ [6], using the experimental values for $\lvert V_{cd}\rvert f_{D}$ and $\lvert V_{cs}\rvert f_{D_{s}}$.

At variance with most of the existing LQCD calculations (see, e.g., Ref. [3]), that provide only the value of the vector form factor at zero 4-momentum transfer, we have evaluated both the vector and scalar form factors in the whole accessible range of values of $q^{2}$ in the experiments, i.e. from $q^{2}=0$ up to $q^{2}_{\rm max}=(M_{D}-M_{\pi(K)})^{2}$. In our calculations quark momenta are injected on the lattice using non-periodic boundary conditions [7, 8] and the matrix elements of both vector and scalar currents are determined for a plenty of kinematical conditions, in which parent and child mesons are either moving or at rest.

The data coming from different kinematical conditions exhibit a remarkable breaking of Lorentz symmetry due to hypercubic effects (see Refs. [9, 10] for preliminary results) for both $D\to\pi$ and $D\to K$ semileptonic form factors. We show also evidence that the hypercubic artifacts may be largely affected by the difference between the parent and the child meson masses. This may represent an important warning in the case of the determination of the form factors governing semileptonic $B$-meson decays.

We present the subtraction of the hypercubic artifacts and the determination of the Lorentz-invariant semileptonic vector and scalar form factors after the combined extrapolations to the physical pion mass and to the continuum limit. A set of synthetic data points, representing our results for $f_{+}^{D\pi(K)}(q^{2})$ and $f_{0}^{D\pi(K)}(q^{2})$ for several selected values of $q^{2}$, is provided (see Tables 4 and 5) and the corresponding covariance matrix is available upon request. Our results are compared with the momentum dependence of the experimental data from BELLE [11], BABAR [12, 13], CLEO [14] and BESIII [15], and a good overall agreement is obtained. Some deviations are nevertheless visible at high values of $q^{2}$ for both $D\to\pi\ell\nu$ and $D\to K\ell\nu$ decays.

At zero 4-momentum transfer the results of our study are

[TABLE]

Using the updated experimental values for $\lvert V_{cd}\rvert f_{+}^{D\to\pi}(0)=0.1426(19)$ and $\lvert V_{cs}\rvert f_{+}^{D\to K}(0)=0.7226(34)$, obtained by the Heavy Flavor Averaging Group (HFAG) in Ref. [16], the following values of the CKM matrix elements $|V_{cd}|$ and $|V_{cs}|$ are derived:

[TABLE]

where the errors are from the lattice calculation and from experiment, respectively, showing that the dominant error is the theoretical one.

Including the determination of $|V_{cb}|$ from $B$-meson decays [17], the unitarity test of the second-row of the CKM matrix yields

[TABLE]

which confirms unitarity at the level of several percent.

The physics described by our simulations and used throughout this work corresponds to the isospin symmetric limit of QCD, where $m_{u}=m_{d}=m_{ud}$, and the quark electric charges are neglected. Therefore, isospin breaking and electromagnetic corrections have to be added separately in phenomenological analyses.

The paper is organized as follows. In Section 2 we describe the simulation details. In Section 3 we present the computation of the vector and scalar form factors $f^{D\pi(K)}_{+}(q^{2})$ and $f^{D\pi(K)}_{0}(q^{2})$ using the matrix elements of the weak vector and scalar quark currents relevant for the $D\to\pi(K)$ transition, obtained from 2- and 3-point correlation functions. In Section 4 the evidence of Lorentz symmetry breaking in the momentum dependence of the form factors is presented and discussed. In Section 5 we describe the strategy adopted in order to extract the physical, Lorentz invariant, vector and scalar form factors. This is based on a global fit of all the data corresponding to the time and spatial components of the vector current and to the scalar current for all the lattice ensembles, studying simultaneously the dependence on $q^{2}$, the light-quark mass $m_{\ell}$ and the lattice spacing $a$, and using a phenomenological Ansatz to describe the hypercubic effects. In Section 6 the results for $f_{+}^{D\pi(K)}(q^{2})$ and $f_{0}^{D\pi(K)}(q^{2})$ from the global fit are shown and compared with experimental data. The extraction of the CKM matrix elements $|V_{cd}|$ and $|V_{cs}|$ is also presented together with the test the unitarity of the second row of the CKM matrix. Finally our conclusions are summarized in Sec. 7.

2 Simulation details

The gauge ensembles used in this work have been generated by the ETMC with $N_{f}=2+1+1$ dynamical quarks, which include in the sea, besides two light mass-degenerate quarks, also the strange and the charm quarks [4, 5]. The ensembles are the same adopted in Refs. [18] to determine the up, down, strange and charm quark masses, using the experimental value of the pion decay constant $f_{\pi}$ to set the lattice scale111With respect to Ref. [18] the number of independent gauge configurations adopted for the ensemble D15.48 has been increased to $90$ to improve the statistics.. They have been used also to determine the leptonic decay constants $f_{K}$, $f_{D}$ and $f_{D_{s}}$ in Ref. [6] and the vector and scalar form factors of the semileptonic $K\to\pi\ell\nu$ decay in Ref. [19].

The gauge fields are simulated using the Iwasaki gluon action [20], while sea quarks are implemented with the Wilson Twisted Mass Action at maximal twist [21, 22, 23]. In order to avoid the mixing of strange and charm quarks in the valence sector we have adopted the non-unitary setup described in Ref. [24], in which the valence strange quarks are regularized as Osterwalder-Seiler (OS) fermions [25], while the valence up and down quarks have the same action of the sea. The use of different lattice regularisations for the valence and sea heavy quarks avoids completely the effects of the strange and charm mixing without modifying the renormalization pattern of operators in massless schemes and produces only a modification of discretization effects. Moreover, since we work at maximal twist, an automatic ${\cal{O}}(a)$-improvement [23, 24] is guaranteed also for our non-unitary setup.

The QCD simulations have been carried out at three different values of the inverse bare lattice coupling $\beta$, to allow for a controlled extrapolation to the continuum limit, and at different lattice volumes. For each gauge ensemble we have used a number of gauge configurations corresponding to a separation of 20 trajectories to avoid autocorrelations. We have simulated quark masses in the range from $\simeq 3\,m_{ud}$ to $\simeq 12\,m_{ud}$ in the light sector, from $\simeq 0.7\,m_{s}$ to $\simeq 1.2\,m_{s}$ in the strange sector, and from $\simeq 0.7\,m_{c}$ to $\simeq 1.2\,m_{c}$ in the charm sector, where $m_{ud}$, $m_{s}$ and $m_{c}$ are the physical values of the average up/down, strange and charm quark masses respectively, as determined in Ref. [18]. The lattice spacings are found to be $a=\{0.0885\,(36),0.0815\,(30),0.0619\,(18)\}\,\mathrm{fm}$ at $\beta=\{1.90,1.95,2.10\}$ respectively, the lattice volume goes from $\simeq 2$ to $\simeq 3$ fm and the pion masses, extrapolated to the continuum and infinite volume limits, range from $\simeq 210$ to $\simeq 450\,\,\mathrm{MeV}$ (see Ref. [18] for further details).

In our study of the semileptonic $D\to\pi\ell\nu$ and $D\to K\ell\nu$ form factors we make use of the input parameters (values of quark masses and lattice spacings) obtained from the eight branches of the analysis carried out in Ref. [18]. The various branches differ by:   i) the choice of the scaling variable, which was taken to be either the Sommer parameter $r_{0}/a$ [26] or the mass of a fictitious pseudoscalar (PS) meson made of two strange-like quarks $aM_{s^{\prime}s^{\prime}}$;   ii) the fitting procedures, which were based either on Chiral Perturbation Theory (ChPT) or on a polynomial expansion in the light quark mass (for the motivations see the discussion in Section 3.1 of Ref. [18]);   iii) the choice between two methods, denoted as M1 and M2 which differ by ${\cal{O}}(a^{2})$ effects (see, e.g., Ref. [27]), used to determine non-perturbatively the values of the mass renormalization constant (RC) ${\cal{Z}}_{m}=1/{\cal{Z}}_{P}$ [18]. Throughout this work the results corresponding to the various branches of the analysis are combined to form our averages and errors according to Eq. (28) of Ref. [18].

The basic simulation parameters and the masses of the $\pi$, $K$ and $D$ mesons corresponding to each ensemble used in this work are collected in Table 1.

3 Lattice calculation of the vector and scalar matrix elements

The matrix element of the weak vector current $\widehat{V}_{\mu}$ between an initial $D$-meson state and a $\pi$($K$)-meson final state decomposes, as required by the Lorentz symmetry, into the two form factors $f_{+}(q^{2})$ and $f_{-}(q^{2})$:

[TABLE]

where $P=\pi(K)$ can be either the pion or the kaon and the 4-momentum transfer $q$ is given by $q\equiv p_{D}-p_{P}$. The scalar form factor $f_{0}$ is then defined as

[TABLE]

so that the kinematic identity $f_{+}(0)=f_{0}(0)$ is satisfied by definition. The scalar form factor is proportional to the 4-divergence of $\braket{\widehat{V}_{\mu}}$ so that, thanks to the Ward-Takahashi identity (WTI), $f_{0}$ can be determined from the matrix element of the scalar density $\widehat{S}$ between the $D$-meson and the $\pi$($K$)-meson states:

[TABLE]

where $x=\ell(s)$.

From Eqs. (5) and (7) it turns out that the determination of the scalar and vector form factors can be carried out by computing the matrix elements $\braket{\widehat{V}_{\mu}}$ and $\braket{\widehat{S}}$. These two quantities can be extracted from the large (Euclidean) time distance behavior of a convenient combination of 2- and 3-point correlation functions in lattice QCD, which are defined as

[TABLE]

where $t^{\prime}$ is the time distance between the source and the sink, $t$ is the time distance between the insertion of the current $\widehat{\Gamma}=\{\widehat{V}_{\mu},~{}\widehat{S}\}$ and the source, while $P_{5}^{D}=i\,\bar{c}\gamma_{5}u$ and $P_{5}^{\pi(K)}=i\,\bar{d}(\bar{s})\gamma_{5}u$ are the interpolating fields of the $D$ and $\pi(K)$ mesons. The Wilson parameters $r$ of the two valence quarks in any PS meson are always chosen to have opposite values, $r_{c}=r_{s}=r_{d}=-r_{u}$, so that the squared PS mass differs from its continuum counterpart only by terms of ${\cal{O}}(a^{2}m\Lambda_{QCD})$ [23].

The statistical accuracy of the correlators (8-9) is significantly improved by using the all-to-all quark propagators evaluated with the so-called “one-end” stochastic method [28], which includes spatial stochastic sources at a single time slice chosen randomly (see Ref. [29] where the degenerate case of the electromagnetic pion form factor is discussed in details). Statistical errors on the quantities directly extracted from the correlators are always evaluated using the jackknife procedure, while cross-correlations are taken into account by the use of the eight bootstrap samplings (with ${\cal{O}}(100)$ events each) corresponding to the eight analyses of Ref. [18] (see Section 2).

In the case of charm quarks it is a common procedure to adopt Gaussian-smeared interpolating fields [30] for both the source and the sink in order to suppress faster the contribution of the excited states, leading to an improved projection onto the ground state at relatively small time distances. For the values of the smearing parameters we set $k_{G}=4$ and $N_{G}=30$. In addition, we apply APE-smearing to the gauge links [31] in the interpolating fields with parameters $\alpha_{APE}=0.5$ and $N_{APE}=20$. The Gaussian smearing is applied as well also for the light and strange quarks. The values of $M_{\pi}$ and $M_{K}$ reported in Table 1 are consistent within the statistical errors with the corresponding ones listed in the Table 1 of Ref. [19], computed using local interpolating fields.

As is well known, at large time distances 2- and 3-point correlation functions behave as

[TABLE]

where $Z_{D}$ and $Z_{P}$ are the matrix elements $\braket{0\lvert\,P_{5}^{D}(0)\,\rvert\,D(\vec{p}_{D})}$ and $\braket{0\lvert\,P_{5}^{P}(0)\,\rvert\,P(\vec{p}_{P})}$, which depend on the meson momenta $\vec{p}_{D}$ and $\vec{p}_{P}$ because of the use of smeared interpolating fields, while $E_{D}$ and $E_{P}$ are the energies of the $D$ and $P$ mesons. These energies and matrix elements can be extracted directly by fitting the exponential behavior, given by the r.h.s of Eq. (10), of the corresponding 2-point correlation functions. The time intervals $[t_{\rm min},\,t_{\rm max}]$ adopted for the fit (10) are listed in Table 2. We have checked that the extracted values of $E_{D(P)}$ are nicely reproduced (within the statistical errors) by the continuum-like dispersive relation $E_{D(P)}^{\rm\,disp}=\sqrt{M_{D(P)}^{2}+|\vec{p}_{D(P)}|^{2}}$, where $M_{D(P)}$ is the meson mass extracted from the 2-point correlator corresponding to the meson at rest. We decided to use for the analysis the energy values $E_{D(P)}^{\rm\,disp}$ instead of those directly extracted from the fit.

As for the 3-point correlators (11) the usual choice of the time distance $t^{\prime}$ between the source and the sink is to maximize it, i.e. to put $t^{\prime}=T/2$. Since we are using smeared interpolating fields, it is convenient to choose smaller values of $t^{\prime}$, which allow to decrease significantly the statistical noise satisfying at the same time the dominance of the ground-state signals. We have optimized the choice of the values of $t^{\prime}$ for the various lattice spacings and volumes, which can be read off in the last column of Table 2.

In our lattice setup we employ maximally twisted fermions and therefore the vector and scalar currents, $\widehat{V}_{\mu}$ and $\widehat{S}$, renormalize multiplicatively [23]. Introducing the local bare currents $V_{\mu}=\bar{c}\gamma_{\mu}q$ and $S=\bar{c}q$, where $q=d(s)$, and keeping the same value for the Wilson parameters of the initial and final quarks, one has

[TABLE]

where ${\cal{Z}}_{V}$ and ${\cal{Z}}_{P}$ are the renormalization constants (RCs) of the vector and pseudoscalar densities for standard Wilson fermions, respectively. The twisted quark masses renormalize multiplicatively with a RC ${\cal{Z}}_{m}$ given by ${\cal{Z}}_{m}=1/{\cal{Z}}_{P}$ [23], which means that the product $(m_{c}-m_{q})\braket{\widehat{S}}$ is renormalization group invariant. Therefore, according to Eq. (7), the scalar form factor $f_{0}(q^{2})$ is related to the (bare) matrix element $\braket{S}$ by

[TABLE]

where $\mu_{q}$ and $\mu_{c}$ are the bare light (strange) and charm quark masses, respectively.

As in the case of the semileptonic $K_{\ell 3}$ decay [19], the matrix elements $\langle S\rangle$ and $\langle\widehat{V}_{\mu}\rangle$ (see Eq. (5)) can be extracted from the time dependence of the ratios $R_{\mu}$ $(\mu=0,1,2,3)$ and $R_{S}$ of 2- and 3-point correlation functions, defined as in Eqs. (8-9) but using the local bare currents $V_{\mu}$ and $S$, namely

[TABLE]

where the correlation functions $\widetilde{C}_{2}^{D(P)}(t)$ are given by

[TABLE]

which at large time distances behave as

[TABLE]

i.e. without the backward signal. Note that the denominator of Eq. (15) is nothing but the numerator evaluated in the mass-degenerate limit for the valence quarks in the current. Such mass-degenerate quarks have the same lattice regularisation of the corresponding ones in the numerator, so that the RC ${\cal{Z}}_{V}$ is the same for the two terms in the ratio and therefore it cancels out.

At large time distances one has:

[TABLE]

so that, up to lattice artifacts, the renormalized quantity $|\braket{\widehat{V}_{\mu}}|^{2}$ and the bare one $|\braket{S}|^{2}$ can be extracted directly from the plateau of $R_{\mu}$ and $R_{S}$, independently of the matrix elements $Z_{D}$ and $Z_{\pi(K)}$ of the interpolating fields. In the r.h.s. of Eq. (19) we have used the fact that, due to the charge conservation, $\braket{P(p_{P})}{\widehat{V}_{\mu}}{P(p_{P})}=2\,p_{P\mu}$ and $\braket{D(p_{D})}{\widehat{V}_{\mu}}{D(p_{D})}=2\,p_{D\mu}$. Taking the square root of $R_{\mu}$ and $R_{S}$ we can get the absolute value of the matrix elements $\braket{\widehat{V}_{\mu}}$ and $\braket{S}$, while their sign can be easily inferred from those of the correlators $C_{V_{\mu}}^{DP}(t,t^{\prime},\vec{p}_{D},\vec{p}_{P})$ and $C_{S}^{DP}(t,t^{\prime},\vec{p}_{D},\vec{p}_{P})$ in the relevant time regions.

When a spatial component of the momentum of either the parent or the child meson is vanishing, the corresponding matrix element $\braket{P(p_{P})}{V_{i}}{P(p_{P})}$ (or $\braket{D(p_{D})}{V_{i}}{D(p_{D})}$) ($i=1,2,3$) is also vanishing and the correlator $C^{PP}_{V_{i}}$ (or $C^{DD}_{V_{i}}$), cannot be used in the denominator of the r.h.s of Eq. (15). In these cases the quantity $2p_{Pi}\,/\!\braket{P(p_{P})}{V_{i}}{P(p_{P})}$ (or $2p_{Di}\,/\!\braket{D(p_{D})}{V_{i}}{D(p_{D})}$) is replaced by $2E_{P}\,/C^{PP}_{V_{0}}$ (or $2E_{D}\,/C^{DD}_{V_{0}}$).

When both the pion (kaon) and the D-meson are at rest, only two ratios, $R_{S}$ and $R_{0}$, can be constructed, providing two independent determinations of the scalar form factor $f_{0}(q^{2})$ at the kinematical end-point $q^{2}_{\rm max}=\left(M_{D}-M_{P}\right)^{2}$, which may differ by lattice artifacts.

In order to inject momenta on the lattice we use the same procedure adopted in Ref. [19] for the $K_{\ell 3}$ decays. In particular for the valence quark fields we impose twisted boundary conditions (BC’s) [7, 8, 32] in the spatial directions and anti-periodic BC’s in time. The sea dynamical quarks, on the contrary, have been simulated with periodic BC’s in space and anti-periodic ones in time. The choice of using twisted BC’s for the valence quark fields is crucial in order to remove the strong limitations, resulting from the use of periodic BC’s, to the accessible kinematical regions of momentum-dependent quantities like, e.g., form factors. Furthermore we remark that, as shown in Refs. [33, 34], for physical quantities which do not involve final state interactions (like, e.g., meson masses, decay constants and semileptonic form factors), the use of different BC’s for valence and sea quarks produces only finite size effects (FSEs) which are exponentially small.

The quark 3-momentum is then given by

[TABLE]

where $n_{x,y,z}$ are integers and the $\theta$ values, adopted for the different gauge ensembles and democratically distributed along the three spatial directions, are collected in Table 3. They have been chosen in order to obtain momenta with values ranging from $\approx 150\,\mathrm{MeV}$ to $\approx 650\,\mathrm{MeV}$ for all the various lattice spacings and volumes222The correlators used in this work have been calculated within the PRACE project PRA067 “First Lattice QCD study of B-physics with four flavors of dynamical quarks”. The values of the quark momentum were not chosen having in mind the investigation of hypercubic effects in the semileptonic form factors. In particular the use of spatially symmetric values of the quark momentum (see Table 3) is not ideal for such a purpose..

The 3-point correlation functions $C^{DP}_{V_{\mu}}(t,t^{\prime},\vec{p}_{D},\vec{p}_{P})$ and $C^{DP}_{S}(t,t^{\prime},\vec{p}_{D},\vec{p}_{P})$ have been simulated imposing periodic BC’s to the spectator light quark and partially twisted BC’s (21) to the initial $c$ and final $u(s)$ quarks. With this choice the $\pi$, $K$ and $D$ meson (spatial) momenta are given by $\vec{p}_{D}=(2\pi/L)\left(\theta_{D},\theta_{D},\theta_{D}\right)$ and $\vec{p}_{\pi(K)}=(2\pi/L)\left(\theta_{\pi(K)},\theta_{\pi(K)},\theta_{\pi(K)}\right)$, where $\theta_{D}$ and $\theta_{\pi(K)}$ can assume for each gauge ensemble the values of the parameter $\theta$ given in Table 3.

As described in Ref. [23] the use of two kinematics with opposite spatial momenta of the initial and final mesons, given by opposite signs of the corresponding $\theta$, allows for an ${\cal{O}}(a)$ improvement on the matrix elements $\braket{\widehat{V}_{\mu}}$ and $\braket{S}$ performing the following average:

[TABLE]

Furthermore since we are using democratically distributed momenta in the three spatial directions, the matrix elements of the spatial components of the vector current $\braket{\widehat{V}_{i}}_{\rm imp}$ are equal to each other. Therefore, in order to improve the statistics, we average them to get

[TABLE]

The quality of the plateau for the matrix elements $\braket{\widehat{V}_{0}}_{\rm imp}$, $\braket{\widehat{V}_{\rm sp}}_{\rm imp}$ and $\braket{S}_{\rm imp}$ is illustrated in Fig. 1 in the case of the $D\to\pi$ transition. The time intervals adopted for fitting Eqs. (19-20) are symmetric around $t^{\prime}/2$ (see Table 2 for the values of $t^{\prime}$ for each specific gauge ensemble) and equal to $[t^{\prime}/2-2,~{}t^{\prime}/2+2]$. These values are compatible with the dominance of the $\pi$, $K$ and $D$ mesons ground-state observed along the time intervals of Table 2 for the two-point correlation functions.

Thus, from the 2- and 3-point lattice correlators we are able to extract the three ${\cal{O}}(a)$-improved matrix elements $\braket{\widehat{V}_{0}}_{\rm imp}$, $\braket{\widehat{V}_{\rm sp}}_{\rm imp}$ and $\braket{S}_{\rm imp}$. The standard procedure for determining the scalar and vector form factors $f_{0}(q^{2})$ and $f_{+}(q^{2})$ is to assume the following Lorentz-covariant decomposition

[TABLE]

with $P_{\mu}\equiv(p_{D}+p_{P})_{\mu}$. Explicitly one has

[TABLE]

which represent a redundant mathematical system consisting of three equations depending on just two form factors.

We then determine $f_{0}(q^{2})$ and $f_{+}(q^{2})$ by minimizing the $\chi^{2}$-variable constructed using Eqs. (28-30). In the next Section we present and discuss the result of this determination in which, as anticipated, we found evidence for Lorentz symmetry breaking terms.

4 Form factors and hypercubic effects

After a small interpolation of our lattice data to the physical values of the strange and charm quark masses, $m_{s}^{phys}(2\,\rm{GeV})=99.6\,(4.3)\,$MeV and $m_{c}^{phys}(2\,\rm{GeV})=1.176\,(3.9)\,$GeV taken from Ref. [18], we determine the vector and scalar form factors $f_{+}^{D\pi(K)}$ and $f_{0}^{D\pi(K)}$ for each gauge ensemble and for each choice of parent and child meson momenta. The momentum dependencies of the semileptonic form factors $f_{+}^{D\pi(K)}$ and $f_{0}^{D\pi(K)}$ are illustrated in the upper (lower) panels of Fig. 2, where different markers and colors correspond to different values of the child meson momentum for the ETMC ensemble A60.24, i.e. at fixed values of the parent and child meson masses as well as of the lattice spacing and volume. Therefore, if the Lorentz-covariant decomposition (26-27) were adequate to describe our data, the extracted form factors would depend only on the squared 4-momentum transfer $q^{2}$. This is not the case and an extra dependence on the value of the child meson momentum is clearly visible in Fig. 2 beyond the statistical uncertainties.

As is well known, the lattice breaks Lorentz symmetry and is invariant only under discrete rotations by multiple of $90^{\circ}$ in each direction of the Euclidean space-time. Therefore, the matrix elements (22-25), and consequently the form factors extracted from Eqs. (28-30), may depend also on hypercubic invariants333Hypercubic symmetry is also broken on our lattices because of the different length of the temporal and spatial dimensions. This effect however is expected to be subdominant and will be neglected in what follows.. Hypercubic effects are known to affect lattice calculations and they have been discussed for instance in Refs. [35, 36]. It is however the first time that these effects are observed in the $D$-meson semileptonic form factors. In Refs. [37] and [38, 39] $N_{f}=2+1$ results for the $f_{+}^{D\pi(K)}$ form factor have been obtained by FNAL/MILC and HPQCD collaborations, respectively, using only a limited number of kinematic conditions, restricted in particular to the parent $D$-meson at rest. Also the ETMC reported $N_{f}=2$ results for the $D$-meson semileptonic form factors in Ref. [40], but the kinematics were limited to the Breit-frame ($\vec{p}_{D}=-\vec{p}_{P}$). Recently in Ref. [41] both the $D\to\pi$ and $D\to K$ semileptonic transitions have been investigated using $N_{f}=2+1$ domain-wall fermions assuming the $D$-meson at rest, while in Ref. [42] the FNAL/MILC collaboration has addressed the determination of the semileptonic form factors, using $N_{f}=2+1+1$ MILC ensembles with HISQ fermions and tuning properly the child meson momentum to reach directly $q^{2}=0$, but working only in the reference frame where the $D$-meson is at rest. We argue that all these choices may obscure the presence of hypercubic effects in the lattice data.

The behavior observed in Fig. 2 might be (at least partially) related to finite size effects (FSEs) (see, e.g., Ref. [43] for the case of $K_{\ell 3}$ decays). The possible impact of FSEs has been investigated by comparing the results corresponding to the two gauge ensembles, A40.24 and A40.32, which share the same pion mass and lattice spacing, but have different lattice sizes, $L=24a$ and $L=32a$, respectively, as illustrated in Fig. 3.

It can clearly be seen that for the $D\to\pi(K)$ semileptonic vector and scalar form factors FSEs appear to be negligible within the current statistical uncertainties, except for the slope of the $D\to\pi$ scalar form factor (upper right panel in Fig. 3). Hypercubic effects for the two gauge ensembles A40.24 and A40.32 are found to be comparable and they do not appear to depend on the lattice size $L$.

Thus, since in our setup all the current matrix elements are ${\cal{O}}(a)$-improved, the breaking of the Lorentz invariance is expected to be produced by ${\cal{O}}(a^{2})$ hypercubic effects, whose subtraction will be discussed in the next Section in order to get the Lorentz-invariant semileptonic vector and scalar form factors.

Before closing the Section, it is worth noting that no evidence of hypercubic effects within the current statistical uncertainties was found in the case of the $K\to\pi\ell\nu$ semileptonic form factors analyzed in Ref. [19], where the same gauge configurations and the same parent and child momenta were adopted444In Ref. [19] the vector and scalar form factors for the $K_{\ell 3}$ decays have been constructed using local interpolating fields for both the pion and the kaon. We have checked that no hypercubic effects are visible also in the case of smeared interpolating fields.. This suggests that the hypercubic artifacts may be governed by the difference between the parent and the child meson masses. Such an indication is confirmed by the results given in Fig. 4, where the transition between two charmed PS mesons with masses close to the $D$-meson one has been considered. The momentum dependencies of the corresponding form factors show no evidence of hypercubic effects within the statistical uncertainties.

The dependence of hypercubic artifacts upon the mass difference between the parent and the child mesons is clearly a very important issue, which warrants further investigations. It may represent an important warning in the case of the determination of the form factors governing semileptonic $B$-meson decays into lighter mesons.

5 Subtraction of the hypercubic effects

As shown in the previous Section the form factors $f_{+}$ and $f_{0}$ for both the $D\to\pi$ and $D\to K$ decays exhibit a sizeable Lorentz-symmetry breaking due to hypercubic effects generated at finite lattice spacing. Let’s start by introducing Euclidean 4-momenta defined, in the case of the 4-momentum transfer $q_{\mu}=(q_{0},~{}\vec{q})$, as

[TABLE]

so that $\sum_{\mu}q_{\mu}^{E}q_{\mu}^{E}=-q^{2}$. Hypercubic invariants can be constructed using the two momenta $q_{\mu}^{E}$ and $P_{\mu}^{E}=(p_{D}+p_{P})_{\mu}^{E}$ as

[TABLE]

For $n+m=2$ there are three invariants, which are Lorentz invariants, namely $q^{2}$, $q\cdot P$ and $P^{2}$. The three invariants can be rewritten in terms of $q^{2}$ and the parent and child meson masses. For $n+m=4$ there are five hypercubic invariants: $q^{[4]}$, $q^{[3]}P^{[1]}$, $q^{[2]}P^{[2]}$, $q^{[1]}P^{[3]}$, $P^{[4]}$, where $q^{[4]}$ ($P^{[4]}$) stands for $q^{[4]}P^{[0]}$ ($q^{[0]}P^{[4]}$).

As already shown in Fig. 2, the form factors $f_{+,0}$ calculated using Eqs. (28-30) do not depend on $q^{2}$ only. A possible way to describe the observed hypercubic effects is to address them directly on the vector and scalar matrix elements. We start by considering the following decomposition of the vector matrix elements:

[TABLE]

where $\braket{\widehat{V}_{\mu}^{E}}_{\rm Lor}$ is the Lorentz-covariant term

[TABLE]

while $\braket{\widehat{V}_{\mu}^{E}}_{\rm hyp}$ is given by

[TABLE]

and the quantities $H_{i}$ ($i=1,2,3,4$) are additional hypercubic form factors. Note that in the Lorentz-covariant term (34) we have explicitly considered that the form factors $f_{+,0}$ can be affected by discretization errors of order $\mathcal{O}(a^{2})$, which are unrelated to hypercubic effects and may depend on $q^{2}$ as well as on the parent and child meson masses.

Eq. (35) is the most general structure, up to order $\mathcal{O}(a^{2})$, that transforms properly under hypercubic rotations and is built with third powers of the components of the two momenta $q_{\mu}^{E}$ and $P_{\mu}^{E}$. The Lorentz-invariance breaking effects are encoded in the four structures $\left(q_{\mu}^{E}\right)^{3}$, $\left(q_{\mu}^{E}\right)^{2}P_{\mu}^{E}$, $q_{\mu}^{E}\left(P_{\mu}^{E}\right)^{2}$ and $\left(P_{\mu}^{E}\right)^{3}$ as well as in the hypercubic form factors $H_{i}$, which, we assume, depend only on $q^{2}$ and the parent and child meson masses. Note that the decomposition (33-35) implies that the form factors $f_{+,0}$ calculated using Eqs. (28-30) do depend not only on $q^{2}$, but also on the five hypercubic invariants $q^{[4]}$, $q^{[3]}P^{[1]}$, $q^{[2]}P^{[2]}$, $q^{[1]}P^{[3]}$, $P^{[4]}$.

For the $H_{i}$ form factors we adopt a simple polynomial form in terms of the $z$ variable [44, 45]

[TABLE]

where $z$ is defined as

[TABLE]

with $t_{+}$ and $t_{0}$ given by

[TABLE]

In Eq. (36) the coefficients $d_{0,1,2}^{i}$ will be treated as free parameters.

Let’s now turn to the matrix elements of the scalar density. We consider the following decomposition of the scalar density matrix elements:

[TABLE]

where $\braket{S}_{\rm Lor}$ is the Lorentz-invariant term

[TABLE]

while $\braket{S}_{\rm hyp}$ is given by

[TABLE]

and the quantities $\widetilde{H}_{i}$ ($i=1,2,3,4,5$) are additional hypercubic form factors.

The Ward-Takahashi Identity (WTI) relates the 4-divergence of the vector current to the scalar density. Let’s introduce the WTI breaking term $\Delta_{WTI}^{\rm hyp}$ defined as

[TABLE]

which implies

[TABLE]

The quantity $\Delta_{WTI}^{\rm hyp}$ can be evaluated directly using the matrix elements $\braket{\widehat{V}_{\mu}^{E}}$ and $\braket{S}$. Its dependence on the parent child momenta is illustrated in Fig. 5 in the case of the gauge ensemble A30.32.

It can be seen that the WTI-violating term $\Delta_{WTI}^{\rm hyp}$ is small, but not vanishing and it cannot depend only on the Lorentz-invariant $q^{2}$. Instead the lattice data suggest a simple, approximate linear dependence on the hypercubic invariant $q^{[4]}$ (see right panels in Fig. 5). This implies a quite simplified structure for $\Delta_{WTI}^{\rm hyp}$ in Eq. (43) and consequently for $\braket{S}^{\rm hyp}$ in Eq. (41), namely: $\widetilde{H}_{i}+H_{i}=0$ for $i=2,3,4$, $\widetilde{H}_{5}=0$ and $\widetilde{H}_{1}+H_{1}=H_{S}\neq 0$. In other words one has $\Delta_{WTI}^{\rm hyp}=a^{2}q^{[4]}H_{S}$ and

[TABLE]

where for the hypercubic form factor $H_{S}$ we adopt the simple Ansatz

[TABLE]

with $d_{0,1}^{S}$ being free parameters.

The structure of the hypercubic artifacts is thus given by Eqs. (35) and (44) in terms of five hypercubic form factors $H_{i}$ ($i=1,2,3,4,S$). The latter cannot be determined by analyzing the matrix elements of vector and scalar currents calculated separately for each gauge ensemble in the present work. A simultaneous, global fit of all the data (more than one thousand lattice points corresponding to the time and spatial components of the vector and scalar currents for the 15 ETMC gauge ensembles of Table 1) has to be performed by considering the dependencies on $q^{2}$, $m_{\ell}$ and $a^{2}$ of the form factors $f_{+,0}$ as well as the $q^{2}$ and $m_{\ell}$ dependencies of the five hypercubic form factors $H_{i}$ ($i=1,2,3,4,S$).

For the form factors $f_{+,0}(q^{2},a^{2})$ we have adopted the modified z-expansion of Ref. [46], viz.

[TABLE]

where we assume for the coefficients $c_{+(0)}^{D\to\pi(K)}(a^{2})$ a simple linear dependence on $a^{2}$ and $z_{0}\equiv z(q^{2}=0$), so that the condition $f_{+}(0,a^{2})=f_{0}(0,a^{2})=f(0,a^{2})$ is explicitly fulfilled at finite lattice spacing. In the r.h.s. of Eqs. (46-47) the terms at second order in the $z$-variable are constrained by the analyticity requirements described in Ref. [46].

As for the functions $\mathcal{P}_{+,0}$, in the case of $D\to\pi$ transition we adopt the single-pole expressions

[TABLE]

while for the $D\to K$ channel we use

[TABLE]

In the case of the $D\to\pi$ pole factors (48) and (49) the quantities $M_{V}$ and $M_{S}$ represent the vector and scalar pole masses, respectively. They are treated as free parameters in the fitting procedure. In the case of the $D\to K$ decays the data are fitted equally well even excluding the pole term in the scalar form factor and therefore we choose $\mathcal{P}_{0}^{D\to K}(q^{2})=1$. Conversely the physical vector meson $D_{s}^{*}$ has a mass below the cut threshold $\sqrt{t_{+}}=(M_{D_{s}}+M_{K})$. Consequently the pole factor (50), including a simple discretization effect proportional to $a^{2}$, is introduced to guarantee the applicability of the $z$-expansion.

In Eq. (48) the quantity $K_{FSE}^{0}(L)$ takes into account the FSEs observed in Fig. 3 by adopting the following phenomenological form

[TABLE]

where $C_{FSE}^{0}$ is a free parameter and $\xi_{\ell}=2Bm_{\ell}/(16\pi^{2}f^{2})$, with $B$ and $f$ being the SU(2) low-energy constants entering the LO chiral Lagrangian and determined in Ref. [18].

For the vector form factor at zero 4-momentum transfer, $f^{D\to\pi(K)}(0,a^{2})$, we use the following Ansatz

[TABLE]

where the coefficients $F_{+}$, $b_{1}$, $b_{2}$ and $D$ are treated as free parameters in the fitting procedure, while $A^{\pi(K)}$ is the chiral-log coefficient predicted by the hard pion SU(2) Chiral Perturbation Theory (ChPT) [47], given by

[TABLE]

where for the coupling constant $\widehat{g}$ we adopt the value $\widehat{g}=0.61$ [17].

Using the ingredients described above we have performed the global, combined fit of all the data for the matrix elements $\braket{\widehat{V}_{0}}$, $\braket{\widehat{V}_{sp}}$ and $\braket{S}$, which amount to a total of $1110$ data points for both the $D\to\pi$ and $D\to K$ transitions. The total number of free parameters is $24$ ($19$) in the case of $D\to\pi(K)$ channel, namely:

• •

in Eq. (53) 3 parameters for $D\to\pi$ ($F_{+}$, $b_{1}$ and $D$) and 4 parameters for $D\to K$ ($F_{+}$, $b_{1}$, $b_{2}$ and $D$);

• •

in Eq. (46) 2 parameters for $c_{+}^{D\to\pi(K)}$ (i.e., $c_{+}=c_{+}^{0}+c_{+}^{1}a^{2}$);

• •

in Eq. (47) 2 parameters for $c_{0}^{D\to\pi(K)}$ (i.e., $c_{0}=c_{0}^{0}+c_{0}^{1}a^{2}$);

• •

1 parameter ($M_{V}$) in Eq. (48) for $D\to\pi$ and 1 parameter ($P_{+}$) in Eq. (50) for $D\to K$;

• •

2 parameters in Eqs. (49) and (52) ($M_{S}$ and $C_{FSE}^{0}$) only for $D\to\pi$,;

• •

2 parameters in Eq. (45) for the hypercubic form factor $H_{S}$;

• •

in Eq. (36) 3 parameters for each of the four hypercubic form factors $H_{1}$, $H_{2}$, $H_{3}$ and $H_{4}$ for $D\to\pi$ and 2 parameters (i.e., $d_{2}^{i}=0$ ) for $D\to K$.

The quality of the fit is quite good obtaining $\chi^{2}/\rm{d.o.f.}\simeq 1.2$ for both the $D\to\pi$ and $D\to K$ transitions. We have tried to include extra terms in Eqs. (46-47) either proportional to $z^{2}$ (including the analyticity requirement of Ref. [46] through an appropriate term proportional to $z^{3}$) or proportional to the light-quark mass $m_{\ell}$ in the coefficients $c_{+,0}^{D\to\pi(K)}$. Since the differences in the results for both the hypercubic corrections and the form factors $f_{+,0}$ are negligible with respect to the other errors and the values of the new parameters turn out to be consistent with [math], such extended fits are not used for estimating systematic uncertainties.

From the global, combined fit we obtain both the momentum dependence of the Lorentz-invariant form factors $f_{+,0}$ and of the five hypercubic form factors $H_{i}$ ($i=1,2,3,4,S$). The momentum dependence of $f_{+,0}$, extrapolated to the physical pion mass and to the continuum and infinite volume limits, will be discussed and compared to the experimental data in Section 6. Here we compute the hypercubic form factors $H_{i}$ coming from the global fit in order to check the quality of the subtraction of the hypercubic effects for each gauge ensemble, i.e. at finite lattice spacing and volume and for the unphysical pion masses given in Table 1.

In Fig. 6 we show the same form factors given in Fig. 2 after the hypercubic contributions determined by the global fit have been subtracted from the matrix elements $\braket{\widehat{V}_{\mu}}$ and $\braket{S}$ using Eqs. (33) and (39). It can be seen that the hypercubic effects are properly removed and both the scalar and the vector form factors depend now only on the $4-$momentum transfer $q^{2}$ within the statistical uncertainties.

In the limiting case where the parent and the child mesons are the same, Eq. (35) reduces to a simpler expression, namely

[TABLE]

because only even terms under the exchange of the initial and final PS mesons survive. In Section 4 we have noted that within the statistical uncertainties there is no evidence of hypercubic effects when the initial and final meson have the same masses (see Fig. 4). This might be an indication that the hypercubic form factors $H_{2}$ and $H_{4}$ can be neglected. Thus, we have repeated the global fitting procedure assuming $H_{2}=H_{4}=0$, which reduces the number of free parameters to $18$ and $15$ for the $D\to\pi$ and $D\to K$ transitions, respectively.

The differences in the results for both the hypercubic corrections and the form factors $f_{+,0}$, obtained including ($H_{2}\neq H_{4}\neq 0$) or excluding ($H_{2}=H_{4}=0$) the two hypercubic form factors $H_{2}$ and $H_{4}$, are found to be negligible within the current statistical uncertainties. Therefore, in what follows we adopt the fitting procedure, in which we assume $H_{2}=H_{4}=0$, as our reference fit for estimating the uncertainties due to various sources of systematic errors as well as for obtaining the results for the form factors $f_{+,0}^{D\to\pi(K)}(q^{2})$, extrapolated to the physical pion mass and to the continuum and infinite volume limits, which will be discussed and compared to the experimental data in Section 6.

We stress again that an important feature of our analysis with respect to previous studies of the semileptonic $D\to\pi(K)$ form factors is the use of a plenty of kinematical conditions corresponding to parent and child mesons either moving or at rest. Using only a limited number of kinematical conditions, for instance the Breit-frame in which $\vec{p_{D}}=-\vec{p}_{\pi(K)}$ or the $D-$meson at rest, the presence of the hypercubic effects may not be manifest.

This point is illustrated in Fig. 7, which shows the subset of our data for the scalar $D\to\pi$ (left panel) and $D\to K$ (right panel) form factor $f_{0}$ corresponding only to the D-meson at rest both before and after the subtraction of the hypercubic effects determined in the global fitting procedure. Lorentz-symmetry breaking is not manifest in the limited set of data points with $\vec{p}_{D}=0$, but it is not negligible. This holds for the scalar form factor $f_{0}$, while in the case of the vector form factor $f_{+}$ we find that Lorentz-symmetry breaking effects are less pronounced in the subset of data corresponding to the D-meson at rest. We stress that the differences between the data with and without hypercubic effects are a ${\cal{O}}(a^{2})$ effect proportional to hypercubic invariants. Thus, any analysis of the data without the subtraction of hypercubic effects, based directly on parameterizations like Eqs. (46-47), where only discretization effects unrelated to hypercubic invariants are considered, is in principle inadequate and may lead to different results in the continuum limit.

6 Results from the global fit and comparison with experimental data

The momentum dependencies of the physical Lorentz-invariant vector and scalar form factors, extrapolated to the physical pion mass and to the continuum and infinite volume limits, are shown in Fig. 8 for both the $D\to\pi$ and $D\to K$ transitions. Our results exhibit a remarkable precision in the full range of values of $q^{2}$ covered by the experiments (i.e., $0\leq q^{2}\leq q^{2}_{\rm{max}}=(M_{D}-M_{\pi(K)})^{2}\simeq 3.0(1.9)$ GeV2).

Our results for the vector form factors $f^{D\pi}_{+}(q^{2})$ and $f^{DK}_{+}(q^{2})$ can be compared with the corresponding values determined by BELLE, BABAR, CLEO and BESIII collaborations in Refs. [11, 12, 13, 14, 15], where the partial decay rates have been measured (see also Refs. [48, 49] for a summary of the experimental results). The agreement is good except at high values of $q^{2}$, where some deviations are visible.

In Fig. 9 our main results for the vector and scalar form factors are compared with those obtained by choosing only the kinematical configurations corresponding to the $D$-meson rest frame and by performing the extrapolations to the physical pion mass and to the continuum and infinite volume limits without including the hypercubic terms (35) and (41). In other words, the continuum extrapolation is based only on the discretization terms contained in Eqs. (46-47). It can be seen that the neglect of hypercubic effects in the analysis and the use of a limited subset of data lead to some distortions of the extrapolated form factors, which are more pronounced in the case of the scalar form factor. Such distortions are found to be comparable with present global uncertainties within one standard-deviation. They may become more relevant as the precision of the data will be increased in the future.

In Table 4 we provide a set of synthetic data points for the vector and scalar $D\to\pi$ form factors, $f^{D\pi}_{+}(q^{2})$ and $f^{D\pi}_{0}(q^{2})$, with the corresponding total uncertainties, calculated at eight selected values of $q^{2}$ between [math] and $q^{2}_{\rm{max}}=(M_{D}-M_{\pi})^{2}$.

The errors in Table 4 take into account the uncertainties induced by:

• •

the statistical noise and the fitting procedure;

• •

the errors in the determinations of the input parameters, namely the values of the average $u/d$ quark mass $m_{ud}$ the value of the charm quark mass $m_{c}$, the lattice spacing $a$ and the $SU(2)$ ChPT LECs $f$ and $B_{0}$, determined in Ref. [18];

• •

the chiral extrapolation, evaluated combining the results obtained using the SU(2) ChPT fit on all our lattice data and a fit with $b_{2}=0$ in Eq. (53) applied only to the data with $M_{\pi}<390$ MeV555In this case the total number of data values reduces to $814$, since the results for the gauge ensembles A80.24, A100.24, B85.24 and B75.32 (see Table 1) are excluded from the analysis.;

• •

the FSEs, evaluated by comparing the results obtained with and without the FSE factor (52);

• •

the discretization effects, calculated by comparing our main results with those obtained including in Eqs. (46-47) extra terms proportional to $(a\Lambda_{\rm{QCD}})^{4}$. Using a value for $\Lambda_{QCD}$ equal to $\simeq 0.35$ GeV, we expect that the values of the coefficients of the extra terms are in a natural range of order ${\cal{O}}(1)$. Therefore, we adopt for the coefficients of the extra terms a (conservative) prior distribution equal to $0\pm 3$.

Similarly, in Table 5 we provide a set of synthetic data points for the vector and scalar $D\to K$ form factors, $f^{DK}_{+}(q^{2})$ and $f^{DK}_{0}(q^{2})$, with the corresponding total uncertainties for eight selected values of $q^{2}$ between [math] and $q^{2}_{\rm{max}}=(M_{D}-M_{K})^{2}$. Note that, at variance with the case of the $D\to\pi$ transition, the uncertainty related to FSEs is not considered, because the data for the $D\to K$ transition do not show any visible FSE.

In order to allow a direct use of the synthetic data points without using our bootstrap samples, we have calculated the covariance matrix among the synthetic data points contained either in Table 4 or in Table 5. The corresponding covariance matrices are available upon request to allow to fit our synthetic data with any functional form, that can be adopted for describing the momentum dependence of the form factors.

For a direct use of our results for the form factors $f_{+,0}^{D\to\pi(K)}$ we provide in the Appendix the values of the parameters of the z-expansions of our global fit after the extrapolations to the physical pion point and to the continuum and infinite volume limits, including the corresponding covariance matrices.

From Tables 4 and 5 our results at zero 4-momentum transfer are

[TABLE]

which are consistent within the errors with the FLAG [3] averages $f_{+}^{D\to\pi}(0)=0.666~{}(29)$, based on the result of Ref. [38], and $f_{+}^{D\to K}(0)=0.747~{}(19)$ from Ref. [39].

Using the experimental values

[TABLE]

given by HFAG in Ref. [16], we obtain our results for the CKM matrix elements $|V_{cd}|$ and $|V_{cs}|$:

[TABLE]

where the errors are from the lattice calculation and from the experiments respectively, showing that the dominant error is the theoretical one. Our results (58) can be compared with the determinations of $|V_{cd}|$ and $|V_{cs}|$ based on the $D$ and $D_{s}$ leptonic decay constants, $f_{D}=207.4(3.8)$ MeV and $f_{D_{s}}=247.2(4.1)$ MeV, obtained in Ref. [6] using the same ETMC gauge configurations. Using the experimental values of $f_{D}|V_{cd}|=46.06(1.11)$ MeV and $f_{D_{s}}|V_{cs}|=250.66(4.48)$ MeV from Ref. [50] one gets

[TABLE]

where again the errors are from the lattice calculation and from the experiments, respectively. At variance with the semileptonic case, the theoretical uncertainties of $|V_{cd}|$ and $|V_{cs}|$, obtained from the leptonic decays, are comparable to (or even smaller than) the experimental ones.

An alternative way to extract the CKM matrix elements $|V_{cd}|$ and $|V_{cs}|$ is to combine directly the momentum dependence of the semileptonic form factors obtained from lattice QCD simulations with the experimental $q^{2}$-bins of the differential $D\to\pi(K)\ell\nu_{\ell}$ decay rates. The application of such a strategy using the form factors determined in this work will be presented elsewhere [Riggio:2017zwh].

The results (58-61) are presented in Fig. 10 as ellipses in the $|V_{cd}|-|V_{cs}|$ plane corresponding to a $68\%$ probability contour. In Fig. 10 also the ellipses corresponding to the leptonic and semileptonic FLAG averages [3] for $|V_{cd}|$ and $|V_{cs}|$ are shown, as well as the constraint imposed by the second-row unitarity, indicated by a dashed line.

Using $|V_{cb}|=0.0360(9)$ from Ref. [17] we can perform the check of the unitarity of the second row of the CKM matrix. We find

[TABLE]

which test the second-row unitarity at the level of several percent for both semileptonic and leptonic modes.

7 Conclusions

We have presented the first lattice $N_{f}=2+1+1$ determination of the vector and scalar form factors of the $D\to\pi\ell\nu$ and $D\to K\ell\nu$ semileptonic decays, which are relevant for the extraction of the CKM matrix elements $\lvert V_{cd}\rvert$ and $\lvert V_{cs}\rvert$ from experimental data.

Our analysis is based on the gauge configurations produced by ETMC with $N_{f}=2+1+1$ flavors of dynamical quarks at three different values of the lattice spacing with pion masses as small as 210 MeV. Quark momenta are injected on the lattice using non-periodic boundary conditions. The matrix elements of both vector and scalar currents are determined for a plenty of kinematical conditions in which parent and child mesons are either moving or at rest.

Lorentz symmetry breaking due to hypercubic effects is clearly observed in the data and included in the decomposition of the current matrix elements in terms of additional form factors. We found evidence that the hypercubic artifacts may be governed by the difference between the parent and the child meson masses. This represents a very important issue, which warrants further investigations, since it might become particularly relevant in the case of the determination of the form factors governing semileptonic $B$-meson decays into lighter mesons.

In the present work the values of the quark momentum were not chosen having in mind the investigation of hypercubic effects in the semileptonic form factors. In particular the use of spatially symmetric values of the quark momentum is not ideal for such a purpose. We have planned to perform new calculations of the semileptonic form factors removing the above constraint and optimizing the choice of the non-periodic boundary conditions. Nevertheless, we stress that the main structure of the hypercubic effects on the matrix elements of the vector and scalar currents has been understood in the present work.

After the extrapolations to the physical pion mass and to the continuum limit we determine the vector and scalar form factors in the whole kinematical region from $q^{2}=0$ up to $q^{2}_{\rm max}=(M_{D}-M_{\pi(K)})^{2}$ accessible in the experiments, obtaining a good overall agreement with experiments. Some deviations are visible at high values of $q^{2}$ for both $D\to\pi\ell\nu$ and $D\to K\ell\nu$ decays. A set of synthetic data points, representing our results for $f_{+}^{D\pi(K)}(q^{2})$ and $f_{0}^{D\pi(K)}(q^{2})$ for several selected values of $q^{2}$, is provided, including also the covariance matrix for the data at different values of $q^{2}$.

At zero 4-momentum transfer we get

[TABLE]

Using the experimental values (57) for $\lvert V_{cd}\rvert f_{+}^{D\to\pi}(0)$ and $\lvert V_{cs}\rvert f_{+}^{D\to K}(0)$ from HFAG [16] we determine

[TABLE]

Including also the determination of $|V_{cb}|$ from $B$-meson decays [17], the test of the second-row of the CKM matrix is

[TABLE]

Acknowledgements

We warmly thank N. Carrasco and F. Sanfilippo for their valuable contribution to the initial stages of the present work. We thank our colleagues of the ETMC for fruitful discussions. We gratefully acknowledge the CPU time provided by PRACE under the project PRA067 “First Lattice QCD study of B-physics with four flavors of dynamical quarks” and by CINECA under the specific initiative INFN-LQCD123 on the BG/Q system Fermi at CINECA (Italy). V.L., S.S. and C.T. thank MIUR (Italy) for partial support under the contract PRIN 2010-2011 and PRIN 2015.

Appendix: The z-expansion of the physical vector and scalar form factors

After the extrapolations to the physical pion point and to the continuum and infinite volume limits, the z-expansions of the vector and scalar form factors, adopted in this work, are written in the case of the $D\to\pi$ transition as

[TABLE]

The values of the five parameters $f^{D\to\pi}(0)$, $c_{+}^{D\to\pi}$, $P_{V}$, $c_{0}^{D\to\pi}$ and $P_{S}$ are collected in Table 6, with the corresponding covariance matrix given in Table 7.

Analogously in the case of the $D\to K$ transition the z-expansions of the vector and scalar form factors read as

[TABLE]

where the values of the three parameters $f^{D\to K}(0)$, $c_{+}^{D\to K}$ and $c_{0}^{D\to K}$ are collected in Table 8, with the corresponding covariance matrix given in Table 9.