The Adler $D$-function for ${\cal N}=1$ SQCD regularized by higher covariant derivatives in the three-loop approximation

We calculate the Adler $D$-function for ${\cal N}=1$ SQCD in the three-loop approximation using the higher covariant derivative regularization and the NSVZ-like subtraction scheme. The recently formulated all-order relation between the Adler function and the anomalous dimension of the matter superfields defined in terms of the bare coupling constant is first considered and generalized to the case of an arbitrary representation for the chiral matter superfields. The correctness of this all-order relation is explicitly verified at the three-loop level. The special renormalization scheme in which this all-order relation remains valid for the $D$-function and the anomalous dimension defined in terms of the renormalized coupling constant is constructed in the case of using the higher derivative regularization. The analytic expression for the Adler function for ${\cal N}=1$ SQCD is found in this scheme to the order $O(\alpha_s^2)$. The problem of scheme-dependence of the $D$-function and the NSVZ-like equation is briefly discussed.


Introduction
The Adler D-function [1] is closely related to the normalized cross-section of electron-positron annihilation to hadrons, namely to the ratio R(s) = σ(e + e − → hadrons) σ 0 (e + e − → µ + µ − ) = 12πIm Π(s), where σ 0 (e + e − → µ + µ − ) = 4πα 2 * /3s with α * being the QED fine-structure constant. Π(s) is the perturbatively evaluated expression for the hadronic vacuum polarization function, which plays an important role in studies of contributions of the strong interaction to various physical quantities. In particular, it is needed for determining strong interaction contributions to the theoretical expression for the muon anomalous magnetic moment (for the most recent discussion see [2]). The D-function is connected with R(s) by the dispersion relation (2) and can be used to compare the theoretical QCD predictions with the available experimental data for R(s) [3]. In the region where perturbation theory is applicable a theoretical expression for the D-function is defined by a sum of diagrams with two external lines of the Abelian gauge field, whereas the quarks and gluons are propagating in the internal loops. Massless perturbative theoretical expressions for the D-function are known in several gauge models. In QCD it was analytically evaluated at the three-loop level (to the order O(α 2 s )) in Refs. [4,5] and numerically in Ref. [6]. At the level O(α 3 s ) the D-function has been analytically calculated in [7]. This result was confirmed in Refs. [8,9]. At present, the QCD analytical expression for the D-function is known to the order O(α 4 s ) [10,11]. In Ref. [12] the α 2 s corrections to the D-function were analytically evaluated in a theoretical model of strong interactions which contains QCD supplemented with multiplets of coloured scalar fields.
A more consistent theoretical model of strong interactions containing coloured scalar fields is N = 1 SQCD, which is continuing to attract attention of both theoreticians and experimentalists (for a review see, e.g., [13]). In this model the O(α s )-corrections to the D-function have been evaluated in [14,15].
An all-loop exact formula relating the D-function of N = 1 SQCD to the anomalous dimension of the matter superfields has been proposed and proved in [16,17] in the case of using the higher covariant derivative regularization [18,19,20,21]. This equation is valid for the D-function and the anomalous dimension defined in terms of the bare coupling constant and has the form where N f flavours of the matter superfields with electric charges q α lie in the fundamental representation of SU (N ). The proof was based on the method used in [22] for the derivation of the NSVZ β-function [23,24,25,26] defined in terms of the bare coupling constant in the N = 1 supersymmetric quantum electrodynamics regularized by higher covariant derivatives by direct summation of Feynman supergraphs. The key idea of this method is based on the following observation made in [27,28]: For supersymmetric gauge theories regularized by higher derivatives, quantum contributions to the two-point Green function of the gauge superfield factorize into integrals of total and even double total derivatives in the momentum space in the limit of the vanishing external momentum. This allows performing explicit integration with respect to one loop momentum in a given contribution to the β-function in the n-loop approximation and relating it to a contribution to the anomalous dimension in the (n − 1)-loop approximation. The idea was fully realized in two different ways in Refs. [22,29], where it was proved to all loops that the β-function is indeed given by integrals of double total derivatives with respect to a loop momentum. Due to this property, it is related to the anomalous dimension of the matter superfields through the NSVZ relation if both renormalization group (RG) functions are defined in terms of the bare coupling constant. 1 These results were verified by explicit three-loop calculations in Ref. [31]. Note that with dimensional reduction the factorization into double total derivatives does not take place [32], and the RG functions defined in terms of the bare coupling constant do not satisfy the NSVZ relation [33]. The formula (3) itself is an analogue of the NSVZ β-function for the N = 1 SYM theory with matter, where tr( The relation (3) even looks very similar to the exact NSVZ β-function for the N = 1 supersymmetric electrodynamics [34,35], but it is necessary to remember that the anomalous dimension on the right hand side is calculated in the non-Abelian theory.
Eq. (4) has been derived by various methods, including instanton calculus [23,25,36], analysis of the anomaly supermultiplet [24,26,37], and the non-renormalization theorem for the topological term [38]. Explicit calculations in the ordinary perturbation theory, i.e. involving direct evaluation of contributions of multiloop Feynman diagrams, performed with dimensional reduction [39,40] accompanied by the DR subtraction scheme, yielded the β-function in the three- [41,42,43,44] and even four-loop approximations [45]. The result coincided with Eq. (4) only in the two-loop approximation. In higher loops Eq. (4) is satisfied only after a specially tuned finite renormalization of the coupling constant, which should be done in each order of perturbation theory [42,43,46]. Note that the very fact of its existence is quite non-trivial, as was noticed in [42], because the NSVZ relation imposes certain scheme-independent restrictions on the divergences [47,48]. The general equations describing the scheme dependence of the NSVZ relation have been derived in [49,48].
At present, in the case of using dimensional reduction there is no general prescription leading to the subtraction scheme in which the NSVZ relation would hold to all orders (the NSVZ subtraction scheme). Moreover, dimensional reduction is mathematically inconsistent [40], and an attempt to remove the inconsistencies causes supersymmetry breaking by quantum corrections in higher loops [50,51]. However, at present, it is not clear in which order this occurs. By contrast, the regularization by higher covariant derivatives appeared to be a remarkable computational tool for perturbative calculations in supersymmetric gauge theories, see e.g. [31,52,53,54,55]. It was first introduced in Refs. [18,19], and afterwards generalized to supersymmetric theories in [20,21]. It is mathematically consistent and preserves both gauge symmetry and supersymmetry. The use of the regularization by higher derivatives for multiloop calculations allowed constructing the NSVZ subtraction scheme in N = 1 SQED [56]. To formulate the corresponding renormalization prescription, it is necessary to fix a value x 0 of x = ln(Λ/µ), where Λ is the dimensionful parameter of the regularized theory and µ is a renormalization scale, and require that the renormalization constants satisfy the condition Then one can show that in this subtraction scheme the β-function and the anomalous dimension defined in terms of the renormalized coupling constant coincide with the corresponding RG functions defined in terms of the bare coupling constant in all orders of perturbation theory [56]. Given that for the latter RG functions the NSVZ relation holds, its validity directly follows for the former ones. The non-Abelian analogue of the boundary conditions (5) has been written in [57]. It was verified by an explicit three-loop calculation for the terms quartic in the Yukawa couplings in [55]. This calculation is non-trivial, because the considered part of the three-loop β-function is scheme-dependent. The result exactly confirms the prescription for the NSVZ scheme given in [57]. The subtraction scheme in which the photino mass renormalization in softly broken N = 1 SQED satisfies the NSVZ-like relation [58,59,60] can also be constructed in a similar way [61]. However, a subtraction scheme in which the exact expression for the Adler D-function is valid for the RG functions defined in terms of the renormalized coupling constant has not yet been found. Certainly, it would be interesting to construct such a scheme and to demonstrate how it works by an explicit calculation. Note that a non-trivial check can be obtained only starting from the three-loop approximation. The three-loop calculation is also a useful verification of Eq. (3), because its proof in Ref. [17], although exhaustive, is rather technical and complicated. This paper is devoted to the calculation of the three-loop Adler D-function and the two-loop anomalous dimension for N = 1 SQCD regularized by higher covariant derivatives and to the verification of the relation (3) in the corresponding approximation. By an explicit calculation we prove the validity of Eq. (3) for the RG functions defined in terms of the bare coupling constant at three loops. Next, we construct the subtraction scheme in which the Adler D-function and the anomalous dimension defined in terms of the renormalized coupling constant satisfy Eq. (3) to all orders and illustrate it with a three-loop calculation. Our result differs from the one for N = 1 SQED obtained in [56] only in that we have to impose the condition similar to (5) on the renormalization constant for the strong coupling constant as well.
Let us comment on the version of the higher covariant derivative regularization we use in this work. Introduction of higher covariant derivatives into the action significantly modifies the vertices, which makes it hard to use this kind of regularization in practical calculations. But this circumstance is mitigated by the structure of quantum corrections in supersymmetric gauge theories regularized by higher derivatives, which we mentioned earlier. Moreover, calculations can be simplified even more, if one uses a version of regularization by higher derivatives breaking the BRST symmetry (a two-loop calculation with such a regularization was done in [52]). As a result, the Slavnov-Taylor identities have to be restored by the introduction of additional counterterms into the action (the corresponding renormalization procedure was constructed for both non-supersymmetric and supersymmetric Yang-Mills theories in [62,63] and [64,65], respectively). On the other hand, it can sometimes be useful to sacrifice the convenience of calculations in favour of the BRST invariance. For example, a recent computation in [54] with a BRST-preserving version of the higher derivative regularization in the Yang-Mills theory revealed the non-renormalization of the triple vertices with two ghost legs and one leg of the quantum gauge superfield at one loop, which later transformed into a non-renormalization theorem proved in all loops in [57], with the proof to a great extent relying on the Slavnov-Taylor identities. In this work we regularize the theory by higher covariant derivatives without breaking the BRST invariance, which considerably complicates the vertices but, nevertheless, leads to calculable expressions.
The paper is organized as follows: in Sect. 2 we describe N = 1 SQCD (with the matter superfields in an arbitrary representation of the gauge group) and specify the version of the higher covariant derivative regularization which we will use throughout the calculations. In Sect. 3 we provide the three-loop expression for the Adler D-function written in terms of loop integrals. In Sect. 4 the result for the anomalous dimension in the two-loop approximation is presented as a sum of integrals and the validity of Eq. (3) is verified in the considered approximation. In Sect. 5 we construct the renormalization prescription fixing the subtraction scheme in which the relation (3) is valid to all loops for the RG functions defined in terms of the renormalized coupling constant. This prescription is verified by an explicit three-loop calculation in Sect. 6 after obtaining the final expressions for γ(α s0 ) and D(α s0 ) by evaluating the integrals for the simplest choice of the higher derivative regulator function. Feynman amplitudes for two-and three-loop diagrams contributing to the Adler D-function can be found in Appendix A, for oneand two-loop diagrams contributing to the anomalous dimension, in Appendix B. Appendix C is devoted to the evaluation of the integrals arising in the three-loop calculation.
2 N = 1 SQCD: the action, regularization, and renormalization We will consider the N = 1 supersymmetric non-Abelian gauge theory with an arbitrary simple gauge group G interacting with an external Abelian gauge field. The chiral matter superfields φ α and φ α for each value of α = 1, . . . , N f are assumed to lie in the representation R +R. The bare action of the theory in the massless limit can be written as This theory is invariant under the gauge transformations of the G × U (1) group. V is the non-Abelian gauge superfield and g 0 is the corresponding bare coupling constant.
where t A are generators of the fundamental representation of the group G normalized by the condition tr(t A t B ) = δ AB /2. In the terms containing the matter superfields V = g 0 V A T A , where T A are generators of the representation R. The Abelian external classical gauge superfield is denoted by V and e 0 is the bare coupling constant corresponding to the U (1) group. Each of the chiral matter superfields φ α belongs to the representation R of the group G and has the charge q α with respect to U (1). The superfields φ α belong to the representationR and have the U (1) charge −q α . The regularization is implemented by adding to the action the term containing higher covariant derivatives, 3 where Λ is the dimensionful parameter of the regularized theory playing the role of the ultraviolet cut-off. The function R(y), such that R(0) = 1, is the regulator which should rapidly grow at infinity, and we introduced the covariant derivatives in the chiral representation Note that regularizing terms of this type preserve BRST invariance [66,67,68] of the total action after the gauge-fixing procedure. We perform gauge fixing by adding the term where ξ 0 is the bare gauge parameter. Due to the Slavnov-Taylor identities quantum corrections to the two-point Green function of the gauge superfield V are transversal, so that the gauge fixing term is not renormalized. This implies that the renormalization of the gauge parameter ξ is related to the renormalization of the coupling constant g and of the superfield V . According to [54], the one-loop running of the gauge parameter is described by the equation where ξ and g are the renormalized gauge parameter and the renormalized coupling constant, respectively. Choosing, for simplicity, the Feynman gauge ξ = 1, in the approximation considered in this paper we can make in the gauge fixing term the replacement ξ 0 g 2 0 → g 2 . Then the propagator of the gauge superfields V A is given by the expression The action for the Faddeev-Popov ghosts has the form with the subscript Adj implying that the ghosts lie in the adjoint representation of the non-Abelian group. After adding the term S Λ to the action S only the one-loop divergences remain [69]. They should be regularized by the help of the Pauli-Villars method [70]. The one-loop divergences coming from the matter loop can be regularized by introducing a single set of anticommuting Pauli-Villars superfields Φ α and Φ α with the action where M = aΛ with a being a dimensionless parameter independent of the coupling constant. Following Ref. [54], to cancel the divergences coming from the loops of the gauge superfield V and of the Faddeev-Popov ghosts in the one-loop approximation, we introduce an additional set of (commuting) Pauli-Villars superfields with the action where M ϕ = a ϕ Λ and a ϕ is a parameter analogous to a which does not depend on the coupling constant. Due to the gauge invariance of the action under the U (1) gauge transformations, quantum corrections to the two-point Green function of V are transversal, Here α 0 = e 2 0 /4π is the bare coupling constant corresponding to the U (1) group, α s0 = g 2 0 /4π is the bare coupling constant corresponding to the group G, and is the transversal projection operator. Similarly, quantum corrections to the two-point Green function of the non-Abelian quantum gauge superfield are also transversal due to the BRST symmetry. This fact can be proved in the standard way by using the Slavnov-Taylor identities [71,72]. Quantum contributions to the two-point Green function of the chiral matter superfields enter the effective action as The renormalized coupling constants α(α 0 , α s0 , Λ/µ), α s (α s0 , Λ/µ) and the renormalization constants for the chiral matter superfields Z i j (α s , Λ/µ), such that are defined in the standard way. 4 Instead of the functions α(α 0 , α s0 , Λ/µ) and α s (α s0 , Λ/µ) one can equivalently use the renormalization constants respectively.
In the one-loop approximation the coupling constant α s is related to the corresponding bare constant by the equation where b 11 and b 12 are constants whose particular values depend on the subtraction scheme. Following [16,17], the Adler D-function in terms of the bare coupling constant is defined as This function can be expressed in terms of the two-point Green function of the Abelian gauge superfield V , The equality follows from the fact that the derivative of the function d −1 expressed in terms of the renormalized coupling constants with respect to ln Λ, with α and α s fixed, vanishes in the limit p → 0. This limit is needed in order to get rid of the finite terms proportional to (p/Λ) k , where k is a positive integer.
Similarly, the anomalous dimension of each of the matter superfields φ α (which is the same for all α = 1, . . . , N f ) defined in terms of the bare coupling constant can be also expressed in terms of the two-point Green function of the chiral matter superfields, In this paper we argue that for the theory (6) in the three-loop approximation the D-function and the anomalous dimension of the matter superfields are related by the equation where dim(R) is the dimension of the representation R. This equation is the generalization of Eq. (3) to the considered case. Really, if the chiral matter superfields φ α belong to the fundamental representation of the SU (N ) group, then so that Substituting these expressions into Eq. (24), we obtain Eq. (3) derived in [16,17].

The three-loop Adler D-function
According to Eq. (22) the Adler function defined in terms of the bare coupling constant is related to the two-point Green function of the superfield V . That is why the D-function is contributed to by superdiagrams with two external lines of the Abelian gauge superfield (which certainly do not contain its internal lines). There is a large number of such diagrams in the three-loop approximation. That is why we do not draw all of them in this paper. Instead of this, in Fig. 1 we present only graphs without external lines. Attaching two external lines of the Abelian gauge superfield V in all possible ways to the loops of the matter superfields in these diagrams, we obtain the total three-loop contribution to the Adler D-function.
In Fig. 1 the solid lines represent propagators of the matter superfields and the corresponding Pauli-Villars superfields, the wavy lines correspond to propagators of the non-Abelian gauge superfield V , and the dashed lines denote propagators of the Faddeev-Popov ghosts. The dotted lines represent propagators of the Pauli-Villars superfields ϕ i . They are considered separately from the Pauli-Villars superfields Φ α and Φ α , because ϕ i do not interact with the external Abelian gauge superfield V unlike Φ α and Φ α .
To obtain the D-function in the considered approximation, it is necessary to attach two external lines of the superfield V to the solid lines in these graphs and calculate the expression (22) finding the function d −1 − α −1 0 by the help of Eq. (15). The results for contributions of the graphs presented in Fig. 1 to Eq. (22) are collected in Appendix A. All these expressions are given by integrals of double total derivatives. Their sum can be written as where I i (q) with i = 1, . . . , 4 are given by the differences of the massless contributions of the usual matter superfields φ, φ and the massive contributions of the Pauli-Villars superfields Φ, Φ,  Figure 1: Two-and three-loop graphs contributing to the Adler D-function. Note that, for simplicity, some graphs vanishing in the considered (Feynman) gauge are not presented in this figure.
The contribution I 5 (q) comes from the graphs (7), (8) and (9), which have two loops of the matter superfields, and has a different structure. The explicit expressions for the integrals I 1 (q, m), . . . , I 4 (q, m), I 5 (q), and the one-loop contribution I 0 (q) have the form .
The function f (k/Λ) entering the expression for I 4 (q, m) is related to the one-loop contribution to the polarization operator of the quantum gauge superfield and has been calculated in Ref. [73]. It is given by Eq. (76) in Appendix A. From Eq. (27) we see that in the considered approximation the Adler D-function is given by integrals of double total derivatives with respect to the momentum q µ of the matter loop. This exactly agrees with the results of Refs. [16,17]. Such a structure of loop integrals allows taking the integral over d 4 q. Note that although being integrals of total derivatives, the contributions to the D-function do not vanish due to singularities of the integrands. When evaluating each of the integrals, one should surround singular points by spheres of infinitely small radii and transform the volume integral into a sum of surface integrals over these spheres and over the infinitely large sphere. The last one can be discarded if the integrand decreases at infinity rapidly enough, where the function F (q 2 ) is assumed to be non-singular and rapidly decreasing at infinity. Using this equation it is possible to take all integrals over d 4 q in Eq. (27). Then, after some transformations, the result for the Adler D-function can be written in the form Note that in Eq. (36) the differentiation with respect to ln Λ should precede the integration. It is necessary in order to get rid of the infrared divergences and to obtain well-defined integrals. This differentiation should be done with the renormalized coupling constant α s fixed, as is required by Eq. (22).

The two-loop anomalous dimension of the matter superfields
According to Eq. (24) the three-loop Adler function is related to the two-loop anomalous dimension of the matter superfields. Therefore, to verify this equation it is necessary to calculate the two-loop contribution to the two-point Green function of the matter superfields. The corresponding superdiagrams are presented in Fig. 2. Note that they can be obtained from the graphs in Fig. 1 by cutting the matter loops in as many ways as possible. As we are calculating their contribution to the effective action, only one particle irreducible diagrams resulting from this procedure should be kept. The expressions for quantum corrections received by G(α s0 , Λ/q) from these diagrams are collected in Appendix B. Using them, the anomalous dimension defined in terms of the bare coupling constant α s0 is calculated by the help of Eq. (23). Namely, it is necessary to differentiate the sum of all contributions to ln G with respect to ln Λ and take the limit of the vanishing external momentum.  Figure 2: Diagrams contributing to the two-loop anomalous dimension. Each graph in this figure can be obtained from the corresponding graph in Fig. 1 by cutting a line of the matter superfields. For simplicity, some evidently vanishing diagrams are omitted exactly as in Fig. 1.
where the function f (k/Λ) is the same as in (36) and is given by Eq. (76). Comparing Eqs. (36) and (37) one can easily verify that the NSVZ-like relation (24) is really valid in the considered approximation for the RG functions defined in terms of the bare coupling constant α s0 , The equality takes place, because both sides of this equation are given by the same integrals, so it is valid independently of the subtraction scheme. This follows from the fact [56] that the RG functions defined in terms of the bare coupling constant are scheme-independent for a fixed regularization. This can be confirmed by an explicit calculation in the considered approximation, for the simplest form of the regulator R(y) (see Eq. (7)). Let us choose it in the form where n is a positive integer. Then it is possible to find explicit expressions for the RG functions entering Eq. (38) calculating the integrals in Eq. (37). Details of this calculation are described in Appendix C, and the result has the form where Therefore, from Eq. (38) we conclude that for the regulator (39) the three-loop Adler D-function defined in terms of the bare coupling constant is given by the expression We see that finite constants defining the subtraction scheme (such as b 11 and b 12 in Eq. (20)) do not enter Eqs. (40) and (42). Consequently, these functions are scheme independent (for a fixed regularization) in agreement with the general statement proved in [56].

The NSVZ-like subtraction scheme
The RG functions entering the relation (24) are defined in terms of the bare coupling constant and, therefore, do not depend on the subtraction scheme for a fixed regularization [56]. On the other hand, RG functions are usually defined as functions of renormalized couplings [74] and depend on the renormalization prescription. This implies that the relation (24) is not valid for the Adler D-function and the anomalous dimension defined in terms of the renormalized coupling constant in an arbitrary subtraction scheme. We will show that the scheme in which (24) holds for the RG functions defined in terms of the renormalized coupling constant (the NSVZ-like scheme) can be constructed by imposing simple conditions on the renormalization constants. Namely, let us fix a value x 0 of ln(Λ/µ) and impose the requirement that the renormalization constants satisfy the equations In the particular case x 0 = 0 this is an analogue of minimal subtractions (see, e.g., [75]) in the sense that only terms proportional to ( ln Λ/µ) n , with n ≥ 1 being a positive integer, survive in the renormalization constants. We will call this scheme HDMS (i.e. Higher Derivatives + Minimal Subtractions). First, we prove that such a scheme exists. For convenience, let us introduce a new variable x = ln(Λ/µ). Then, we consider some renormalization prescription and find the functions α 0 (α, α s , x), α s0 (α s , x) and Z(α s , x) i j for it. In general, the conditions (43) are not satisfied, so that Evidently, for a fixed finite value of x 0 the functions a(α, α s ), b(α s ) and g(α s ) i j are finite. Then we choose new renormalized coupling constants α ′ , α ′ s and new renormalization constants for the chiral matter superfields Z ′ i j such that By construction, the renormalized Green functions expressed in terms of the renormalized coupling constants α and α s are finite. Certainly, they remain finite after the finite renormalization (45). Furthermore, from Eqs. (44) and (45) we conclude that the renormalization constants (43). This implies that the scheme fixed by the conditions (43) exists. Therefore, it only remains to show that in this scheme (24) is valid for the RG functions defined in terms of the renormalized coupling constant.
The Adler D-function and the anomalous dimension are defined in terms of the renormalized coupling constant by the equations From (47) we obtain Let us set x = x 0 in this equation. Then, the second term in this expression vanishes, because according to the first equation in (43) Moreover, the third equation in (43) gives α s (α s0 , x 0 ) = α s0 . Therefore, from Eq. (48) we obtain This implies that both definitions of the anomalous dimension give the same function if the boundary conditions (43) are satisfied. Similar statement can be proved for the D-function. For this purpose we rewrite the definition (46) in the equivalent form As earlier, the differentiation gives Setting x = x 0 , we see that the last two terms vanish due to the boundary conditions (43), while α s (α s0 , x 0 ) = α s0 and α(α 0 , α s0 , x 0 ) = α 0 . Consequently, we obtain Therefore, under the boundary conditions (43) the RG functions defined in terms of the renormalized coupling constants coincide with the ones defined in terms of the bare coupling constants. Taking into account that the latter functions satisfy the NSVZ-like relation in the case of using the higher covariant derivative regularization, we conclude that the former ones also satisfy it in the subtraction scheme specified by the prescription (43). This implies that Eq. (43) produces the NSVZ-like scheme in all orders of perturbation theory with the higher covariant derivative regularization.

The NSVZ-like scheme in the three-loop approximation
Let us verify the general statements considered in the previous section by the explicit threeloop calculation. We will start with Eq. (40) and integrate the RG equation (23). Solving this equation for ln Z and taking into account Eq. (20), which describes the one-loop running of the coupling constant α s , we obtain where g 1 , g 21 , g 22 , and g 23 are arbitrary finite constants. 5 Existence of these arbitrary constants follows from arbitrariness of choosing a subtraction scheme. Fixing values of these constants and the other similar constants, one fixes the subtraction scheme. The finite constants do not enter the RG functions defined in terms of the bare coupling constant, but they are present in the scheme dependent expressions for the RG functions defined in terms of the renormalized coupling constant. To calculate the anomalous dimension defined as a function of the renormalized coupling, one needs to rewrite the equation (54) in terms of α s0 , differentiate it with respect to ln µ, and express the result in terms of α s . This gives The presence of the finite constants g 1 , b 11 and b 12 confirms the scheme-dependence of this RG function.
Eq. (55) can be compared with the result obtained with dimensional reduction in the DRscheme, see [42] and references therein. In our notation it is written as We see that the scheme-independent terms proportional to C(R) 2 coincide, while the other terms coincide if we choose the values of finite constants satisfying the equations This implies that our result agrees with the result of [42], because the difference in the regularization can always be compensated by a proper choice of the subtraction scheme.
Integrating the RG equation (21) (keeping α s fixed) one can relate the renormalized coupling constant α (corresponding to the U (1) group) to the bare coupling constant α 0 , As we have already mentioned, the arbitrary constants d 1 , d 2 , d 31 , d 32 , and d 33 , which appeared here as constants of integration, together with g 1 , g 21 , g 22 , g 23 , b 11 , and b 12 , fix the subtraction scheme in the considered approximation. From Eq. (58) we can obtain the renormalization constant Z α = α/α 0 . Expressing the right hand side of (58) in terms of α s0 , differentiating both sides with respect to ln µ, and then returning back to α s we obtain the Adler D-function defined in terms of the renormalized coupling constant, Due to the presence of the finite constants this expression is scheme-dependent. 6 However, from Eqs. (55) and (59) we see that the terms proportional to C(R) 2 i j in the two-loop anomalous dimension and the terms proportional to tr C(R) 2 in the D-function are scheme independent. The scheme-independence of such terms in various RG functions was discussed in Refs. [76,77,78,47] in the Abelian case and in Ref. [48] in the non-Abelian case. Now, let us fix the values of the finite constants by imposing the boundary conditions (43). First, we fix a value x 0 and consider the equation Z αs (α s , x 0 ) = 1. The expression for Z αs = α s /α s0 can be found from Eq. (20) and is written as  (54) we obtain the equation g 1 = −x 0 . Therefore, Consequently, the Adler D-function defined in terms of the renormalized coupling constant in the NSVZ-like scheme obtained with the higher covariant derivative regularization is For the particular case of the SU (N ) gauge group and the matter superfields in the fundamental (φ) and antifundamental ( φ) representations we obtain so that in the NSVZ-like scheme the anomalous dimension γ(α s ) i j ≡ γ(α s )δ j i and the Dfunction take the form Comparing Eqs. (40) with (55) and (42) with (59) we conclude that Eq. (61) ensures the identical equality between γ and γ as well as between D and D in the considered approximation.
This exactly confirms the general arguments presented in the previous section. Moreover, it is evident now that in the subtraction scheme (43) the RG functions defined in terms of the renormalized coupling constant satisfy the NSVZ-like relation in the considered approximation.

Conclusion
In the three-loop approximation for N = 1 SQCD regularized by higher covariant derivatives we have verified the validity of the NSVZ-like equation (3) (and its more general version (24)), which relates the Adler D-function to the anomalous dimension of the matter superfields defined in terms of the bare coupling constant. Moreover, we have obtained an explicit expression for the Adler function in the NSVZ-like scheme to the order O(α 2 s ). In full accordance with Ref. [17] the three-loop D-function defined in terms of the bare coupling constant is given by integrals of double total derivatives. After integration with respect to the momentum of the matter loop they coincide with the integrals giving the two-loop anomalous dimension. Note that doing the calculations we use a version of the higher covariant derivative regularization that preserves the BRST invariance of the action. Although it leads to complicated integrals, for the higher derivative regulator (39) the RG functions in the considered approximation can be explicitly calculated. Note that the RG functions defined in terms of the bare coupling constant satisfy the NSVZ-like relation independently of the renormalization prescription with the higher covariant derivative regularization. However, for the RG functions defined in terms of the renormalized coupling constant this relation is valid only in a certain subtraction scheme. In this paper we present the prescription which gives this scheme in all loops. It has also been verified by the explicit three-loop calculation, which confirmed its correctness in the considered approximation.
It would be also interesting to calculate O(α 2 s ) contribution to the Adler D-function in the DR-scheme. Having in mind the results of [42,43], one may suggest that the NSVZ-like equation (66) is not satisfied in the DR-scheme, and a specially tuned finite renormalization is needed to obtain the RG functions for which it is valid.
Each of them corresponds to the sum of all superdiagrams obtained from a considered supergraph by attaching two external lines of the Abelian gauge superfield. All these expressions are written in the Euclidean space after the Wick rotation. Note that, for convenience, some contributions are summed.
; (74) d d ln Λ graph(10) + graph(11) + . . . + graph (21) where the function f (k/Λ), related to the one-loop polarization operator, has the form B Supergraphs contributing to the two-point Green function of the matter superfields For the sake of completeness, in this section we present expressions for contributions of superdiagrams presented in Fig. 2 to the function G i j (α s0 , Λ/q). They are written in the Euclidean space after the Wick rotation. For convenience, certain groups of diagrams are summed.

C Calculation of the integrals
Let us calculate the anomalous dimension defined in terms of the bare coupling constant given by Eq. (37) for the regulator function (39). It is important that the differentiation with respect to ln Λ should be done at a fixed value of the renormalized coupling constant g (or, equivalently, α s = g 2 /4π) before the integrations. That is why, first, we rewrite the expression for the anomalous dimension in terms of α s using Eq. (20). In the considered approximation this gives Some integrals in this equation can be calculated using the results of Refs. [79,80] for the simplest regulator R(y) = 1 + y n , where n is a positive integer, The last equality follows from the fact that the expression in the square brackets is not singular in the limit l → 0. (Note that Eq. (99) is valid for an arbitrary function R(x) rapidly growing at infinity, such that R(0) = 1.) Therefore, we obtain Collecting the results for all integrals, the anomalous dimension defined in terms of the bare coupling constant can be presented as Rewriting the right hand side in terms of the bare coupling constant α s0 we finally obtain π 2 N f T (R)C(R) i j ( ln a + 1)