Analytic expressions for Debye functions and the heat capacity of a solid

Analytic expressions for the $N$-dimensional Debye function are obtained by the method of brackets. The new expressions are suitable for the analysis of the asymptotic behavior of this function, both in the high and low temperature limits.


Introduction
The N dimensional Debye functions play an important role in study of a variety of problems in statistical physics and solid state physics, especially in calculations of heat capacity of solids. This function appeared first in a model proposed by Debye [5] describing the heat capacity of a crystalline solid, which with some variations it is still used today.
These functions, up to recent times only known through their integral representation, have created enough interest in their evaluation for arbitrary values of N and absolute temperature T . Debye functions can be expressed as [1, ch. 27]: for |X| < 2π and N ≥ 1. Here B k are the Bernoulli numbers. Numerical computation of these functions appear in [19,20]. The first analytical expression for such functions, other than integral representations, may be found in [6]. The work presented here gives expressions for D N (X) in terms of the polylogarithm functions.
The method of brackets [11] is used in the current work to evaluate the Debye functions. This method gives (new) analytic expressions for them, reproducing the results presented in [6], as well as some new expressions.
The method of brackets was developed in the context of calculations of multidimensional definite integrals appearing in evaluation of Feynman diagrams [10,14,15,16]. It consists in a small number of heuristic rules that yield the evaluation of a wide range of integrals. These rules admit an easy implementation in a computer algebra system. The reader will find more details in [2,9,11,12,13].
The content of the paper is described next. Section 2 introduces the method of brackets, Section 3 uses this method to evaluate Debye functions and their analytic expressions. In particular, expressions that are free of integral representations are presented here, recovering those presented by [6]. These results are then used to study the asymptotic behaviour of these functions in limiting values of the temperature. Section 4 uses these representations to evaluate the internal energy and heat capacity in solids. The emphasis here is on the new expression for Debye functions to show that the manipulation of them simplifies the computation of the limits T → 0 and T → ∞ for temperature presented in [6].

Basic of Method of Brackets (MoB)
The method of brackets is a generalized version of the Negative Dimensional Integration Method (NDIM) [3,4,7,18,21], a technique developed to evaluate Feynman diagrams. In quantum field theories, Feynman diagrams correspond to multi-variable integrals that represent physical processes.
This method evaluates definite integrals in one or several dimension over the interval [0, ∞]. The procedure introduces the notion of a bracket and converts the integrand in a series of brackets. The method contains a small number of heuristic rules which transform the evaluation of an integral into the solution of a small linear system of equations. A summary of these rules is presented below. More details may be found in [2,9,11,13].
Rule 0. For a ∈ C, the bracket associated to a is the divergent integral Rule 1. The expansion of an arbitrary function. The use of the method of brackets requires to replace components of the integrand by their his corresponding power series, that is, it is required to represent an arbitrary function f (x) as: where C(n) are the coefficients in the expansion, α and β are arbitrary (complex) exponents and φ n is defined by: .
For multidimensional integrals one needs expansions in several variables, such as The notation φ 12 is frequently used for φ n1 φ n2 .
Rule 2. Polynomial expansion. An expression of the form (A 1 + · · · + A r ) µ often appears in the evaluation of integrals. The expansion This rule has been established in [11].
Rule 3: Eliminating integration symbols. Once the first two rules are applied, the integral is converted into a bracket series using the definition of bracket.
Rule 4: Finding solutions. The result of applying the previous rules to an integral is that its value is represented by a bracket series J. The rule to evaluate this series is given in the special case when number of sums and brackets is the same (this is the so-called index zero case): the bracket series is .φ nr C(n 1 , · · · , n r ) a 11 n 1 +· · ·+a 1r n r +c 1 · · · a r1 n 1 +· · ·+a rr n r +c r .
The coefficient C(n 1 , ..., n r ) depends on the parameters of the integral and the index of the sum {n i } , i = 1, ..., r. The value of this multiple sum is declared to be where A = {a ij } and the values {n * i } (i = 1, ..., r) are the solutions of the linear system obtained by the vanishing of the brackets: If the matrix A is not invertible and the number of sums is larger than the number of brackets, there is an extension of the procedure described here to evaluate the integral. Details may be found in [11,13].
The Debye functions is defined by: The following extension is considered here: Here N is zero or a positive integer, X and α are positive parameters. The parameter α is introduced here to find alternative expressions for these extensions.
3.1. A bracket series for D N (α, X). The computation of a bracket series for D N (α, X) is described next. The first step is the expansion of the denominator in the integrand to obtain: The expansion of the exponential function is and replacing in (3.3) produces The change of variables y = t/(X − t) converts the last integral to [0, ∞) as (y + 1) N +n3+2 dy, and the desired bracket series of bracket is The final bracket series for D N (α, X) is × 1 + n 1 + n 2 N + n 3 + 2 + n 4 + n 5 N + n 3 + n 4 + 1 .
An expression for the integral (3.2) is now obtained from (3.6). The method of brackets yields four different series: The influence of the parameter α is discussed first, because in addition to parameter X, it allows to discriminate different series. The four solutions S j are power series in α or 1/α, so that S 1 and S 4 are expansions in α → ∞ and S 2 and S 3 are expansions in α → 0. The same situation occurs with respect to the parameter X. Each series represents the integral (3.2). Their analysis is described next.
(1) The series S 4 must be neglected because the term with n 1 = 0 diverges.
(2) The series S 3 is naturally truncated at n 2 = 0. Since this index is associated to the powers of X −1 , it represents an asymptotic approximation for case X >> 1. A detailed study including condition α → 1 yields: where ζ(s) is the Riemann zeta function. (3) The series S 1 and S 2 are both convergent as power series in X. Both are expressions for D N (α, X), but it turns out that they are equivalent.

