Rational Relation Conjectures

Dustin Steinhauer, Randy Qian^{∗}

February 12, 2003

Mathematics Department Princeton University Princeton, NJ 08544

Abstract

We will look computationally for simple rational factor relationships between Mahler measures of specific polynomials and L-functions using techniques involving the LLL basis reduction algorithm. We will also look for rational relationships between ζ(2n + 1) and an experimentally postulated, rapidly converging formula which is of significant interest to us.

1 First Section: Using Continued Fractions to “Guess the Rational Number”

Throughout this paper we will be attempting to form rational linear dependencies between irrational (and often transcendental) quantities. A natural place to start is with the problem of determining a single rational number from a computer output. More precisely, if we are given the first several digits of the decimal expansion of a rational number, it would be useful to be able to derive the rational expression (i.e., p/q with integers p and q such that (p, q) = 1) to which it corresponds. By definition, any finite decimal expansion represents a fraction of the form x/10^{a},where x is some integer and a is the number of digits past the decimal point that we are given, so we want to find a rational expression “better” than that one. In particular, we are interested in having fewer digits in the denominator.

We note that we can take the continued fraction expansion of a rational number and it will clearly terminate since any finite continued fraction corresponds to a rational number and vice-versa. That a finite continued fraction equals a rational number is clear from simply carrying out the finite number of arithmetic operations involved. If we begin with a rational number 0 < p_{1}/p_{0} < 1, we can carry out Euclid’s Algorithm to generate the continued fraction:

∗E-mail: dsteinha@princeton.edu, rqian@princeton.edu

\frac{p_{1}}{p_{0}} = \frac{1}{\frac{p_{0}}{p_{1}}} = \frac{1}{q_{0} + \frac{p_{2}}{p_{1}}}

where p_{0} = p_{1}q_{0} + p_{2} and p_{1} > p_{2} ≥ 0 is guaranteed by long division property of the integers. So we now expand \frac{p_{2}}{p_{1}} in the same manner as above and repeat the results it generates. At an arbitrary step in the process, we have:

\frac{p_{n+1}}{p_{n}} = \frac{1}{\frac{p_{n}}{p_{n+1}}} = \frac{1}{q_{n} + \frac{p_{n+2}}{p_{n+1}}}

so we have p_{n} > p_{n+1} ≥ 0 for any n. This means that the p_{n} form a strictly decreasing sequence bounded by 0. Since p_{n} is always an integer, we will take at most p_{0} steps to reach 0, where our continued fraction terminates, and is therefore finite.

We can also take the continued fraction expansion for the finite decimal expansion in question, and given enough digits, its continued fraction expansion will bear the same terms up to a certain point as the rational number. So to find a rational number from a decimal expansion, we could cut off its continued fraction expansion at some point and unravel it into a normal fraction. It is clear then that we should look for an uncharacteristically large term Q in the continued fraction since it would correspond to a very small correction to the value of the continued fraction. If we examine how a large partial quotient Q affects a much smaller previous partial quotient q, we see that … \frac{1}{q+1/Q} ≈ \frac{1}{q} …, so q alone might suffice since it is only slightly adjusted by 1/Q. We can then reasonably assume that the value of the continued fraction with all the terms before that point would be very close to the rational expression we want, and unravel that to get our normal fraction.

For example, suppose we are given the digits .005869134058 … and told that these are the first 12 digits of a rational number. What is our “best guess?” The partial quotients of this number are

[0, 170, 2, 1, 1, 1, 1, 2, 1, 4, 5, 1, 9758, 7, 4, …].

It is reasonable to suspect that from some point on, a continued fraction representation of this decimal will represent the decimal’s deviation from the rational number rather than the rational number itself. If the denominator of the rational number is much smaller than 10^{12} then we have a better guess than the obvious guess if we cut off the continued fraction expansion just before the first large partial quotient. In fact, the continued fraction

[0, 170, 2, 1, 1, 1, 1, 2, 1, 4, 6]

equals a rational number with a particularly small denominator, \frac{1379}{234958}.

2 Second Section: Using the LLL Algorithm to “Guess the Rational Numbers”

The method for approximations described above works well when we want to solve numerical equations that are essentially first degree polynomials, namely, we looked for a rational solution to x in ax + b = 0 given sufficiently many digits of a/b. For higher order polynomial, a more general and different technique using a lattice reduction algorithm called LLL can be used. We will rely on mathematica’s implementation of LLL in the form of its LatticeReduce function which takes a matrix as its input. It takes the basis formed by the rows of the matrix and finds a new basis of “shorter” vectors in the form of another matrix whose rows spans the same space and are linear combinations of the original basis vectors [W1].

