Fourloop renormalization of QCD with a reducible fermion representation of the gauge group: anomalous dimensions and renormalization constants
HJE
Fourloop renormalization of QCD with a reducible
K.G. Chetyrkin 0 1 2 4
M.F. Zoller 0 1 3
0 Luruper Chaussee 149 , 22761 Hamburg , Germany
1 WolfgangGaedeStra e 1 , 76131 Karlsruhe , Germany
2 II Institut fur Theoretische Physik, Universitat Hamburg
3 Institut fur Physik, Universitat Zurich , UZH
4 Institut fur Theoretische Teilchenphysik, Karlsruhe Institute of Technology , KIT
We present analytical results at fourloop level for the renormalization constants and anomalous dimensions of an extended QCD model with one coupling constant and an arbitrary number of fermion representations. One example of such a model is the QCD plus gluinos sector of a supersymmetric theory where the gluinos are Majorana fermions in the adjoint representation of the gauge group.
Perturbative QCD; Renormalization Group

elds
2
2
1 Introduction
2
Introduction
2.2.1
Direct fourloop calculation in the Feynman gauge with massive
HJEP06(217)4
tadpoles
Indirect fourloop calculation using threeloop massless propagators
computation of the gauge group factors
The behaviour of Green's functions w.r.t. a shift of the renormalization scale is described
by the anomalous dimensions of the elds and parameters of the theory, which enter the
Renormalization Group Equations (RGE). For QCD the full set of fourloop
renormalization constants and anomalous dimensions was presented in [2]. The results for the fourloop
QCD
function [3, 4] and the fourloop quark mass and eld anomalous dimensions had
already been available [5{7].1
In this paper we consider a model with a nonabelian gauge group, one coupling
constant and a reducible fermion representation, i.e. any number of irreducible fermion
representations. The function for the coupling this model was computed in an earlier work [1].
Here we provide the remaining Renormalization Group (RG) functions in full dependence
on the gauge parameter .
Apart from completing the set of renormalization constants and the RGE of the theory,
which is important in itself, the gauge boson and ghost propagator anomalous dimensions
serve another purpose. These quantities are essential ingredients in comparing the
momentum dependence of the corresponding propagators derived in nonperturbative calculations
on the lattice, with perturbative results (see e.g. [11{18]).
This paper is structured as follows:
rst, we will give the notation and de nitions
for the model and the computed RG functions We will also repeat how the special case of
QCD plus Majorana gluinos in the adjoint representation of the gauge group can be derived
from our more general results. Then we will present analytical results for the fourloop
1Recently, the veloop QCD
function has been obtained for QCD colour factors [8] as well as for a
generic gauge group [
9
] (see, also, [10]).
{ 1 {
anomalous dimensions of the gauge boson, ghost and fermion eld as well as the ones for
the ghostgluon vertex, the fermiongluon vertex and the fermion mass in Feynman gauge
for compactness. The renormalization constants and anomalous dimensions for a generic
gauge parameter
can be found in machine readable form as supplementary material to
this article.
of the gauge group is given by
HJEP06(217)4
(2.1)
(2.2)
(2.3)
(2.4)
(2.5)
(2.6)
with the structure constants f abc. We have one quadratic Casimir operator CF;r for each
fermion representation, de ned through
and CA for the adjoint representation. The dimensions of the fermion representations are
given by dF;r and the dimension of the adjoint representation by NA. The traces of the
di erent representations are de ned as
TF;r
ab = Tr T a;rT b;r
= Tiaj;rTjbi;r:
At fourloop level we also encounter higher order invariants in the gauge group factors
which are expressed in terms of symmetric tensors
daR1a2:::an =
1
X
n! perm
Tr T a (1);RT a (2);R : : : T a (n);R ;
LQCD =
+
4
1 a G
G
a
mq;r q;r q;r + gs q;rA=aT a;r q;r ;
with the gluon eld strength tensor
G
a
The index r speci es the fermion representation and the index q the fermion avour, q;r
is the corresponding fermion
eld and mq;r the corresponding fermion mass. The number
of fermion avours in representation r is nf;r for any of the Nrep fermion representations.
The generators T a;r of each fermion representation r ful ll the de ning anticomuting
relation of the Lie Algebra corresponding to the gauge group:
hT a;r; T b;ri = if abcT c;r
Tiak;rTkaj;r = ij CF;r;
{ 2 {
sentation, R = A, where Tbac;A =
i f abc.
where R can be any fermion representation r, noted as R = fF; rg, or the adjoint
repre
An important special case of this model is the QCD plus gluinos sector of a
supersymmetric theory where the gluinos are Majorana fermions in the adjoint representation of the
gauge group. Here we have Nrep = 2 and
nf;1 = nf ;
TF;1 = TF ;
CF;1 = CF ;
ng~
2
;
nf;2 =
TF;2 = CA;
CF;2 = CA;
(2.7)
HJEP06(217)4
should also be compensated by a factor 12 .
the factor 12 in front of the number of gluinos ng~ being a result of the Majorana nature
of the gluinos (see e.g. [19]). This can be understood in the following way: it has been
shown in [20] that one can treat Majorana fermions by rst drawing all possible Feynman
diagrams and choosing an arbitrary orientation (fermion ow) for each fermion line. Then
Feynman rules are applied in the same way as for Dirac spinors, especially one can use the
same propagators p= im for the momentum p along the fermion ow and
p=i m for p against
the fermion ow. Closed fermion loops receive a factor ( 1). One then applies the same
symmetry factors as for scalar or vector particles, e.g. a factor 12 for a loop consisting of
two propagators of Majorana particles. For this work we generate our diagrams using one
Dirac eld
for all fermions, i.e. we produce both possible fermion
ows in loops unless
they lead to the same diagram. The latter case is exactly the one where the symmetry
factor 12 must be applied. The rst case means that the loop was doublecounted which
By adding counterterms to the Lagrangian (2.1) in order to remove all possible UV
divergences we arrive at the bare Lagrangian expressed through renormalized elds, masses
and the coupling constant:
LQCD;B =
a
1 Z(4g)gs2 f abcAb Ac 2
1
2
3
2 1
4 1
1
(2.8)
Z(q;r) i
2
!
mq;rZm(q;r)Z(q;r)
2
were we have already used the fact that Z = Z(2g).
3
Due to the SlavnovTaylor identities all vertex renormalization constants are connected
and can be expressed through the renormalization constant of the coupling constant and
{ 3 {
the renormalization constants of the elds appearing in the respective vertex:
Zgs = Z(3g) Z(2g)
1 3
Zgs =
1
qZ(4g) Z(2g)
3
Zgs = Z(ccg)
1
Zgs = Z(q;r)
1
3
2
;
1
3
3
;
3
2
Z(2c)qZ(2g)
(a) = a2 dda za(1)(a) :
"a + (a) ;
(2.9)
(2.10)
(2.11)
(2.12)
(2.13)
(2.14)
(2.15)
(2.16)
(2.17)
(2.18)
and the fact that
=
3
and the function of the model.
from the corresponding anomalous dimension, a nite and usually more compact quantity,
(2g). Using (2.17) one can reconstruct renormalization constants
{ 4 {
In the MSscheme using regularization in D = 4
2" space time dimensions all
renormalization constants have the form
where a = 16gs2 2 . From the fact that the bare parameter aB = Zaa 2" (with Za = Zg2s ) does
not depend on the renormalization scale
one nds
Given a renormalization constant Z the corresponding anomalous dimension is de ned as
(a; ) =
2 d log Z(a; )
d 2
= a
X
n=1
(n)( ) an :
From the de nition of anomalous dimensions (2.16) it follows that
(a; ) = ("a
(a))
d log Z(a; )
da
3
The 1particleirreducible Feynman diagrams needed for this project were generated with
QGRAF [21]. We compute Z
2
3
3
(2c), Z(2g) and Z(q;r) from the 1PI selfenergies of the elds
Aa , c and q;r as well as Z(ccg) and Z(q;r) from the respective vertex corrections and Zm(q;r)
1 1
from the 1PI corrections to a Green's function with an insertion of one operator
and an external fermion line of type (q; r). We used two di erent methods to calculate
these objects, rst a direct fourloop calculation in Feynman gauge with massive tadpoles
and then an indirect method where fourloop objects are constructed from propagatorlike
threeloop objects to derive the full dependence on the gauge parameter := 1
.
Direct fourloop calculation in the Feynman gauge with massive tadpoles
counterterm M2 ZM(22g) Aa Aa
2
= 0 (Feynman gauge) the topologies of the diagrams were identi ed with the C++
programs Q2E and EXP [22, 23]. In this approach all diagrams were expanded in the
external momenta in order to factor out the momentum dependence of the treelevel
vertex or propagator, e.g. q q
q2g
were projected onto scalar integrals, using e.g. q q
for the gluon selfenergy. Then the tensor integrals
q4
as well as gq2 as projectors for the
gluon selfenergy. After this we set all external momenta to zero since the UV
divergent part of the integral does not depend on
nite external momenta. We then use the
method of introducing the same auxiliary mass parameter M 2 in every propagator
denominator [24, 25]. Subdivergencies / M 2 are cancelled by an unphysical gluon mass
restoring the correct UV divergent part of the diagrams.
This method was e.g. used in [3, 4, 26{30] and is explained in detail in [31].
For the expansions, application of projectors, evaluation of fermion traces and
counterterm insertions in lower loop diagrams we used FORM [
32, 33
]. The scalar tadpole
integrals were computed with the FORMbased package MATAD [34] up to threeloop
order. At four loops we use the C++ version of FIRE 5 [35, 36] in order to reduce the scalar
integrals to Master Integrals which can be found in [4]. Technical details of the reduction
are described in the previous paper [
29
].
2.2.2
Indirect fourloop calculation using threeloop massless propagators
The case of a generic gauge parameter is certainly possible to treat in the same massive
way but calculations then require signi cantly more time and computer resources.2 As a
result we have chosen an alternative massless approach which reduces the evaluation of any
Lloop Zfactor to the calculation of some properly constructed set of (L
1)loop massless
propagators [38{41]. As is wellknown (starting already from L = 2 [42]) calculation of
Lloop massive vacuum diagrams is signi cantly more complicated and timeconsuming
than the one of corresponding (L
1)loop massless propagators.
The approach is easily applicable for any Zfactor except for Z3 [2]. The latter problem
is certainly doable within the massless approach but requires signi cantly more human
e orts in resolving rather sophisticated combinatorics.3
On the other hand, one could
2Nevertheless, it has been done recently along theses lines in [37] for the case of one irreducible fermion
representation.
3Very recently the problem has been successfully solved in two radically di erent ways [8] and [
9
].
{ 5 {
ψ
ψ
(a)
(d)
R2
R3
R3
R1
R2
R1
Aa
μ
R4
R4
¯
ψ
¯
ψ
(b)
ψ
(e)
R2
R1
Aa
μ
R1
R2
R3
R3
¯
ψ
Ab
ν
(c)
ψ
(f)
c
a
R1
R2
R1
¯
ψ
c¯
b
restore the full dependence of Z3 from all other renormalization constants and from the
fact that the charge renormalization constant Zg is gauge invariant [2, 37]. As Zg in QCD
with fermions transforming under arbitrary reducible representation of the gauge group
has been recently found in [1] we have proceeded in this way. For calculation of 3loop
massless propagator we have used the FORM version of MINCER [43].
2.2.3
computation of the gauge group factors
The calculation of the gauge group factors was done with an extended version of the FORM
package COLOR [44] already used and presented in [1]. We take the following steps:
1. For the generation of the diagrams in QGRAF [21] we use one eld A for the adjoint
representation (gauge boson) and one eld
for all the fermion representations.
This has the advantage that we do not produce more Feynman diagrams than in
QCD. Each fermion line in a diagram gets a line number and is treated as a di erent
representation from the other fermion lines. Since we compute diagrams up to
fourloop order we need up to four di erent line representations R1; : : : ; R4 (see gure 1)
and Tiaj;R4 = T4(i,j,a). Each fermion loop gets assigned a factor nf .
2. The modi ed version of COLOR [1, 44] then writes the generators into traces
Tr T a1;R : : : T an;R
= TRfRg(a1,...,an); (R = R1; : : : ; R4)
(2.19)
which are then reduced as outlined in [44] yielding colour factors expressed through
traces TFfRg, the Casimir operators cFfRg and cA, the dimensions of the
representations dFfRg and NA.
{ 6 {
3. Now we change from fermion line numbers R1; : : : ; R4 to four explicit physical fermion
representations r by substituting each of the line numbers R1; : : : ; R4 by the sum over
all representations r = 1; : : : ; 4. An example of the substitution of fR1; : : : ;
R4gcolour factors with those of the physical representaions in a oneloop diagram is
Nf*TF1 ! nf;1TF;1 + nf;2TF;2 + nf;3TF;3 + nf;4TF;4:
(2.20)
At higher orders this subtitution becomes much more involved.4 Diagram (a) from
gure 1 now corresponds to a sum of 44 = 256 diagrams with explicit fermion
representations. This lengthy representation of our results is needed for the
renormalization procedure, since e.g. a one loop counterterm to the gluon selfenergy, computed
from a diagram with only R1, must be applied to all the fermion loops in
gure 1
(a,b,d,e). This is most conveniently achieved if each fermionloop is considered a
sum over all physical fermion representations just as it is considered a some over all
(massless) fermion
avours.5 The factors involving daF1;ra2a3a4 , daF1;ra2a3 , da1a2a3a4 and
A
da1a2a3 appear only at fourloop level and do hence not interfere with lower order
A
diagrams with counterterm insertions. They can be treated directly in the next step.
4. After all subdivergencies are cancelled by adding the lowerloop diagrams with
counterterm insertions we simplify and generalize the notation. The explicit colour factors
are collected in sums of terms built from nf;r, CF;r and TF;r over all physical
representations r, e.g.6
nf;1TF;1 !
X nf;iTF;i
nf;2TF;2
nf;3TF;3
nf;4TF;4:
(2.21)
Since we used the maximum number of di erent fermion representations which can
appear in any diagram the result is valid for any number of fermion
representations Nrep.
3
Results
In this section we give the results for the anomalous dimensions of the QCDlike model
with an arbitrary number of fermion representations as described above to fourloop level.
The number of active fermion avours of representation i is denoted by nf;i. Apart from
the Casimir operators CA and CF;i and the trace TF;i the following invariants appear in our
4For this reason it is convenient to collect all combinations Nfx1*TF1x2*CF1x3*TF2x4*CF2x5*TF3x6*CF3x7
*TF4x8*CF4x9 in a function C(x1,...,x9). The factors C(x1,...,x7) are then substituted by the proper
combinations of nf;1, TF;1, cF;1, etc.
5Since renormalization constants in the MSscheme do not depend on masses all fermion avours can be
treated as massless for their computation.
6For
convenience
we
collect
nfx;11nf;2
x2 nf;3
x3 nfx;44TFy;11TFy;22TFy;33TFy;44CFz1;1CFz2;2CFz3;3CFz4;4 in a function
CR(x1,...,x4,y1,...,y4,z1...,44) which are then substituted by the proper sums of colour
factors over all representations r.
{ 7 {
dabcddabcd
A
NA
dF;r
d(A4A) =
A ; d(F4A);i =
d~(4)
FA;r = daFb;rcddabcd
A ; d~(F4F);ri = daFb;rcddaFb;icd ;
NA
dF;r
daFb;icddabcd
A ; d(F4F);ij = daFb;icddaFb;jcd ;
NA
where r is xed and i; j will be summed over all fermion representations. In this section
we give the results for
= 1 (Feynman gauge), the general case
= (1
) can be found
as supplementary material to this article.
anomalous dimension according to (2.16)
From the gauge boson
eld strength renormalization constant Z3(2g) we compute the
243 CA2 + X nf;iTF;i (4CF;i + 5CA) ;
From the ghost eld strength renormalization constant Z3(2c) we compute
3
3
(2c) (1)
(2c) (2)
=
=
1
2 CA;
2449 CA2 +
5
+CA
X nf;inf;jTF;iTF;j CF;j
8315
972
+
86
4 3 + CA
X nf;iTF;i CF;i
4
+
From the fermion eld strength renormalization constant Z2(q;r) we nd
2
2
2
{ 9 {
(3.8)
(3.9)
(3.10)
(3.11)
(3.12)
CF;r
for the anomalous dimension of a representation r fermion eld.
The fermion eldgauge bosonvertex renormalization constant Z1(q;r) yields
(3.13)
(3.14)
(3.15)
(3.16)
3 + C2 6307
A
+
for each representation r and the ghostgauge bosonvertex renormalization constant Z1(ccg)
yields
(3.18)
(3.19)
(3.20)
(3.21)
(3.22)
(3.23)
(3.24)
Finally, the mass anomalous dimension computed from Zm(q;r) is found to be
m
m
m
3 :
(3.25)
We checked that the well known relations
(a)
(a)
a
a
= 2 1(ccg)(a; )
= 2 1(q;r)(a; )
2 3(2c)(a; )
2 2(q;r)(a; )
3
(2g)(a; );
3
are ful lled with the function from [1]. This is also true if we include the full dependence
on the gauge parameter
in the anomalous dimensions. This dependence cancels
in the function. We provide renormalization constants and anomalous dimensions with
the full gauge dependence as supplementary material to this article. We compared these
fully dependent results with [37] for one fermion representation and nd full agreement.
4
Conclusions
We have presented analytical results for the eld anomalous dimensions 3(2g), 3(2c), 2(q;r),
the vertex anomalous dimensions 1(ccg) and 1(q;r) and the mass anomalous dimension (q;r)
m
in a QCDlike model with arbitrarily many fermion representations and with the full
dependence on the gauge parameter .
Acknowledgments
The work by K.G. Chetykin was supported by the Deutsche Forschungsgemeinschaft
through CH1479/11 and in part by the German Federal Ministry for Education and
Research BMBF through Grant No. 05H2015. The work by M.F. Zoller was supported by
the Swiss National Science Foundation (SNF) under contract BSCGI0 157722.
Open Access.
Attribution License (CCBY 4.0), which permits any use, distribution and reproduction in
any medium, provided the original author(s) and source are credited.
[INSPIRE].
[hepph/9703278] [INSPIRE].
[1] M.F. Zoller, Fourloop QCD
function with di erent fermion representations of the gauge
group, JHEP 10 (2016) 118 [arXiv:1608.08982] [INSPIRE].
[2] K.G. Chetyrkin, Fourloop renormalization of QCD: Full set of renormalization constants
and anomalous dimensions, Nucl. Phys. B 710 (2005) 499 [hepph/0405193] [INSPIRE].
[3] T. van Ritbergen, J.A.M. Vermaseren and S.A. Larin, The Four loop function in quantum
chromodynamics, Phys. Lett. B 400 (1997) 379 [hepph/9701390] [INSPIRE].
[4] M. Czakon, The Fourloop QCD
function and anomalous dimensions, Nucl. Phys. B 710
(2005) 485 [hepph/0411261] [INSPIRE].
[5] J.A.M. Vermaseren, S.A. Larin and T. van Ritbergen, The four loop quark mass anomalous
dimension and the invariant quark mass, Phys. Lett. B 405 (1997) 327 [hepph/9703284]
[6] K.G. Chetyrkin, Quark mass anomalous dimension to O( S4), Phys. Lett. B 404 (1997) 161
[7] K.G. Chetyrkin and A. Retey, Renormalization and running of quark mass and eld in the
regularization invariant and MSbar schemes at three loops and four loops, Nucl. Phys. B
583 (2000) 3 [hepph/9910332] [INSPIRE].
[8] P.A. Baikov, K.G. Chetyrkin and J.H. Kuhn, FiveLoop Running of the QCD coupling
constant, Phys. Rev. Lett. 118 (2017) 082002 [arXiv:1606.08659] [INSPIRE].
YangMills theory with fermions, JHEP 02 (2017) 090 [arXiv:1701.01404] [INSPIRE].
general gauge group, JHEP 07 (2016) 127 [arXiv:1606.08662] [INSPIRE].
[11] H. Suman and K. Schilling, First lattice study of ghost propagators in SU(2) and SU(3) gauge
HJEP06(217)4
theories, Phys. Lett. B 373 (1996) 314 [heplat/9512003] [INSPIRE].
[12] D. Becirevic et al., Asymptotic scaling of the gluon propagator on the lattice, Phys. Rev. D
61 (2000) 114508 [hepph/9910204] [INSPIRE].
Rev. D 60 (1999) 094509 [hepph/9903364] [INSPIRE].
[13] D. Becirevic et al., Asymptotic behavior of the gluon propagator from lattice QCD, Phys.
[14] D. Becirevic et al., Gluon propagator, triple gluon vertex and the QCD coupling constant,
Nucl. Phys. Proc. Suppl. 83 (2000) 159 [heplat/9908056] [INSPIRE].
[15] L. von Smekal, K. Maltman and A. Sternbeck, The strong coupling and its running to four
loops in a minimal MOM scheme, Phys. Lett. B 681 (2009) 336 [arXiv:0903.1696]
[INSPIRE].
[arXiv:1012.3135] [INSPIRE].
[16] ETM collaboration, B. Blossier et al., S from Lattice QCD: progresses and perspectives for
a realistic fullQCD determination of the running Strong coupling, PoS(ICHEP 2010)372
[17] B. Blossier et al., RI/MOM renormalization constants (Nf = 4) and the strong coupling
constant (Nf = 2 + 1 + 1) from twistedmass QCD, PoS(LATTICE 2011)223
[arXiv:1111.3023] [INSPIRE].
[18] V.G. Bornyakov, E.M. Ilgenfritz, C. Litwinski, V.K. Mitrjushkin and M. MullerPreussker,
Landau gauge ghost propagator and running coupling in SU(2) lattice gauge theory, Phys.
Rev. D 92 (2015) 074505 [arXiv:1302.5943] [INSPIRE].
[19] L. Clavelli, P.W. Coulter and L.R. Surguladze, Gluino contribution to the three loop
function in the minimal supersymmetric standard model, Phys. Rev. D 55 (1997) 4268
[hepph/9611355] [INSPIRE].
interactions, Nucl. Phys. B 387 (1992) 467 [INSPIRE].
[20] A. Denner, H. Eck, O. Hahn and J. Kublbeck, Feynman rules for fermion number violating
[INSPIRE].
[INSPIRE].
diagrams, hepph/9905298 [INSPIRE].
[21] P. Nogueira, Automatic Feynman graph generation, J. Comput. Phys. 105 (1993) 279
[22] T. Seidensticker, Automatic application of successive asymptotic expansions of Feynman
[23] R. Harlander, T. Seidensticker and M. Steinhauser, Complete corrections of O(
s) to the
decay of the Z boson into bottom quarks, Phys. Lett. B 426 (1998) 125 [hepph/9712228]
[24] M. Misiak and M. Munz, Two loop mixing of dimension ve avor changing operators, Phys.
Lett. B 344 (1995) 308 [hepph/9409454] [INSPIRE].
selfinteraction in the Standard Model, JHEP 06 (2012) 033 [arXiv:1205.2892] [INSPIRE].
fourloop level, JHEP 02 (2016) 095 [arXiv:1508.03624] [INSPIRE].
function of the Higgs selfcoupling in the SM and vacuum stability, JHEP 06 (2016) 175
[arXiv:1604.00853] [INSPIRE].
[arXiv:1407.6608] [INSPIRE].
[arXiv:0807.3243] [INSPIRE].
loops, JHEP 03 (2017) 020 [arXiv:1701.07068] [INSPIRE].
Feynman integrals in the limit of large momenta and masses, arXiv:1701.08627 [INSPIRE].
NIKHEFH9118 [INSPIRE]. [44] T. van Ritbergen, A.N. Schellekens and J.A.M. Vermaseren, Group theory factors for
[9] F. Herzog , B. Ruijl , T. Ueda , J.A.M. Vermaseren and A. Vogt , The veloop  function of [10] T. Luthe , A. Maier , P. Marquard and Y. Schr oder, Towards the veloop function for a [25] K.G. Chetyrkin, M. Misiak and M. Munz, functions and anomalous dimensions up to three loops , Nucl. Phys. B 518 ( 1998 ) 473 [ hep ph/9711266] [INSPIRE].
[26] Y. Schr oder, Automatic reduction of four loop bubbles , Nucl. Phys. Proc. Suppl . 116 ( 2003 ) 402 [ hep ph/0211288] [INSPIRE].
[29] M.F. Zoller , TopYukawa e ects on the function of the strong coupling in the SM at [30] K.G. Chetyrkin and M.F. Zoller , Leading QCDinduced fourloop contributions to the [31] M.F. Zoller , Threeloop function for the Higgs selfcoupling , PoS(LL2014) 014 [32] J.A.M. Vermaseren , New features of FORM, mathph/0010025 [INSPIRE].
[33] M. Tentyukov and J.A.M. Vermaseren , The Multithreaded version of FORM, Comput . Phys. [36] A.V. Smirnov , FIRE5: a C+ + implementation of Feynman Integral REduction, Comput . [40] K.G. Chetyrkin , Combinatorics of R, R 1 and R operations and asymptotic expansions of