3.2.
Analysis of the expressions obtained above.

S 1 as solution.
Rearranging the defining series produces a hypergeometric representation: The previous expression may be written as where hypergeometric function 1 F 1 is the Kummer function. Now use (3.14) where γ(n, Z) is the incomplete Gamma function defined by the integral representation In the important special case of n ∈ N, the function γ(n, Z) can be written as a finite sum and then The formula (3.14) is now transformed to (−Xn 1 ) k k! , and the series S 1 can be written as After some algebraic manipulations, the previous expression is written as The polylogarithm function [17], defined by the series is now used to obtain an expression for the Debye function D N (α, X) in the form (3.19) This formula was first presented in [6]. In addition to this representation, the method of brackets produces a new expression for the Debye function using S 2 .

The series S 2 . A new solution.
As in the computation of S 1 , the series defining S 2 can be written as a sum of values of the incomplete Gamma function: (1 + n 1 ) N +1 γ (N + 1, (1 + n 1 ) X) , and using (3.17), this becomes Proceeding as in the previous case, the Debye function D N (α, X) is now In summary, the method of brackets has produced two equivalent formulations of the representation of the Debye function given in [6]. The first one in (3.2), reproducing the solution presented in [6] and a second expression given in (3.22). This is a new representation for D N (α, X).

Application : Debye Model and heat capacity in solids
An important topic in solid state physics is the determination of heat capacity using quantum treatments [5,8]. The integral expression (3.1) is associated to this problem through a model proposed by Debye [5]. According to this model, the internal energy in solids is given as a function of the absolute temperature T , by with the usual notation for Debye functions, i.e D 3 ΘD T = D 3 1, ΘD T . Here k B is the Boltzmann constant, Θ D is called the Debye temperature and N the number of particles in the system.
Using (3.19), and with the notation u = Θ D /T , the Debye function is The analysis of (4.2) as T → 0 is not easy to obtain directly from here. On the other hand, the new expression (3.22) permits such an analysis. To describe this procedure and with the same notation as before, start with This is described next.

Asymptotic limits.
The classical approach to study limiting behavior of these functions is to conduct approximations for the integral representations, valid in some specific limits (high and low temperatures). These limits can now be studied directly from analytical expressions presented here. The new formulae presented here permit the analysis of limiting high and low temperatures. This study reproduces the results of [19,20]: • At T → ∞, In the analysis of this last formula, the behavior of the polylogaritmic function Li n e −ΘD/T ≪ 1 as T → 0 is used. This can be seen from the power series expansion These are in agreement with the result cited in the literature [19,20].

4.2.
Heat capacity. This is computed using c V = ∂U ∂T V . The limiting behaviors are • For T → ∞ • For T → 0 The analytical expressions for the Debye functions presented here produce results valid for arbitrary temperatures. Using (4.1), and with the notation u = Θ D /T , the value c V is given by

Conclusions
Analytic expressions for the Debye functions have been produced using the method of brackets. These expression differ from the classical integral representations and they involve sums of the polylogarithm function. One of the results presented here reproduces formulas developed in [6].
The new expressions obtained here provide an efficient way to evaluate high and low temperature behavior.