We note that the powers of x in a polynomial p(x) = a_{n}x^{n}+…+a_{1}x+a_{0} can be thought of as a set of basis vectors that span the space of polynomials of degree n or less. The polynomial is then the vector (a_{n}, …, a_{0}). If we are given a single root to a certain precision, we use a trick to pack the powers of this root into the rows of a matrix of proportionate and much larger integers. If we then reduce the basis spanned by the powers of this root, we find a linear combination which takes this basis to a much “shorter” value, namely the magnitudes of the vectors in the new basis are much smaller in magnitude than the a_{i}. Since the linear transformations are also encoded within the new basis of the matrix, they themselves are shorter, and these small valued linear combinations precisely describes a “nice” polynomial since it would map the approximate root very close to 0.

The trick is to take an identity matrix and append a final row which contains all of the powers of the root up to the degree we wish to test and multiply them all by 10^{n}, where n is the number of significant digits past the decimal point that we have calculated. We run the basis reduction on the matrix and we note that the entries in each row of the resultant matrix, save its last entry, correspond to the coefficients of some potential polynomial. The size of the last term can be seen as a measure of how well the polynomial fits to our root.

Example 1:

To illustrate this discussion, we will try to find a “nice” polynomial which has the root a = 5^{1/3} + 2^{1/4}. We know that this expression will require at most a degree 12 polynomial because we adjoin the numbers 5^{1/3} and 2^{1/4} to the ground field of the rationals. This constitutes a degree 12 extension, and hence there is a linear relationship between the 13 vectors of a degree 12 polynomial.

We can calculate 5^{1/3} and 2^{1/4} to 6 digit precision easily, and we find that a = 2.89918. We multiply this value by 10^{6} to get an integer and pack a matrix with the powers of this value adjoined to the identity matrix like this:

( I | (10^{6} ∗ a)^{0} ; ⋮ ; (10^{6} ∗ a)^{12} )

where I corresponds to the 13 X 13 identity matrix. After lattice reduction is performed, we end up with a matrix that looks like:

( −3 2 2 −2 1 1 −1 −2 −1 −2 −2 1 0 ; ⋮ | 2 ⋮ )

All of our final column entries are “short” compared to the original vectors so it is fair to just consider the first row since it is as good a candidate as any other row. Again, this is a linear combination of our original rows. Since we packed the left side with the identity matrix, the operations performed to reach this result are obviously just row entries to the left of the final column. So our polynomial is then

P (x) = −3 + 2x + 2x^{2} − 2x^{3} + x^{4} + x^{5} − x^{6} − 2x^{7} − x^{8} − 2x^{9} − 2x^{10} + x^{11}.

Solving P (x) = 0 indeed yields a solution which has the desired solution to at least the precision which we have specified.

3 Third Section: Mahler Measure

Let P (x_{1}, x_{2}, …, x_{n}) be a polynomial in n variables. Define the logarithm of its Mahler measure to be:

m(P) = ∫_{0}^{1} ⋯ ∫_{0}^{1} log |P (e^{2πix_{1}}, …, e^{2πix_{n}})| dx_{1} ⋯ dx_{n}.

and we often denote the actual Mahler measure by

M(P) = e^{m(P)}.

Lemma 1: Let f(x) be a monic polynomial in 1 variable,

f(x) = x^{n} + … + c_{1}x + c_{0} = ∏_{j=1}^{n}(x − a_{j})

where the a_{j} are the roots of the polynomial. If the roots a_{1}, …, a_{k} are on or inside the unit circle, then

M(f) = ∏_{j=k+1}^{n} |a_{j}|.

Proof: Let r > 1 be such that |a_{j}| < r for j ≤ k and |a_{j}| > r for j > k. Jensen’s formula states that

∫_{0}^{1} log |f(re^{2πix})| dx = log |f(0)| + ∑_{j=1}^{k} log |\frac{r}{a_{j}} | .

Note that the left side is m(f) as r → 1^{+}. Exponentiating both sides and taking the limit as r → 1^{+} we see that

M(f) = |f(0)| ∏_{j=1}^{k} |\frac{1}{a_{j}} | .

But f(0) = c_{0} = (−1)^{n}a_{1} ⋯ a_{n}, so the roots inside the circle cancel, and

M(f) = ∏_{|a_{j}|>1} |a_{j}|.

A cyclotomic polynomial has measure 1. It is not known what the next largest possible measure of a polynomial is, but Lehmer’s conjecture states that the polynomial

