Renormalization and ultraviolet sensitivity of gauge vertices in universal extra dimensions

When computing radiative corrections in models with compactified extra dimen- sions, one has to sum over the entire tower of Kaluza-Klein excitations inside the loops. The loop corrections generate a difference between the coupling strength of a zero-mode gauge boson and the coupling strength of its Kaluza-Klein excitation, although both originate from the same higher-dimensional gauge interaction. Furthermore, this dis- crepancy will in general depend on the cutoff scale and assumptions about the UV completion of the extra-dimensional theory. In this article, these effects are studied in detail within the context of the minimal universal extra dimension model (MUED). The broad features of the cutoff scale dependence can be captured through the so- lution of the functional flow equation in five-dimensional space. However, an explicit diagrammatic calculation reveals some modifications due to the compactification of the extra dimension. Nevertheless, when imposing a physical renormalization condition, one finds that the UV sensitivity of the effective Kaluza-Klein gauge-boson vertex is relatively small and not very important for most phenomenological purposes. Similar conclusions should hold in a larger class of extra-dimensional models besides MUED.


Introduction
In models with additional spatial dimensions, the Standard Model (SM) gauge dynamics are the low-energy remnant of a gauge theory in the full higher-dimensional space. Let us consider one flat extra dimension compactified on a circle with radius R and a non-Abelian gauge group G with gauge coupling g (5) in the five-dimensional (5D) theory. At tree level, there is a simple geometric relation between g (5) and the coupling g in the 4D compactified description, which is given by g (5) = g/ √ πR. At loop level, however, the 4D and 5D couplings develop a non-trivial dependence on the renormalization and cutoff scales [1][2][3]. This is intimately tied to the fact that extradimensional theories are superficially non-renormalizable. Moreover, in the compactified 4D picture, the couplings of the zero-mode gauge bosons and the Kaluza-Klein (KK) gauge bosons require different counterterms at one-loop order [4,5], although they both stem from the same 5D gauge coupling. This implies that these two 4D couplings, while identical at tree-level, differ by loop-induced correction terms. In fact, one also has to entertain the possibility that this correction may depend on assumptions about the nature of the UV completion of the 5D theory.
While the question of scale sensitivity of couplings in extra-dimensional theories has a long history [1][2][3][6][7][8][9], a detailed analysis of the difference between the effective zero-mode gauge bosons and their KK excitations is so far missing. Such an analysis is the topic of this paper, together with an extensive discussion of the influence of the UV cutoff scale on both vertices.
As a concrete example, the decay G 1 → Q 1q in the minimal universal extra dimension model (MUED) is considered, where G 1 and Q 1 are a level-1 KK gluon and quark, respectively, while q is a SM quark of the same flavor. This decay is kinematically allowed due to the loop-induced mass corrections to the KK-gluon and -quark [10,11]. The complete set of one-loop QCD corrections in MUED can be grouped into two parts: (a) diagrams involving only zero and level-1 modes in the loops, and (b) diagrams involving KK modes beyond level 1. Category (a) contains all IR divergent contributions, which cancel against corresponding real emission corrections. It is equivalent to the corrections within a lowenergy effective theory dubbed the two-site coloron model [5]. Category (b) contains the dependence on the cutoff scale and the renormalization procedure. In order to provide a physical renormalization condition for the G 1 Q 1q vertex, it must be compared to the gqq vertex in MUED, which will be matched to the SM below the compactification scale. Here g is the SM gluon field.
While the 5D gauge coupling does not run within the usual framework of dimensional regularization, its dependence on a cutoff scale can be described by the functional flow equation [12,13]. It is very illustrative to compare the flow equation analysis to the explicit summation over the KK tower, as will be shown below.
After a brief review of the conventions and notation for the MUED and coloron model in section 2, the calculation of the next-to-leading order (NLO) corrections to the decay G 1 → Q 1q in the coloron model is described in section 3. This set of corrections corresponds to the category (a) above. The contribution from level-n KK modes (n > 1) in the loops is discussed in section 4. Special emphasis is given to the dependence on the cutoff scale in different approaches (explicit Feynman diagram calculation, large-n approximation, functional flow equation). Finally, section 5 illustrates the impact of the cutoff dependence on physical quantities, such as the decay G 1 → Q 1q , before the conclusions are presented in section 6.