l(x) = x^{10} + x^{9} − x^{7} − x^{6} − x^{5} − x^{4} − x^{3} + x + 1

has the next smallest possible measure. M(l) = 1.1762808183 … which is the value of its sole root outside the unit circle. l(x) is a called a Salem polynomial and its root outside the unit circle is called a Salem number. More precisely, a Salem number is a positive real number greater than one whose minimal polynomial has all of its other roots on or inside the unit circle with at least one pair of conjugate roots on the circle.

Certain results concerning the Mahler measure of polynomials exist that relate their values to other quantities [B]. Many of these relations that have been proven were discovered numerically and subsequently shown rigorously. Furthermore, there are quite a large number of results that have been numerically verified to high accuracy but not proven.

For example, let P (x, y) = x + y + 1. A numerical calculation shows m(P) = .323066 . . .. On the other hand, let χ be the odd Dirichlet character (mod 3), that is, χ(n) = 0 if n ≡ 0 (mod 3), χ(n) = 1 if n ≡ 1 (mod 3), and χ(n) = −1 if n ≡ 2 (mod 3). Note

\frac{m(P)}{L(χ, 2)} = .413497

which was observed to be very close to \frac{3\sqrt{3}}{4π}. Furthermore this can be written more succinctly using the functional equation for L(s, χ), namely,

L(s, χ) = 2^{s}π^{s−1}k^{−s+1/2} cos(\frac{π}{2} s)Γ(1 − s)L(1 − s, χ).

If s = −1 then clearly L(−1, χ) = 0, as the cosine term vanishes and no other term has a singularity. If we then let s = 2, the simple pole of the Γ function cancels the simple zero of L(−1, χ), and m(P) = L′(χ, −1).

There are plenty of unproven conjectures concerning rational combinations of Mahler measure of polynomials. For example, if

P (x, y) = (x + y + 1)(x + 1)(y + 1) + xy
Q(x, y) = x^{3} + (x^{2} + 2x − 1)y + y^{2}

then it appears that 5m(P) = 7m(Q). Furthermore, if the definition of Mahler measure is naturally extended to include Laurent polynomials, we may further observe that

M(x + x^{−1} + y + y^{−1} + 5) = [M(x + x^{−1} + y + y^{−1} + 1)]^{6}

to high precision. This example suggests the following question: Let

P_{j}(x, x^{−1}, y, y^{−1}) = x + x^{−1} + y + y^{−1} + j.

Are there other values of j that produce relationships between polynomials in this family? We will also attempt to test other families of polynomials such as

P_{I}(x, y) = ∑_{I}(xy)^{I}

where I is some set of indices we can vary. For example, we may work with (xy)^{I} = x^{i_{1}}y^{i_{2}} and i_{1} + i_{2} = n for a given integer n. Our numerical mathods, however, were not sufficient to generate many digits of precision; furthermore, they turned out to be too slow to test large numbers of polynomials.

The issue still remains to prove these relations once they are discovered numerically. The number of cases in which a proof has been achieved is far lower than the number of conjectures that have been made using computers. A conjectured relation can be a useful starting point for a proof, however. As was mentioned, in the case of P (x, y) = x + y + 1, its relation to an L-function was indeed discovered numerically. As it turns out, we can prove that close numerical agreement in this case is not a coincidence.

Theorem 2: Let P (x, y) = x + y + 1 and χ be the odd character (mod 3). Then

m(P) = \frac{3\sqrt{3}}{4π}L(χ, 2)

Proof: Define the function log^{+}(x) to be log(x) if x > 1 and 0 otherwise. If we apply Jensen’s formula to a 1-variable monic polynomial with roots a_{1}, …, a_{n}, then we can rewrite the formula for its Mahler measure as

m(P) = ∑_{j=1}^{n} log^{+} |a_{j}|

as only roots outside the disk will be counted. When viewed as a function of y alone, P has one zero at −x − 1. Then

\frac{1}{4π^{2}} ∫_{0}^{2π} dθ ∫_{0}^{2π} log^{+} |e^{iθ} + e^{iϕ} + 1| dϕ

= \frac{1}{2π} ∫_{0}^{2π} log^{+} | − e^{iθ} − 1| dθ

and the problem is reduced to evaluating the integral

\frac{1}{2π} ∫_{−π}^{π} log^{+} |e^{iθ} + 1| dθ.

A moment’s reflection verifies that |e^{iθ} + 1| ≥ 1 if |θ| ≤ 2π/3 so we consider

\frac{1}{2π} ∫_{−2π/3}^{2π/3} log |e^{iθ} + 1| dθ.

|e^{iθ} + 1| = \sqrt{2 + 2 cos θ} and since cos θ = 2 cos^{2}(θ/2) − 1 we get

\frac{1}{π} ∫_{−2π/3}^{2π/3} log (2 cos (\frac{θ}{2})) dθ.

After noting that cos is even and changing variables we get

\frac{2}{π} ∫_{0}^{π/3} log(2 cos θ) dθ
= \frac{2}{3} log 2 + \frac{2}{π} ∫_{0}^{π/3} log(cos θ) dθ.

This final integral is complicated to compute. A math package can provide an answer which we verify by differentiation.

Lemma 2:

∫ log(cos x) dx = x log(cos x) − x log(1 + e^{2ix}) + \frac{ix^{2}}{2} + \frac{i}{2}Li_{2}(−e^{2ix}).

Proof: An integration by parts gives

∫ log(cos x) dx = x log(cos x) + ∫ x tan x dx

so we want to verify that the last 3 terms are the integral of x tan x. The first of these 3 terms differentiates to

− log(1 + e^{2ix}) − \frac{2ixe^{2ix}}{1 + e^{2ix}} .

The dilogarithm cancels the first term, as

\frac{d}{dx} \frac{i}{2} ∑_{k≥1}(−1)^{k} \frac{e^{2ikx}}{k^{2}} = ∑_{k≥1}(−1)^{k+1} \frac{e^{2ikx}}{k} = log(1 + e^{2ix}).

We are left with

\frac{d}{dx} \frac{ix^{2}}{2} − \frac{2ixe^{2ix}}{1 + e^{2ix}} = ix (1 − \frac{2e^{2ix}}{1 + e^{2ix}} ) = ix\frac{e^{−ix} − e^{ix}}{e^{−ix} + e^{ix}} = x tan x.

Continuing with the proof of the theorem, note that ∫_{0}^{π/3} log(cos x) dx is clearly real so we can ignore the imaginary parts of its antiderivative. Recalling that Re(log z) = log |z| and again using the fact that |e^{iθ} + 1| = 2 cos(θ/2), we find that

\frac{2}{π} (θ log(cos θ) − θ log(1 + e^{2iθ}))|_{π/3} = \frac{2}{3} log (\frac{1}{2})

which cancels the \frac{2}{3} log 2 from before. We have succeeded in showing that

m(P) = Re(\frac{i}{π} Li_{2}(e^{−πik/3})) = \frac{1}{π} ∑_{k≥1} \frac{sin(πk/3)}{k^{2}}.

Since sin(π(k + 6)/3) = sin(πk/3) we obtain the sum

\frac{\sqrt{3}}{2π} (\frac{1}{1^{2}} + \frac{1}{2^{2}} − \frac{1}{4^{2}} − \frac{1}{5^{2}} + …) .

Observe that

(\frac{1}{1^{2}} + \frac{1}{2^{2}} − \frac{1}{4^{2}} − \frac{1}{5^{2}} + …) − (\frac{1}{1^{2}} − \frac{1}{2^{2}} + \frac{1}{4^{2}} − \frac{1}{5^{2}} + …)
= 2 (\frac{1}{2^{2}} − \frac{1}{4^{2}} + \frac{1}{8^{2}} − \frac{1}{10^{2}} + …)
= \frac{1}{2} (\frac{1}{1^{2}} − \frac{1}{2^{2}} + \frac{1}{4^{2}} − \frac{1}{5^{2}} + …)

so we have

m(P) = \frac{3\sqrt{3}}{4π} (\frac{1}{1^{2}} − \frac{1}{2^{2}} + \frac{1}{4^{2}} − \frac{1}{5^{2}} + …) = \frac{3\sqrt{3}}{4π}L(χ, 2)

and the proof is complete.

In general, given a character χ and a number s such that Re(s) > 1, L(s, χ) can be calculated from its Dirichlet series. It would be desirable nevertheless to have a better way to calculate L(s, χ). Such a method can be obtained from the functional equation for L(s, χ). If χ is an odd character (mod k), then as was mentioned before,

L(s, χ) = 2^{s}π^{s−1}k^{−s+1/2} cos (\frac{π}{2} s) Γ(1 − s)L(1 − s, χ).

Though this formula describes the analytic continuation of L(s, χ) to the complex plane, for our purposes it can help us compute L(n, χ) for n a positive integer > 2. We will return to this point later.

4 Fourth Section: ζ(2n + 1) Conjectures

It has been conjectured that ζ(2n+1) can be expressed as a formula involving linear rational combinations of S_{±}(2n + 1) and π^{2n+1} [W2]. Here we define S_{±}(n) to be