QCD sector of MUED and Coloron Model
In this section, a brief overview of the QCD sector of MUED and the two-site coloron model are given, to clarify the notation and conventions.
The Lagrangian for the strongly interacting MUED fields can be written as Here L gauge contains the kinetic term of the 5D gluon field. For practical purposes, one also needs to introduce a gauge-fixing term, L gf , and corresponding ghost Lagrangian, L ghost . In Feynman gauge these contributions read Here x 5 is the coordinate of the extra dimension, which is compactified on a circle with radius R. Furthermore, G a M N is the field strength tensor of the 5D gluon field V a M , while c a denotes the ghost field. Capital Latin indices run over all five dimensions (M, N = 0, 1, 2, 3, 5), whereas Greek indices (µ, ν, ...) run over the usual uncompactified four dimensions. The superscripts a, b, ... refer to SU(3) indices in the adjoint representation, whereas g (5) s refers to the 5D QCD coupling, which is related to the 4D coupling according to g πR. The coupling of the gluon field to quarks is described by where T a are the SU(3) generators in the fundamental representation, and Ψ/Ψ refer to the quark fields that are doublets/singlets under SU(2), respectively. Upon compactification and orbifolding, the 5D fields decompose into towers of 4D KK excitations: c a (x, c a n (x) cos Ψ (x, with P ± = (1 ± γ 5 )/2. Notice that the Dirac fermions Q n = Q n,L + Q n,R are constructed from components with the same SU(2)×U(1) quantum numbers, while the left and right components of q = q L + q R come from the 5D doublet and singlet respectively. * The explicit Lagrangian terms in the 4D compactified form can be found for example in Ref. [11].
The two-site coloron model described in Ref. [5] has the same degrees of freedom as zero and level-1 modes of the strongly interacting sector of MUED. It is based on a SU(3)×SU(3) symmetry broken by a non-linear sigma model, which ensures that gauge symmetry is conserved. This is in contrast to a naive truncation of MUED at the KK level 1, which violates the gauge symmetry of the 5D gluon field [14]. A detailed description of the model can be found in Ref. [5].
For the purposes of this work it suffices to point out that the only difference in the Feynman rules of the coloron model and MUED occurs for quadruple couplings of four level-1 bosons. Such couplings do not contribute to the calculation of level-1 KK gluon decays at LO and NLO. Therefore, the result for these decays in the coloron model is identical to what one would obtain in the naively truncated MUED † .

KK-gluon decay in Coloron Model
To account for the IR physics of the corrections to the G 1 Q 1q vertex we calculate the one-loop decay width Γ[G 1 →q Q 1 ] within the coloron model. As already mentioned above, a naïvely truncated MUED model violates gauge invariance [14], whereas the coloron model can be used as a gauge-invariant low-energy approximation of MUED.  Fig. 1 represents the only tree-level Feynman describing the coloron decay. The squared matrix element at leading order is given by

Diagram (A) in
where M G 1 and M Q 1 are the masses of the coloron and the heavy quark, respectively. When matching the coloron model to MUED at leading order, one obtains the constraint M G 1 = M Q 1 = R −1 , so that (11) vanishes. We therefore assume the radiatively induced mass splitting as in MUED, amounting to [11] where Λ is the high-energy cutoff scale of MUED.
Since the mass splitting is small in comparison to the mass scale of G 1 and Q 1 , it is sufficient to perform the calculation of the NLO corrected matrix element only to lowest order in M 2 G 1 − M 2 Q 1 , i.e. the mass splitting can be ignored inside the one-loop integrals.
The phase space integration for the two-body process can be performed analytically and leads to the decay width The renormalization is performed by using the on-shell scheme for the wave-function renormalization of the physical states and MS renormalization for the strong coupling constant. An example loop diagram is shown in Fig. 1

(B).
For the external states the wave function renormalization constants are introduced for the left-and right-handed (massless) SM quarks, the doublet and singlet KK-quarks and the coloron, respectively. As usual, their values are determined through the residues of the renormalized propagators, leading to and We use Σ L,R (p 2 ) and Σ C (p 2 ) to denote the left/right-handed quark self-energies and the transverse coloron self-energy, respectively. Furthermore we renormalize the coupling constant in the MS scheme, leading to a redefinition of the coupling according to where the beta functions can be written as Note that for our analysis we take all SM quarks to be massless, including the top, such that n q = 6. The scale-dependence of the coupling is determined by fixing its value to the SM value of the strong coupling at µ = R −1 , even though their counterterms are differing. Beneath the threshold µ = R −1 we assume the coupling to effectively not run at all. This is achieved through the logarithm term in (18), which cancels the leading-order running of α s (µ) for µ < R −1 . A similar decoupling prescription was used previously in Ref. [17] in the context of supersymmetry.
To remove the infrared divergencies we employ phase-space slicing with two cutoffs as it is described in Ref. [15]. The soft-collinear divergencies from the real emission diagrams, see Fig. 1 (C) for an example, exactly cancel the soft poles of the one-loop functions. To analytically extract and cancel the poles the matrix element and phase space of the threebody decay are expanded for small emission angles θ and small gluon energies E 3 . This approximation is used in the region bounded by the soft-collinear conditions with soft and collinear cutoff parameters δ s and δ c respectively. Both the matrix element and the phase space of the three-body decay factorize in the limit of small gluon energies into the two-body result and an eikonal factor Figure 2: Full LO and NLO decay width for the level-1 KK-gluon decaying into a regular quark and the level-1 KK-quark within the Coloron model. Both curves exhibit the same Λ dependence, since the κ-factor is constant in our approximation. We assume R −1 = 1 TeV .
The remaining angular integrals can be performed analytically in dimensional regularization and are found in the literature, e.g. in Ref. [16].
Hard-collinear divergencies only arise from the diagram with a final-state gluon being emitted from the SM quark leg. The phase space and the three-body matrix element factorize in the limit of a small angle θ between quark and gluon into the Born contribution and a divergent Altarelli-Parisi splitting kernel. In dimensional regularization one finds where A c 1 and A c 0 are numerical constants found e.g. in [15]. The three-body phase space integration in the hard, non-collinear regime can be performed analytically, within the bounds set by the soft-collinear conditions. After carrying out the integration over the three-body matrix element the cancellation of cutoff parameters between the three and two-body processes is then checked explicitly. Due to the approxima- is a constant, independent of any of the masses, and turns out to be κ = 1.238. In Fig. 2 we plot the dependence of the LO and NLO decay width on the cutoff scale Λ, that was introduced through the radiative mass splitting between the KK-quark and gluon.
The calculation has been performed using publicly available computing tools sublimated by in-house routines. To generate the necessary diagrams and amplitudes we incorporated the MUED Feynman rules into FeynArts 3 [18]. The color, Dirac and Lorentz algebra was performed with FeynCalc [19].
To simplify the treatment of tensor loop integrals, the one-loop amplitude was contracted with the Born amplitude and the sum over the spins of external particles carried out before any tensor reduction. As a result, most tensor structures in the numerator of the loop integrand can be cancelled against propagator denominators. For the remaining tensor integrals, Passarino-Veltman reduction has been used [20]. One thus arrives at a final result in terms of standard one-loop basis functions. A similar procedure has been employed for the calculations in the following sections. The IR-divergent basis integrals necessary for the Coloron decay were taken from Ref. [21].

Contribution of higher KK modes to KK-gluon decay
When calculating the contributions from modes with n ≥ 2 to the decay G 1 → Q 1q , one needs to sum up diagrams with KK modes of level n in the loops up to some order N , both in the vertices and in the counterterms. N will later be identified with the cutoff scale ΛR. The wave function renormalization is performed in analogy to the previous section, but now including level-n modes in the loops. Then the remaining UV divergence is absorbed into a MS coupling renormalization of the form from which one then can extract the contribution β n of the level to the overall beta function

Standard Model Vertices
While we are ultimately interested in the decay of the KK-gluon, we will first study the impact of level-n KK modes on the SM QCD vertices, in order to understand the renormalization procedure. There are three basic QCD vertices, the quark-gluon vertex, the three-gluon vertex and the four-gluon vertex. Example diagrams for each are shown in Fig. 3. All modes within in the loop have the same mode-number starting at n = 0, as a consequence of KK-parity.
The standard result for the SM QCD beta function is reproduced by the diagrams with zero modes (n = 0) at one loop and is given by When allowing KK modes with n ≥ 1 in the loop, an additional contribution to the beta function is generated each time when one crosses a mass threshold n R β n = 7 2 These contributions to the beta function lead to a running of the strong coupling constant that deviates from the Standard Model behavior, described for a renormalization scale be- where the sum over KK-modes extends up to the UV-cutoff scale Λ. It has to be noted that the sum is strictly speaking only defined for an integer cutoff number, i.e. ΛR in (27) is understood as the argument of the Gauss floor function ΛR , rounding down to the largest integer not exceeding ΛR.
Since we take all external legs to be massless there is no dipole operator being generated in the Lagrangian. We find that the vector coupling gqq as well as the three and four-point gluon couplings receive the same Wilson coefficient C SM . Summing over all modes within the loop up to the cutoff scale that coefficient then reads Notice the that the first term in the bracket is independent of the renormalization group running and represents a threshold correction which is not contained in the beta function. For large cutoff scales, and setting µ = Λ, the Wilson coefficient has a leading linear dependence on Λ and can be well approximated by For the computational details of the approximation see appendix A.

KK-Vertices
In the case of two external legs being KK-particles the one-loop corrections to the vertex induce both a vertex-like and dipole operator and the total interaction under consideration becomes where q is the KK-gluon momentum. If only SM modes and the first KK-mode are allowed in the loop, the coupling renormalization of eq. (18) contributes to the beta function with β Coloron . The contributions to the vertex function from higher modes are shown in Fig. 4 and can be straightforwardly calculated for every new mode allowed in the loop. The renormalization procedure is again done level by level and requires the on-shell field renormalizations for the external legs. The beta function resulting from the coupling renormalization at every level n ≥ 2 is identical to the one we found for the SM vertices for n ≥ 1, reported in eq. (26).
For the Wilson coefficient of the vector coupling of eq. (30) we separate the NLO contribution into terms proportional to C A , C F and n q and find 10 + (1 + 5n) log n 2 (n + 1) 2 (31) 15 + 4n + 4n 2 + n 8 + 3n + 2n 2 log n 2 (n + 1) 2 − 7 log Additionally we can calculate the dipole coefficient, which is finite and does not require renormalization: As before we set µ = Λ and can expand the Wilson coefficient in the limit of large Λ and find So both the coefficient of the SM operators and the KK operators show the same leading order behavior, although the log ΛR terms differ between eqs. (29) and (33). Additionally we find for the dipole moment The dipole operator therefore is well-defined even in the limit Λ → ∞.

Flow Equation Analysis
The previous results can be put into perspective through a renormalization flow analysis of the 5D uncompactified theory. Our solution strategy of the one-loop integro-differential equation is a non-local heat kernel expansion and closely follows Ref. [24].
The effective action of a theory in saddle-point approximation can schematically be written as where S is the action and Φ represents a generic field. We also introduced coarse graining through an explicit regulator function R k designed to suppress the infrared modes of the theory. A standard choice is Litim's regulator [22], which is given by in terms of the Heaviside theta function Θ. It is optimized in the sense described in [23]. Taking the derivative with respect to the renormalization time parameter t = log k k 0 yields the one loop flow equation, which can be written symbolically as This is the one-loop approximation can be obtained from the full functional renormalization group equation [12], by fixing the effective action to its UV limit Γ Λ = S. On the right-hand side of this flow equation we can identify a heat kernel h k which allows for a non-local expansion [25,26]. To lowest order and in the absence of gravity the expansion of the flow equation reduces to where Ω M N = [D M , D N ] is the 5D gauge connection of the theory and U = D 2 1+∆ the nonderivative part of the Laplacian, which carries all representations of the theory. Furthermore the expansion depends on coefficient functions given by which in turn can be expressed through the "Q-functionals" that take the form which can be found analytically for the optimized regulator [26].
where the first term does not contain any fields and will therefore be ignored. The nonderivative part U of the Laplacian and the gauge connection Ω M N are in this scenario given by Armed with these tools we can apply the non-local heat kernel expansion of eq. (38) to the right hand side of the QED flow equation and find ‡ This form of the flow equation can now be compared to where the polarization function Π k (z) depends on z = −D 2 . This leads to an integrodifferential equation for the wave function renormalization constant Z k and the polarization function To extract the beta function from this equation we consider the limit of zero momentum transfer, z → 0, and exploit the fact that Π k (0) = 0. In that limit the second term of the equation behaves like We can now find the QED-like contributions to the anomalous dimension of the gluon field as well as the beta function, after identifying ∂ t Z k = 2∂ t (g (5) k ) −1 for the 5D coupling, leading to in terms of the 4D effective coupling g = √ 4π g (5) . Note that the MUED model contains both a SU (2) singlet and doublet quark for every flavor, so that an additional factor 2 must be included when comparing with section 4.1.
Non-Abelian QCD contributions: For QCD in D dimensions one can proceed in a similar fashion. The contribution proportional to C A can be most conveniently extracted from the Yang-Mills part of QCD in MUED. The tree-level action in background field gauge reads where D M denotes the gauge covariant derivative with respect to the background field, while D M contains the full field. c a is the corresponding Faddeev-Popov ghost. The one-loop flow equation, applied to the above action, can be brought into the standard heat kernel form and is taken from Ref. [27]: After applying the non-local heat kernel expansion to the first term of the above flow equation and comparing it to the expected form of the effective action we are again left with When taking the limit z → 0 one finds the non-Abelian QCD contribution to the anomalous dimension of the gluon field as well as to the beta function

Phenomenological analysis
From eqs. (47) and (52) one finds that the cutoff dependence of the effective gauge coupling is given by Since this was obtained in the framework of an uncompactified 5D theory, there is no distinction between couplings of zero modes and higher KK modes. The result in eq. (53) is numerically close, but not identical to the linear ΛR terms in eqs. (29) and (33). The difference stems from threshold corrections that contribute to (29) and (33) at each KK level. Thus, even when we study the dependence of the effective vertex coupling on Λ for ΛR 1, i.e. at scales much beyond the compactification radius, there is a non-vanishing impact of the compactification. In other words, when considering MUED at large scales one does not trivially recover the uncompactified 5D theory. This is illustrated in Fig. 5. In this plot it is assumed that the there is some fundamental 5D coupling at the scale Λ 0 = 20R −1 . The curves illustrate the effect of integrating out all energy scales beyond ΛR −1 in different approximations. The solid blue line corresponds to the result (53) from the flow equation. The solid green curve depicts the result from the explicit calculation of the gqq vertex diagrams, as in eq. (28) with µ = Λ, whereas the solid orange curve is the approximation (29). The dotted green and orange curves are the equivalent cases for the G 1 Q 1q vertex, eqs. (32) and (33). The vertical axis is normalized such that the effective gqq coupling at n = 1 is equal to its SM value.
First of all, the plot demonstrates that the large-ΛR expansions (orange) provide an excellent approximation to the level-by-level calculations (green). Furthermore, it shows that all approaches lead to similar results for the dependence of the effective gauge couplings on the truncation scale. However, there are difference at the few-percent level, which stem from the following two facts: The flow equation is insensitive to threshold corrections from the compactification; and the G 1 Q 1q and gqq vertices have different logarithmic dependencies on Λ. Note that the plot also includes the running from the coloron contribution, see eq. (18).  Figure 5: Running of the SM vertex function q − q − g (solid line) and the KK vertex q − Q 1 − G 1 (dashed line). We compare the full summation up to the cutoff with the asymptotic expansion as well as the results obtained from the FRGE running. The integration constants have been fixed such that the couplings coincide at ΛR = 20.
From a low-energy perspective, on the other hand, the dependence of the gauge vertices on the cutoff Λ is not directly observable. However, if KK gluons and KK quarks should be discovered in the future, one can compare the strength of the SM gqq coupling with the G 1 Q 1q vertex. While the leading linear ΛR dependence drops out in this difference, there is still some sensitivity to the cutoff scale from the logarithmic terms in eqs. (29) and (33): Including the contribution to the running from the lowest modes, see eq. (18), and assuming that the gqq and G 1 Q 1q couplings are identical at the scale µ = Λ, one finds that the observable couplings depends on Λ as follows: For ΛR = 10 . . . 50, this amounts to an effect between 3.5% and 6%. This means that the prediction of the decay rate for G 1 → Q 1q , including O(α s ) corrections, depends on the unknown UV completion of MUED. At the same time, the impact of this UV sensitivity is rather mild, at the level of a few percent, which may be negligible for most practical purposes.
While a more complete analysis of different processes in MUED and larger classes of extradimensional models is beyond the scope of this work, we expect that a similar conclusion can be reached in these cases. This conjecture is based on our observation that the leading UV sensitivity cancels when normalizing the KK gauge boson vertex to the SM gauge boson vertex, as dictated by the renormalization flow equation. The cancellation still works when including the threshold corrections at each KK level, and thus it should hold in any extradimensional extension of the SM. The next-to-leading term, while enhanced by log(ΛR), is nevertheless numerically rather modest.

Conclusions
Models with extra dimensions generally feature infinite KK towers of new states that need to be considered when calculating loop corrections within those models. KK-number conserving operators receive corrections from the full spectrum of these modes. The corrections are in general dependent on a cutoff scale Λ at which the model breaks down and some unspecified UV completion is required to describe the physics. In this article we investigated the numerical impact of this unknown cutoff parameter on physical observables. As a concrete framework we considered QCD in a spacetime with one additional universal extra dimension, which is compactified on a circle with a Z 2 orbifold.
In the first half of this article, the problem was discussed from the viewpoint of the compactified 4D effective theory. We computed the full one-loop QCD corrections to the SM gauge-boson vertex gqq as well as the vertex G 1 Q 1q (or, equivalently, G 1Q1 q), which contains two level-1 KK modes. These include vertex diagrams, on-shell counterterms for the external legs, and MS coupling counterterms generated by loops up to KK level n. The MS scale dependence of the gqq vertex can be described through the regular QCD beta function, where the well-known SM result is supplemented by an extra term each time one of the KK thresholds is crossed. Additionally, the Wilson coefficient of the gq/q vertex receives a finite threshold correction from each KK level, which is not described by the beta function.
For the G 1 Q 1q (G 1Q1 q) vertex we proceeded similarly. For concreteness, we considered the physical process G 1 →qQ 1 . The NLO corrections to this decay exhibit soft divergencies in the virtual vertex contributions, which can be cancelled against the real radiation contributions with the two cutoff phase-space slicing method. This part of the calculation is identical to a 4D coloron model. Overall, the beta function for this vertex differs from the SM gq/q vertex, but the contribution from higher KK modes (n > 1) is identical for the gqq and G 1 Q 1q beta functions. Furthermore, the Wilson coefficient for the G 1 Q 1q vertex receives KK threshold corrections, which differ from the ones found for the SM vertex.
The results for both vertices can be conveniently written in terms of an expansion for large values of the cutoff scale Λ. It turns out that their leading terms, which are linear in ΛR, are identical. The first difference between the two vertices appears at the subleading order log (ΛR). Higher orders beyond the log term are numerically very small and can be neglected.
In the second half, we studied the cutoff dependence within uncompactified 5D QCD, using the functional renormalization group flow equation. To solve the equation and extract the beta function we applied the saddle-point approximation and utilized the non-local heat kernel expansion method. Through this approach we were able to find the flow equation analogue of the one-loop beta function. Its coefficient coincides with the leading order behavior of the running of the coupling found for the explicit diagrammatic calculation in the compactified theory. However, the uncompactified 5D framework is not able to reproduce the contributions from the threshold corrections in the 4D framework.
When comparing the diagrammatic 4D calculations of the gqq and G 1 Q 1q vertices, as well as the 5D flow equation result, one finds that the NLO prediction for the decay G 1 →qQ 1 is indeed sensitive to the choice of the UV cutoff and thus to the unknown high-scale physics. However, the numerical impact of this uncertainty is numerically rather modest, since the leading contributions cancel in the comparison.

A Asymptotic Expansion
The analytic summation using elliptical Jacobi theta functions, that has been employed in Ref. [3,6,28], is not an accurate method for our case, since we are defining our coupling renormalization within the MS scheme. We are however able to recover the same results in the large ΛR limit, as the theta function method does within the on-shell scheme, using a somewhat different approach. We begin by performing the summation up to mode N = ΛR analytically. The only non-trivial sums appearing are of the form which can be carried out using the Riemann ζ function. The prime denotes a derivative with respect to the first argument of the generalized ζ function ζ[s, a] = ∞ n=0 (n + a) −s . One solves the sums involving higher moments analogously, like e.q.: where A = 1.282427 · · · denotes the Glaisher-Kinkelin constant.
For the asymptotic expansion of the first sum, in a region where ΛR becomes large we simply utilize Stirling's Formula For the sums involving higher powers of n, the asymptotic expansion of the derivative of the generalized Riemann ζ function can be performed according to It is worth noting at this point that the Riemann ζ function is related to the elliptical Jacobi theta function through making our order by order recovery of the results derived with the aide of the theta function not entirely surprising.