S_{±}(n) = ∑_{k=1}^{∞} \frac{1}{k^{n}(e^{2πk} ± 1)}

.

Such a formula is know to hold for ζ(3):

ζ(3) = \frac{7}{180}π^{3} − 2S_{−}(3)

Using the LLL algorithm, we will be testing many values of ζ(2n + 1) against this conjecture. We will use a rapidly converging numerical method to find many digits of ζ(2n+1) and use this to generate potential coefficients in the formula using the LLL method. Clearly successful use of the LLL algorithm depends on knowing the quantities involved to high precision. The following example on ζ(3) will illuminate some of the difficulties that can arise in the calculations that interest us.

Example 2:

We will now postulate a formula for ζ(3) and discuss the numerical details involved. To do this, we will guess the coefficients for the linear relation a_{1}ζ(3) + a_{2}π^{3} + a_{3}S_{−}(3) = 0. Several numerical difficulties arise immediately. We wish to be able to find a rapidly convergent expression for ζ(3) to determine enough digits when using LLL. We also want to know how many digits of π we need to achieve a certain degree of precision for π^{3}. (Here, S_{−}(3) can be evaluated directly; the expression for its summands involves an exponential term in the denominator, meaning it converges rapidly. The first 10 terms give a result correct to 35 decimal places.) Clearly, the same problems arise when we discuss formulas for ζ(2n + 1) replacing 3 with 2n + 1, with the added difficulty of calculating digits from new dilogarithm expressions that are likely to appear in the exact formulas.

There are well known methods of calculating many digits of π; we just need to know how many digits are needed to verify k digits of π^{n}, which can be done through the binomial equation. Let p be our approximate value of π, and let π = p + ϵ where ϵ is our error term. So

(p + ϵ)^{n} = p^{n} + np^{n−1}ϵ + E_{0}.

Here, E_{0} is a combinatorial expression, and using Stirling’s Formula,

( \begin{matrix}n \\ k \end{matrix} ) ∼ 2^{n}/n^{1/2} we are able to bound this expression by

E_{0} < 2^{n}/n^{1/2}(p^{n−2}ϵ^{2} + … + ϵ^{n}) ∼ 2^{n}p^{n−2}e^{2}/n^{1/2}.

We can be confident that p is at most 3.5, and because ϵ is small, ϵ^{2} will be the most significant term. Hence we have that E_{0} < 8^{n}ϵ^{2}. If we had 10^{−M} digits initially testing for π^{n} we end up with E_{0} < 8^{n}/(10^{2M}) and so if n < 2M we will get E_{0} < 8^{n}/(10^{2M}) < 10^{−M}. So for any power of π, only twice the number of digits of precision desired are necessary.

To determine a formula for ζ(s), or more generally L(s, χ), we must of course calculate L(s, χ) to sufficiently high precision. For Re(s) > 1 we can use the Dirichlet series directly, but for small values of s this series does not converge rapidly. In fact, the first 100 of ζ(3) yields a results correct only to 4 digits. A rough estimate of the error from the series for the ζ function can be estimated as

|ζ(s) − ∑_{n=1}^{N} \frac{1}{n^{s}} | < ∫_{N}^{∞} \frac{dx}{x^{s}} < \frac{1}{N^{s−1}}

and the same estimate clearly applies to L-functions as well. Thus increasing s makes a tangible difference: while we obtain only 4 correct digits of ζ(3) from the first 100 terms of the sum, we can expect 16 correct digits from the first 100 terms of the sum for ζ(9). For small values of s we will attempt to develop some rapidly converging approximate formulas for ζ(s) and L(s, χ) using their functional equations. (This will be spelled out in detail in the final draft.)

For now we will defer the discussion of calculating ζ(3) and simply refer to a table of values for this expression. We use Mathematica to implement S_{−}(3) in its series form and calculate 100 terms. With these values accurate to 10 decimal places, we form the matrix with the identity matrix appended on the right with S_{−}(30), π^{3}, ζ(3) respectively, all multiplied by 10^{10}, namely:

( I_{3} | \begin{matrix}18713728 \\ 310062766803 \\ 120205569032 \end{matrix} )

After LLL, we obtain:

( −360 7 −180 ; ⋮ | −219 ; ⋮ )

where we include only the top row since its rightmost value is clearly the “shortest” of all the rows. This corresponds to the formula:

−180ζ(3) + 7π^{3} − 360S_{−}(3) = 0

or just

ζ(3) = \frac{7}{180}π^{3} − 2S_{−}(3)

Hence, we have some motivation for continuing this process into ζ(2n+1) for arbitrary n.

In executing the LLL algorithm for larger values of s, the larger the initial integers on the right-most side, the closer to 0 the dot product of the resultant factors of each of our “ingredients” with the vector containing the respective “ingredients”, {π^{s}, S_{−}(s), S_{+}(s), ζ(s)} will be. This is because if we have a precise formula, it will send the values of our “ingredients” closer to 0 the more accuracy we achieve, so we can eliminate false-positives linear combinations, which will diverge further from zero with respect to a precise linear combination as we increase accuracy. For our purposes, we will say that a vector in the resultant reduced-basis row vector is “short” if the right most value is at most 1/100 the value of the other row vectors in our matrix. We will then adjust the accuracy of the calculations needed until some vector in our result satisfies this criterion. We know easy to implement and rapidly convergent formulae for all the expressions except ζ(s) so far, and we will derive an exponentially convergent formula for calculating it.

Intermezzo: Calculating ζ(k).

We want to calculate ζ(k) with k > 1. Define the ξ function to be

ξ(s) = π^{−s/2}Γ(s)ζ(s)

and consider the function ξ(s)/(s − k). ξ has simple poles at 0 and 1 with residues -1 and 1, respectively, and the factor of 1/(s−k) adds a simple pole at k with residue ξ(k). Consider the contour illustrated in Figure 1. Since the Γ function decays rapidly as Im(s) → ∞, the integrals along sides II and IV vanish as T → ∞. Let

A = \frac{1}{2πi} ∫_{I} ξ(s)\frac{ds}{s − k}
B = \frac{1}{2πi} ∫_{III} ξ(s)\frac{ds}{s − k} = \frac{1}{2πi} ∫_{I} ξ(s)\frac{ds}{1 − s − k}

where, in the last step we have applied the made the change of variables s → 1 − s and used the functional equation for ξ, namely ξ(s) = ξ(1 − s). Then

A + B = ξ(k) + \frac{1}{1 − k} + \frac{1}{k} .

So we want to evaluate the integrals A and B. Expanding ξ(s) we get

A = \frac{1}{2πi} ∫_{I} ∑_{n=1}^{∞} n^{−s}π^{s/2}Γ(s/2)\frac{ds}{s − k} = ∑_{n=1}^{∞} \frac{1}{2πi} ∫_{I} \frac{π^{s/2}Γ(s/2)}{s − k}n^{−s} ds.

Letting s = 2w we find that

A = ∑_{n=1}^{∞} \frac{1}{2πi} ∫_{(k/2+1)}(πn^{2})^{−w}Γ(w)\frac{2 dw}{2w − k}

where the integral is carried out along the line Re(w) = \frac{k}{2} + 1.

Call each integral in the sum I_{n}. We will evaluate this integral by considering the contour in Figure 2. The integrals along the horizontal lines do not contribute because of the decay of Γ. The integral along the left side is always assumed to intersect the real axis at a half-integer, so since

Γ(1 − s) = \frac{π}{Γ(s) sin(πs)},

this piece does not contribute either. Let y = πn^{2}. The integrand has a simple pole at k/2 with residue Γ(k/2)y^{−k/2} and since Γ has simple poles at the non-positive integers with residues (−1)^{n}/n!,

I_{n} = Γ(k/2)y^{−k/2} − ∑_{n=0}^{∞} \frac{(−y)^{n}}{n!(n + \frac{k}{2})}.

This sum resembles a series expansion of e^{−y} about the origin, and in fact,

∑_{n=0}^{∞} \frac{(−y)^{n}}{n!(n + \frac{k}{2})} = 2y^{−k/2} ∫_{0}^{√\overline{y}} t^{k−1}e^{−t^{2}}dt

and making the change of variables u = t^{2} we have,

2y^{−k/2} ∫_{0}^{y} u^{\frac{k}{2} −1}e^{−u} \frac{du}{2}.

If we define the incomplete gamma function to be

Γ(s, y) = ∫_{y}^{∞} u^{s−1}e^{−u}du

for all y > 0, then Γ(s, y) is entire, and

I_{n} = y^{−k/2}Γ(k/2, y).

From this we conclude that

A = ∑_{n=1}^{∞} I_{n} = π^{−k/2} ∑_{n=1}^{∞} \frac{Γ(k/2, πn^{2})}{n^{k}}.

By a similar argument we find that

B = π^{\frac{k−1}{2}} ∑_{n=1}^{∞} \frac{Γ(\frac{1−k}{2}, πn^{2})}{n^{1−k}}.

Putting this all together we have

ξ(k) = A + B + \frac{1}{k(k − 1)}

which implies that

ζ(k) = \frac{1}{Γ(k/2)} ∑_{n=1}^{∞} \frac{Γ(\frac{k}{2}, πn^{2})}{n^{k}}

+ \frac{π^{k− \frac{1}{2}}}{Γ(k/2)} ∑_{n=1}^{∞} \frac{Γ(\frac{1−k}{2}, πn^{2})}{n^{1−k}}

+ \frac{π^{k/2}}{k(k − 1)Γ(k/2)}.

The sums in this expression for ζ(k) are rapidly convergent, unlike the Dirichlet series that usually defines t he ζ function. We can find a bound for the summand by repeatedly integrating by parts:

Γ (\frac{k}{2} , πn^{2}) = ∫_{πn^{2}}^{∞} u^{\frac{k}{2} −1}e^{−u}du
= e^{−πn^{2}}(πn^{2})^{\frac{k}{2} −1} + (\frac{k}{2} − 1) ∫_{πn^{2}}^{∞} u^{\frac{k}{2} −2}e^{−u}du
⋮
≤ [(\frac{k}{2} − 1) (\frac{k}{2} − 2) … (\frac{1}{2}) + (\frac{k}{2} − 2) … (\frac{1}{2}) + … + (\frac{1}{2})] e^{−πn^{2}}
≤ \frac{k}{2} [(\frac{k}{2} − 1) (\frac{k}{2} − 2) … (\frac{1}{2})] e^{−πn^{2}}
≤ \frac{k}{2} Γ (\frac{k}{2} + 1) e^{−πn^{2}}

This is faster than exponential decay with respect to n, so our sum will converge very rapidly after just a few terms in A and B. The denominator of B grows with respect to n, but clearly the decay e^{−πn^{2}} overwhelms the denominator.

In order to run the LLL algorithm we need at least 100 digits of precision. Suppose we calculate the first M terms of each sum. Clearly the sum in B decays more slowly in n because of the factor of 1/n^{1−k}. Nevertheless we can bound our error since

∑_{n=M+1}^{∞} \frac{Γ(\frac{k}{2}, πn^{2})}{n^{1−k}} ≤ 2 \frac{k}{2} Γ (\frac{k}{2} + 1) ∫_{M}^{∞} e^{πn^{2}}n^{k−1}dn

≤ k Γ (\frac{k}{2} + 1) ∫_{M^{2}}^{∞} e^{−πu}u^{k/2}du
≤ k [Γ (\frac{k}{2} + 1)]^{2} ∫_{M^{2}}^{∞} e^{−πu}du
≤ k [Γ (\frac{k}{2} + 1)]^{2} e^{−πM^{2}}.

For k = 100, putting M = 15 yields over 175 digits of accuracy. In practice, any M ≤ 20 is sufficient as we have that e^{−π} < 10^{−1} so roughly speaking, e^{−πM^{2}} gives about M^{2} digits of accuracy.

We compute Γ(s, y) for small y (say, ≤ 10) with the sum

Γ(s, y) = Γ(s) − ∑_{n=0}^{∞} \frac{(−1)^{n}y^{n+s}}{n!(n + s)}

where the n! gives rapid convergence. For larger n we use a continued fraction expansion:

Γ(s, y) = \frac{e^{−y}y^{s}}{y + 1 − s − _{\frac{1(1−s)}{y+3−s−_{\frac{2(2−s)}{y+5−s− \frac{3(3−s)}{y+7−s−...}}}}}}

The L-functions share similar rapidly convergent formulae, which we make use of [C]. If χ is an even character (mod f) then

Γ (\frac{s}{2}) L(s, χ) = ∑_{n=1}^{∞} \frac{χ(n)Γ(\frac{s}{2}, \frac{πn^{2}}{f})}{n^{s}} + G(χ)π^{s−1/2}f^{−s} ∑_{n=1}^{∞} \frac{\overline{χ(n)}Γ(\frac{1−s}{2}, \frac{πn^{2}}{f})}{n^{1−s}}

and for an odd character χ,

Γ (\frac{s + 1}{2}) L(s, χ) = ∑_{n=1}^{∞} \frac{χ(n)Γ(\frac{s+1}{2}, \frac{πn^{2}}{f})}{n^{s}}+\frac{G(χ)}{i}π^{s−1/2}f^{−s} ∑_{n=1}^{∞} \frac{\overline{χ(n)}Γ(\frac{2−s}{2}, \frac{πn^{2}}{f})}{n^{1−s}}

where

G(χ) = ∑_{n (mod N)} χ(n) exp (\frac{2πin}{N}) .

For each of the following s ≥ 5, 100 digits is more than enough to obtain a clearly shorter row vector from LLL. It is not really necessary to use fewer digits at this point since we have enough computational power to quickly finish our calculations.

ζ(5) = \frac{1}{294}π^{5} − \frac{72}{35}S_{−}(5) − \frac{2}{35}S_{+}(5)
ζ(7) = \frac{19}{56700}π^{7} − 2S_{−}(7)
ζ(9) = \frac{125}{3704778}π^{9} − \frac{992}{495}S_{−}(9) − \frac{2}{495}S_{+}(9)
ζ(11) = \frac{1453}{425675250}π^{11} − 2S_{−}(11)
ζ(13) = \frac{89}{257432175}π^{13} − \frac{16512}{8255}S_{−}(13) − \frac{2}{8255}S_{+}(13)
ζ(15) = \frac{13687}{390769879500}π^{15} − 2S_{−}(15)
ζ(17) = \frac{397549}{112024529867250}π^{17} − \frac{261632}{130815}S_{−}(17) − \frac{2}{130815}S_{+}(17)
ζ(19) = \frac{7708537}{21438612514068750}π^{19} − 2S_{−}(19)
ζ(21) = \frac{68529640373}{188106381576225253125}π^{21} − \frac{4196352}{2098175}S_{−}(21) − \frac{2}{2098175}S_{+}(21)
ζ(23) = \frac{4472029801}{1211517431782539131250}π^{23} − 2S_{−}(23)
ζ(25) = \frac{2598638688071}{6948173623016040171631875}π^{25} − \frac{67100672}{33550335}S_{−}(25) − \frac{2}{33550335}S_{+}(25)
ζ(27) = \frac{149780635937}{3952575621190533915703125}π^{27} − 2S_{−}(27)

Note the obvious pattern if s is of the form 4n + 3. If s is of the form 4n + 1, then the S_{−} and S_{+} terms are of the form

−\frac{2}{d}((d + 1)S_{−}(s) + S_{+}(s)).  (4.1)

In testing these formulae for accuracy, we found that every increase in the precision of machine calculations provided successively smaller error terms. We have good reason then to suspect that these are indeed precise expressions for ζ(2n + 1) for n ≥ 1.

5 Fifth Section: Mahler Measure Results

We tried to compute relations between the Mahler Measures of families of 2-variable polynomials and L-functions. Using our rapidly convergent expression for the L-function, we can quickly evaluate one part of each conjectured Mahler Measure relation, and all that remains is to compute the Mahler Measure of each polynomial involved. This proved to be somewhat difficult using the default methods in Mathematica. The precise expression we hoped to integrate is:

m(P) = ∫_{0}^{1} ∫_{0}^{1} log |P (e^{2πit}, e^{2πiu})| dt du.

and we see immediately that difficulties arise when P (e^{2πit}, e^{2πiu}) = 0 since log |P | blows up at such points. To, get by this problem, we will assume the contributions at these singularities are small due to the fast convergence of ∫_{0}^{1} log(x^{n}) dx for any n ≥ 0. We can then approximate such points by replacing the values around P (e^{2πit}, e^{2πiu}) = 0 by some small ϵ:

m(P) = ∫_{0}^{1} ∫_{0}^{1} log | max[P (e^{2πit}, e^{2πiu}), ϵ]| dt du.

Using this method computationally, we are able to achieve about 5 digits of precision for most polynomials of interest no matter how small we pick ϵ or how finely we partition the domain of integration.

We have verified several relations from Boyd’s paper to this precision. However, we do not have enough knowledge and insight about the Mahler Measure to systematically pick a family of polynomials to test. Moreover, some refinement of technique is required before we can establish a more convincing relationship past 5 digits.

6 Bibliography

References

[B] D. W. Boyd Mahler’s Measure and Special Values of L-Functions, Experimental Mathematics 7:1,37-82 (1998).

[C] H. Cohen. High precision computation of Hardy-Littlewood constants,

http://www.math.u-bordeaux.fr/ cohen/hardylw.dvi. (2003).

[W1] E. Weisstein. LLL Algorithm. Mathworld,

http://mathworld.wolfram.com/LLLAlgorithm.html. (2003).

[W2] E. Weisstein. Riemann Zeta Function. Mathworld,

http://mathworld.wolfram.com/RiemannZetaFunction.html. (2003).

7 Acknowledgements

We would like to thank Florin Spinu for his immense help in the completion of this paper by patiently answering our questions and providing us with notes and various reference materials. We also thank Steven J. Miller for recommending the technique we used in calculating Mahler Measure.