The evaluation of Knopfmacher's curious limit

Daniel Lichtblau
August, 2000

Abstract

We begin with a function expressed as a certain infinite product. It is a  twice-mutated  variation of another product that has its origins in counting  irreducible factors of univariate polynomials over Galois fields.  Knopfmacher's limit is taken as we approach 1 from below in this product. We derive and execute an algorithm that finds  a good approximation to this limit using moderate computational resources. We  also investigate the coefficients of the power series for the logarithm of  this product; these will be shown to exhibit size fluctuations that render  the straightforward power series inadequate for the purpose of estimating the  limit in question.

Statement of the problem

The problem was posed by Arnold Knopfmacher in a note to the Usenet news  group comp.soft-sys.math.mathematica,  January 1999 [Knopfmacher 1999a] (herein we make  small changes in notation). We wish to obtain a numerical estimate (say 8 decimal digits) of the limit as x tends to 1 from below of the function

p[x]=(((∏_(k=2)^∞) (1-((x^(m[k]))/(k+1))))/(1-x))

where m[k]=k-(k/(d[k])) and d[k] is the smallest prime factor of k.

There are similar formulas in [Knopfmacher and Warlimont 1995] analyzing  probabilities related to numbers of irreducible factors of distinct degrees  in univariate polynomials over Galois fields. According to [Knopfmacher 1999b] his coauthor studied a limiting case of such a formula. It gave rise  to a limit seemingly harder than the one above. Knopfmacher extracted the  present problem as a simpler case to try first, and subsequently posted it to  the news group. Thus this example might be described as a double mutation of a problem that is of independent interest elsewhere in the realm of number  theory. It later turned out that the original problem was the more tractable but by then the mutation had acquired a life of its own. In this report we  show how to compute a good approximation to the above limit using a blend of  theory and the computational capabilities of Mathematica [Wolfram 1999] (Mathematica (TM) is a registered trademark of Wolfram Research, Incorporated).

In Mathematica this function may be written as

d[k_] := Divisors[k][[2]]
m[k_]:=k-(k/(d[k]))
p[x_]:=(((∏_(k=2)^∞)(1-((x^(m[k]))/(k+1))))/(1-x))

Crude bounds

First we deduce crude bounds for the limit. Among other things this will  demonstrate the existence of a lim inf and lim sup.

Proposition 1:
(a) A lim inf for p[x] as x->1 is given by (e^EulerGamma).
(b) A lim sup is given by 2 (e^EulerGamma).

Proof: For 0<x<1 note that

(((∏_(k=2)^∞)(1-((x^(⌈(k/2)⌉))/(k+1))))/(1-x))<p[x]<(((∏_(k=2)^∞)(1-((x^(k-1))/(k+1))))/(1-x))

Everything is positive so the inequalities are preserved on taking logarithms:

-Log[1-x]+(∑_(k=2)^∞)Log[1-((x^(⌈(k/2)⌉))/(k+1))]< Log[p[x]]<
	-Log[1-x]+(∑_(k=2)^∞)Log[1-((x^(k-1))/(k+1))]

We tackle the second inequality in order to prove (b).

Log[p[x]]<-Log[1-x]+(∑_(k=2)^∞)Log[1-((x^(k-1))/(k+1))]=(∑_(k=1)^∞)((x^k)/k)-(∑_(k=1)^∞)(∑_(j=1)^∞)((((((x^k)/(k+2))))^j)/j)=(∑_(k=1)^∞)((2 (x^k))/((k^2)+2 k))-(∑_(k=1)^∞)(∑_(j=2)^∞)((x^(j k))/(j (((k+2))^j)))

The first sum in the limit as x->1 is readily computed to be 3/2.

We evaluate this second sum at x=1 (justification of the interchange of sum and limit will be deferred to  the next section). We first reverse the order of summation.

innersum = Expand[Sum[1/(j*k^j), {k,3,Infinity}]]

We now sum this.

fullsum=Sum[innersum, {j,2,Infinity}]
approx=N[fullsum]

logupbnd =3/2-fullsum
upbnd = N[Exp[logupbnd]]

So we see that an upper bound for p[x] as x->1 is about 3.56.

We now prove (a). To get a lower bound we now minorize the exponents.

Log[p[x]]>-Log[1-x]+(∑_(k=2)^∞)Log[1-((x^(⌈(k/2)⌉))/(k+1))] =(∑_(k=1)^∞)((x^k)/k)-(∑_(k=2)^∞)(∑_(j=1)^∞)((((((x^(⌈(k/2)⌉))/(k+2))))^j)/j)=

((2 x)/3)+ (∑_(k=2)^∞)((1/k)-(1/(2 k))-(1/(2 k+1))) (x^k) -(∑_(k=2)^∞)(∑_(j=2)^∞)((x^(j ⌈(k/2)⌉))/(j (((k+1))^j))) = 

((2 x)/3)+ (∑_(k=2)^∞)((x^k)/(4 (k^2)+2 k))-(∑_(k=2)^∞)(∑_(j=2)^∞)((x^(j ⌈(k/2)⌉))/(j (((k+1))^j)))

We will evaluate this at x=1 to obtain EulerGamma. Note that the single summation can readily be done  by hand using basic calculus methods (one obtains 5/6-Log[2]; we omit the details). So the limiting value of the first two terms is 3/2-Log[2]. As it is important for the actual limit estimate later we indicate how  to do the double summation without recourse to software; we will also show  how it may be evaluated using Mathematica. We evaluate at x=1 so we want to find

Limit[(∑_(k=3)^M)(∑_(j=2)^∞)(1/(j (k^j))),M->∞]

Limit[(∑_(k=3)^M)(∑_(j=2)^∞)(1/(j (k^j))),M->∞]

Now rewrite this as

Limit[(∑_(k=3)^M)(∑_(j=1)^∞)((1/(j (k^j)))-(1/k)),M->∞]

Evaluating the inner sum gives

Limit[(∑_(k=3)^M)(-Log[1-(1/k)]-(1/k)),M->∞]

which in turn is

Limit[((Log[3]-Log[2])+...+(Log[M]-Log[M-1])-(∑_(k=3)^M)(1/k)),M->∞]

Telescoping and simplification gives

Limit[Log[M]-(∑_(k=1)^M)(1/k),M->∞]-Log[2]+(3/2)=-EulerGamma-Log[2]+3/2

Negating this and adding the result for the other terms gives a lower  bound of EulerGamma for our sum. We exponentiate to obtain a lower bound for p[1] that is near 1.78, or half of our upper bound.

N[(E^EulerGamma)]

Approximating the actual limit

Here I will discuss a nice analysis and resulting estimate worked out by Jurgen Tischer of Universidad del Valle, Columbia [Tischer 1999]. First I  will mention a method I had previously tried that is not terribly promising.  In brief we truncate the series for the logarithm, evaluate using exact or  high precision arithmetic at x=1, and that gives an estimate for the limit. Using this approach it is not  difficult to approximate the result to four places without excessive  computational effort. One learns that p[1]≈2.292. There are some problems with this method. First, it appears to require  significant work to get even one more decimal place of precision. Moreover  for truncation after N terms it is not clear how we might get a bound on the error. For odd k we get terms that are roughly (1/(2 k)) in the series and this alone makes an error bound for this approximation  appear to be less than promising. We will return to these matters in the next  section.

To proceed with the approximation we will look at the logarithm in more detail. It may be written as

Log[p[x]]=-Log[1-x]+(∑_(2|k))Log[1-((x^(k/2))/(k+1))]+(∑_(2∤k,3|k))Log[1-((x^(2 k/3))/(k+1))]+...=-(∑_(k=1)^∞)((x^k)/k)+(∑_(2|k))((x^(k/2))/(k+1))+(∑_(2|k))(∑_(j=2)^∞)((x^(j k/2))/(j (((k+1))^j)))+(∑_(2∤k,3|k))((x^(2 k/3))/(k+1))+(∑_(2∤k,3|k))(∑_(j=2)^∞)((x^((2/3)j k))/(j (((k+1))^j))) +...

If we truncate the power series at degree N, say, it can be shown that error from dropping terms from the double sums  (where j≥2) is O((1/N)). But it turns out that we need not rely on this error bound; Tischer derived a rather nice exact formula for these terms, summed to  infinity. After showing his derivation we will then approximate the sum of  the "main" terms, which are the single summations in the formula above.

We begin by rewriting Log[p[x]] (we will call this lp[x] from now on) as

lp[x]=x+((x^2)/2)+(∑_(k=2)^∞)(((x^(k+1))/(k+1))-((x^(m[k]))/(k+1))) -(∑_(k=2)^∞)(∑_(j=2)^∞)((1/j))(((((x^(m[k]))/(k+1))))^j)

Proposition 2:  lp[x]->EulerGamma+Log[2]+(∑_(k=2)^∞)(((x^(k+1))/(k+1))-((x^(m[k]))/(k+1))) as x->1 from below.

Proof: The limiting value of first two terms in (1) is obvious so we focus on the last sum. We showed above how to evaluate  it by hand at x=1, now we will do so using Mathematica. It is not hard to see that for x<1 everything in sight converges nicely and so we can switch summation order  to obtain

(∑_(j=2)^∞)((1/j))(∑_(k=2)^∞)(((((x^(m[k]))/(k+1))))^j)

Fix j and look at the inner sum. Since (k/2)≤m[k]≤k and we approach from x<1 we have

(∑_(k=2)^∞)(((((x^k)/(k+1))))^j)≤(∑_(k=2)^∞)(((((x^(m[k]))/(k+1))))^j)≤(∑_(k=2)^∞)(((((x^(k/2))/(k+1))))^j)

Simple limit arguments and Abel's theorem [Rudin 1976 p. 174] demonstrate  that as x->1 the middle sum is squeezed to

(∑_(k=2)^∞)((((1/(k+1))))^j)

We must justify the interchange of limit with the outer sum. This may be done using the Weierstrass M-test [Rudin 1976 p. 148]. First note that for all x∈[0,1] we have

((1/j))(∑_(k=2)^∞)(((((x^(m[k]))/(k+1))))^j)≤((1/j))(∑_(k=2)^∞)((((1/(k+1))))^j)

and all the summands are nonnegative. Therefore if we know that the right hand side summed over j converges then we  have also justified the interchange of limit and outer sum. Alternatively we might overestimate the tail of the outer sum with an  integral and show that it goes to zero as the lower bound increases.

Our sum at  x=1 is

(∑_(j=2)^∞)((1/j))(∑_(k=2)^∞)((((1/(k+1))))^j)=(3/2)-EulerGamma-Log[2]

as was demonstrated in the previous section. This finishes the proof of proposition 2.

Now we must work with that remaining sum. As all summands are negative we may reorder them; we have

(∑_(k=2)^∞)(((x^(k+1))/(k+1))-((x^(m[k]))/(k+1)))=(∑_(d[k]=2)^∞)(((x^(k+1))/(k+1))-((x^(((k/2))))/(k+1))) +(∑_(d[k]=3)^∞)(((x^(k+1))/(k+1))-((x^((((2 k)/3))))/(k+1)))+(∑_(d[k]=5)^∞)(((x^(k+1))/(k+1))-((x^((((4 k)/5))))/(k+1))) + ...

We are in effect sieving our sum over smallest prime divisors. It will turn out that we can evaluate these subsums and eventually truncate in a way that gives a very accurate result. To start out, we want to evaluate

(∑_(d[k]=Prime[j])^∞)(((x^(k+1))/(k+1))-((x^(((((Prime[j]-1) k)/(Prime[j])))))/(k+1)))

We will need some simple functions below.

q[j_]:=(∏_(k=1)^j)Prime[k]
r[j_]:=(∏_(k=1)^(j-1))(Prime[k]-1)
frac[j_]:=((r[j])/(q[j]))

The terms k for which Prime[j] is the smallest divisor larger than 1 fall into finitely many congruence sets.

For example, when the prime in question is 5, the applicable values for k are {5,25,35,55,65,85,95,...}. This may be partitioned as {5,35,65,95,...} and {25,55,85,...}. In each case, the step size is 30.

Proposition 3: For Prime[j] the step size in each congruence class will be

q[j]=(∏_(m=1)^j)Prime[m]

Proof: This is perhaps well known but we provide a proof for sake of completeness. Clearly if Prime[j] is the smallest divisor of t then it is also the smallest divisor of t+q[j], t+2 q[j], .... Suppose there is some smaller (Overscript[q, ∼])[j] with this property. Then for some k<j we have Prime[k]∤(Overscript[q, ∼])[j]. Hence Prime[k] is relatively prime to both t and (Overscript[q, ∼])[j].

Suppose y and z are distinct members of the set

{t,t+(Overscript[q, ∼])[j],t+2 (Overscript[q, ∼])[j],...,t+(Prime[k]-1) (Overscript[q, ∼])[j]}

with z>y. Then

z-y=((Subscript[k, 1])-(Subscript[k, 2])) (Overscript[q, ∼])[j]

where

(Subscript[k, 1])-(Subscript[k, 2])<Prime[k]

Hence

Prime[k]∤z-y

so y and z occupy distinct congruence classes modulo Prime[k]. Thus a pigeonhole argument shows that Prime[k] must divide one of these, contradicting the assumption that Prime[j] is the smallest prime divisor.

We will also need to know the number of congruence classes for given Prime[j].

Proposition 4: The number of congruence classes for Prime[j] is

r[j]=(∏_(k=1)^(j-1))(Prime[k]-1)

Proof: This too is known but will likewise be proven for sake of completeness.  First a check: For j=1, Prime[j] is 2, the number of congruence classes is 1. For j=2, Prime[j] is 3, the number of congruence classes is again 1. For j=3, Prime[j] is 5 and there are two congruence classes. For j=4, Prime[j] is 7 and there are 8 congruence classes. The general assertion is that in the set 1≤n≤q[j] there are r[j] elements whose smallest prime divisor is Prime[j]. We checked this explicitly above for the first few primes. A general proof will be obtained working directly with frac[j]=r[j]/q[j]. It suffices to show that

frac[j]=(((Prime[j-1]-1) )/(Prime[j])) frac[j-1]

This in essence defines the sieving process. At step j the fraction of integers with smallest prime divisor equal to Prime[j] will be (1/(Prime[j])) times the fraction not already removed for having smaller divisors. In  other words, we have

frac[j]=((1-(∑_(k=1)^(j-1))frac[k])/(Prime[j]))

Thus

frac[j]=(((1-(∑_(k=1)^(j-2))frac[k])-frac[j-1])/(Prime[j]))=(((((1-(∑_(k=1)^(j-2))frac[k]) Prime[j-1])/(Prime[j-1]))-frac[j-1])/(Prime[j]))=((Prime[j-1] frac[j-1]-frac[j-1])/(Prime[j]))

which completes the proof.

Now that we have determined the number of congruence classes for each  prime we will need to show that the relevant limits for each prime are independent of these classes. So we must find the limits of these sums for all congruence class members as x->1. First we form such a sub-subsum as below.

Sum[((x^(k+1))/(k+1))-((x^(((Prime[j]-1) k)/(Prime[j])))/(k+1)),{k,init,∞,step}]

Next we find the limit of this sum as x->1.

s1 = Sum[((x^(k+1))/(k+1))-((x^(((Prime[j]-1) k)/(Prime[j])))/(k+1)),{k,init,∞,step}]

While this may look awkward we see that it has a tractable limit.

lims1 = Limit[s1, x->1]

Actually this can be simplified considerably.

lims1 = FullSimplify[lims1]

We next use an improper high-school simplification because we know it holds on our domain:

lims1 = lims1 /. Log[r_]-Log[s_]->Log[r/s]

In addition to simplicity this limit has the virtue of being independent  of congruence class. That is, it is does not depend on the parameter init. This significantly simplifies the problem of evaluating the full sums.

In terms of the function frac, we have shown

Proposition 5:

Limit[(∑_(d[k]=Prime[j])^∞)(((x^(k+1))/(k+1))-((x^(((((Prime[j]-1) k)/(Prime[j])))))/(k+1))),x->1] = frac[j]Log[1-(1/(Prime[j]))]

We now return to the full sum in (2). Pulling the limit inside the first M terms of the right hand side gives

Limit[(∑_(k=2)^∞)(((x^(k+1))/(k+1))-((x^(m[k]))/(k+1))),x->1]=(∑_(j=1)^M)frac[j] Log[1-(1/(Prime[j]))]+Limit[(∑_(d[k]≥Prime[M+1])^∞)(((x^(k+1))/(k+1))-((x^(((((d[k]-1) k)/(d[k])))))/(k+1))),x->1]

We now want to show the tail gets arbitrarily small as M->∞ in order to interchange limit with infinite sum, as this will yield

Theorem 1:

Limit[(∑_(k=2)^∞)(((x^(k+1))/(k+1))-((x^(m[k]))/(k+1))),x->1]=(∑_(j=1)^∞)frac[j] Log[1-(1/(Prime[j]))]

Proof: The tail is bounded in absolute value by

Limit[(∑_(d[k]≥Prime[M+1])^∞)(((x^(((((d[k]-1) k)/(d[k])))))/(k+1))-((x^(k+1))/(k+1))),x->1]

and we have

(∑_(d[k]≥Prime[M+1])^∞)(((x^(((((d[k]-1) k)/(d[k])))))/(k+1))-((x^(k+1))/(k+1)))<(∑_(d[k]≥Prime[M+1])^∞)(((x^(((((Prime[M]-1) k)/(Prime[M])))))/(k+1))-((x^(k+1))/(k+1)))<(∑_(k=0)^∞)(((x^(((((Prime[M]-1) k)/(Prime[M])))))/(k+1))-((x^(k+1))/(k+1)))

The sum on the right is a special case of (3) above. It has a fairly simple closed form.

tailsum= Sum[((x^(k+1))/(k+1))-((x^(((Prime[M]-1) k)/(Prime[M])))/(k+1)),{k,0,∞}]

We compute the limit as x->1 from below.

Limit[tailsum, x->1, Direction->1]

As this goes to 0 for M→∞ we have shown that the error term in (4) can be made arbitrarity small. This suffices to prove theorem 1.

Now we need to estimate our sum. We have

Corollary: p[x]->2.2921736953 as x->1 from below.

Proof: We will do a careful computation of our estimate and bounded error.  Define

estimate[n_]:=(∑_(j=1)^n)frac[j] Log[1-(1/(Prime[j]))]
error[n_]:=(∑_(j=n+1)^∞)frac[j] Log[1-(1/(Prime[j]))]

We will bound the magnitude of the error term and explicitly compute the  estimate for some large n. To bound this error, first observe that the products frac[j] are conveniently bounded by (constant/(Prime[j])); this follows from a simple telescoping argument since successive factors  are of the form ((Prime[k-1]-1)/(Prime[k])). We will get a reasonable value for constant as part of the computation that estimates the value. For this we require  a faster numeric version of the algorithm for frac[j] and estimate[j]. This will allow us to avoid the incessant recomputation implicit in the  definition of estimate[n]; moreover in using numeric approximation we avoid generating large  rationals. This will introduce some error but as we use Mathematica's significance arithmetic [Keiper 1992] at sufficiently high precision it  will be small compared to the bounds that we also compute. We will overestimate those bounds so arithmetic correct to 20 digits will more than suffice for our purposes.

numestimate[n_,prec_]:= Module[{sum=-1/2*Log[2],frac=1/2,p1=Prime[1],p2},
Do[
p2=Prime[j];
frac *= (((N[p1-1,prec])/p2));
sum += frac*Log[N[1-(1/p2),prec]];
p1=p2,
{j,2,n}];
{frac*p2,sum}]

Timing[{mfact,est} = numestimate[2*10^7+1,25]]

Setting n to be 2* (10^7)+1, an upper bound on the magnitude of error[n] is given by

mfact(∑_(j=n+1)^∞)(1/(Prime[j]))|Log[1-(1/(Prime[j]))]|

This in turn is bounded by

mfact(1+(1/(Prime[n])))(∑_(j=n+1)^∞)(1/((j^2) ((Log[j])^2)))

Here we have overestimated the error by using the asymptotic term j Log[j] from the prime number theorem to underestimate Prime[j] [Ribenboim 1996 p. 249], and we use the basic calculus result that (1+(1/(Prime[n]))) overestimates the ratio of -Log[1-(1/(Prime[j]))] to (1/(Prime[j])) for j>n. We will now ignore this factor because we more than compensate when we  round upward the value of mfact in the fourth decimal place. Thus we have as a bound:

(.0285)(∑_(j=n+1)^∞)(1/((j^2) ((Log[j])^2)))

Using simple calculus we overestimate this sum by an integral:

errormax =(57/2000)(∫_(2 (10^8))^∞)(1/(((j Log[j]))^2))\[DifferentialD]j

N[errormax]

Now we add estimate[2 (10^7)+1] to the exact terms given in proposition 2. Then we will exponentiate the sum to get an approximation of p[x] as x->1 from below. We also multiply this by (e^errormax)-1 to get a bound on the error in our estimate of the limit.

p1 = Exp[EulerGamma+Log[2]+ est]

errbound = N[Exp[errormax]-1]*p1

So we have obtained a result that is correct to about the eleventh decimal place, demonstrating the corollary. In [Tischer 1999] the computations are done with around (10^8) terms. This can achieve a result with approximately one more decimal place of accuracy. In actual fact it was presented as being slightly less accurate because the error bound was looser.

Further remarks about the truncated power series

We will use the Mathematica Series function to work with the truncated power series. First we note a curiousity in the sign pattern of the coefficients.

logp[x_,n_] := Normal[Series[-Log[1-x]+(∑_(k=2)^(2*n))Log[1-((x^(m[k]))/(k+1))], {x,0,n}]]

terms100 = logp[x,100];
coeffs100 = CoefficientList[terms100,x];
signs100 =Map[Sign,coeffs100];
Map[{#[[1]],Length[#]}&,Split[signs100]]

We see an unusual phenomenon wherein signs by-and-large  alternate except for sporadic runs of three consecutive positive coefficients. This holds at least for the first 2000 coefficients. It is an open question as to whether or why this persists.

In attempting to get a numerical approximation from this not-quite-alternating series, at first glance it appeared that one could sum  n terms and hope to obtain O(n) convergence. Further work makes it clear that this cannot be the case, at least not for all values of n. Convergence is slow (though it does in fact seem to converge) for good reason; some of the terms are asymptotically larger than O((1/n)). As we will see, one can bound from above the magnitudes of the coefficients that arise in the power series by O(((log [log[ n]])/n)) for n bigger than 3. Further work will reveal that this bound is the best we can do, and indeed the coefficients that approach it will be not too scarce.

In order to obtain results more readily we will now work with only the  main terms. Specifically, we discard higher-order terms from all but the main logarithmic term (the one that is singular at x=1). The series in question is now

logpmainterms[x_] := (∑_(k=1)^∞)((x^k)/k) -(∑_(k=2)^∞)((x^(m[k]))/(k+1))

where the first sum is for -Log[1-x] and the second is for the main terms from the remaining logarithms.

Thus the remarks below actually apply to a slightly altered problem, but one that is not different in significant aspects from the original (and in  particular, with a bit of work one can show that the derivations of coefficient size bounds will still apply). We write the series for logpmainterms in the symbolic form (∑_(k=1)^∞)(c[k](x^k)), thus implicitly defining a "coefficient function" c[k]. Let us see explicitly how one would describe a general coefficient.  Contributions to the term in (x^k) arise from the (1/k) coefficient in the first sum, as well as from coefficients of the form (1/(j+1)) for all j for which m[j]=k.

Proposition 6: For odd k we have

c[k]=(1/k)-(1/(2 k+1))=((k+1)/(2 (k^2)+k))

Proof: This follows from the following simple observations.

(i) d[2 k]=2, hence m[2 k]=k.

(ii) For j even and j≠2 k, m[j]=(j/2)≠k.

(iii) For j odd, (j/(d[j])) is odd, hence m[j] is even. So m[j]≠k.

We conclude that for odd k there are exactly two contributions to c[k] and they are as in the statement of the proposition.

Corollary: Proposition 6 shows (as previously noted) that the coefficients are in general no smaller than O((1/k)). But they can in fact be bigger. The next result gives an indication of  how large they can be.

Theorem 2:
(a) An upper bound on size of c[n] is (((Subscript[Log, 2])[(Subscript[Log, 2])[n+1]]+1)/n).
(b) There are coefficients that approach bound to within a factor of 3/4. We construct one such explicitly. From the method one will see how to approach the upper bound within any factor less than 1.

Proof: We now investigate coefficients of those terms of even degree. Let us suppose that m[j]=k for some j. Let d[j]=p. Then we know that j-(j/p)=(((p-1) j)/p)=k and thus p k=(p-1) j. Immediately we see that

(i) p-1∣k

(ii) p∣j

(this is in any case obvious because p is given as the smallest prime divisor of j)

(iii) (k/(p-1))=(j/p)

and hence (k/(p-1)) has no prime divisor smaller than p (again because p is the smallest prime divisor of j).

Now we suppose we have all possible j for which m[j]=k. Specifically, suppose

k=(((Subscript[p, 1])-1)/(Subscript[j, 1]))=(((Subscript[p, 2])-1)/(Subscript[j, 2]))=...=(((Subscript[p, s])-1)/(Subscript[j, s]))

where d[(Subscript[j, r])]=(Subscript[p, r]) for all 1≤r≤s and we moreover have ordered these so that

2=(Subscript[p, 1])<(Subscript[p, 2])<...<(Subscript[p, s])≤k+1

and

2 k=(Subscript[j, 1])>(Subscript[j, 2])>...>(Subscript[j, s])≥k+1

As we saw,  (Subscript[p, r])-1∣k for each 1≤r≤s .

Claim 1: (Subscript[p, r])-1 contains all divisors of k that are smaller than (Subscript[p, r]) .

Proof: If not then the equality (Subscript[p, r])k=((Subscript[p, r])-1)(Subscript[j, r]) implies that (Subscript[j, r]) has a divisor (and hence a prime divisor) strictly smaller than (Subscript[p, r]) . This would contradict the fact that (Subscript[p, r])=d[(Subscript[j, r])].

Claim 2: For each 1≤r≤s-1 we have (Subscript[p, r])-1∣(Subscript[p, r+1])-1 .

Proof: (Subscript[p, r])-1 contains all divisors of k that are smaller than (Subscript[p, r]). It cannot contain any other factors; there is no room left over as it is  itself smaller than (Subscript[p, r]). A similar statement holds for the divisors of (Subscript[p, r+1])-1 with respect to (Subscript[p, r+1]). Thus the divisors of (Subscript[p, r])-1 are a proper subset of those of (Subscript[p, r+1])-1.

We obtain the chain (Subscript[p, 2])-1∣(Subscript[p, 3])-1∣...∣(Subscript[p, s])-1 of length s-1 (we do not bother with (Subscript[p, 1]) because it is 2, so (Subscript[p, 1])-1 is 1), and moreover we have inequalities concerning the quotients.  Specifically,  for 1≤r≤s-1 we may write

(Subscript[p, r+1])-1=(Subscript[q, r]) ((Subscript[p, r])-1)

where (Subscript[q, r]) ≥ (Subscript[p, r]). Hence

(Subscript[p, r+1]) ≥ (Subscript[p, r]) ((Subscript[p, r])-1)+1

Now using the fact that (Subscript[p, r])-1≥(2/3)(Subscript[p, r]) we obtain

k+1≥(Subscript[p, s])≥(((((2/3)))^(2^(s-2))))(Subscript[p, 2]^(2^(s-2)))≥(2^(2^(s-2)))

We now take logarithms to base 2 to see (Subscript[Log, 2])[k+1]>(2^(s-2)). Again taking logarithms to base 2, we get

(Subscript[Log, 2])[(Subscript[Log, 2])[k+1]]+2≥s

Thus

|c[k]|=|(1/k)-(∑_(m[j]=k))(1/(j+1))|<(((Subscript[Log, 2])[(Subscript[Log, 2])[k+1]]+1)/k)

where we use the fact that the denominators in the sum all lie between k+1 and 2 k. This completes the proof of part (a).

Now that we have this upper bound for coefficients in the power series, it is natural to ask whether one can derive an asymptotically smaller bound. As stated in part (b), this is not possible. We outline why this is so. First we investigate numerically: we will search for coefficients that approach this bound in ratio. Within the first few thousand terms we can find several that come close.

approxlogp[x_,n_] := Normal[Series[(∑_(k=1)^n)((x^k)/k)-(∑_(k=2)^(2*n))((x^(m[k]))/(k+1)), {x,0,n}]]

Timing[terms4k = approxlogp[x,2^12];]
coeffs4k = CoefficientList[terms4k,x];

Now we look at the sizes of the coefficients. While as we know none can be larger than ((Log[2,Log[2,n]]+1)/n) (for n>3, that is), it appears that some come near to this. For example, in the first 4096 coefficients one is about 0.66 of this maximum and it is roughly three times larger than (1/n).

normedcoeffs = Drop[Abs[coeffs4k],4] * Table[n/(Log[2,Log[2,N[n]]]+1),{n,5,Length[coeffs4k]}];

Max[normedcoeffs]

Let us investigate more closely that large normalized coefficient.

Position[normedcoeffs,Max[normedcoeffs]]

We dropped four coefficients (corresponding to terms of degree zero  through three), so we want to look at the 3067th coefficient in the list, which goes with the term of degree 3066.

(Log[2,Log[2,3067]]+1)/3067// N

Abs[coeffs4k[[3067]]]//N

Abs[coeffs4k[[3067]]*3067]//N

It is claimed above that the best bound we can attain on these coefficients is in fact O(((log [log[ n]])/n)). To see how to approach this bound, it is instuctive to see how we obtain the coefficient for (x^3066).

FactorInteger[3066]

We have a chain of primes

{1+1,2+1,2*3+1,2*3*7+1,2*3*7*73+1}

wherein all factors in the indicated products are themselves prime, and each last factor is at least as large as the previous element in the list (note that 73 is the smallest prime greater than or equal to 43 that works). In general we will look for values of k that have an increasing chain of prime factors (Subscript[q, r]) such that

(Subscript[p, r+1])-1=(Subscript[q, r]) ((Subscript[p, r])-1)

for some ascending set of primes (Subscript[p, r]). For example, to find a term after the 3066th that makes the ratio larger, we look for the smallest prime q larger than 3066 such that 3066 q+1 is also prime. This happens for q=3137, giving k=9618042. The ratio

((Log[2,Log[2,n]]+1)/n)

to the size of the coefficient of (x^9618042) turns out to be near 0.72. Carrying this construction yet another step gives the next q as 9618361, so the value for k is 9618361*9618042 or 92509800069162. The ratio now is about 0.765. This finishes the proof of part (b) and moreover demonstrates how one  might search for coefficients closer in ratio to the upper bound of part (a).

Note that one can get more terms near the upper bound simply by taking for our next prime factor q any prime slightly larger than the smallest one that works, as some of  these may also work. For example, if we use 79 rather than 73 as a factor above (noting that 2*3*7*79+1=3319 is prime), we get a coefficient with ratio just under .66. One thus sees a general means by which to obtain coefficients not far  in ratio below the upper bound given above.

As a final remark, let us recall that a straightforward computation with the truncated power series, taking x=1, gave a result correct to about 4 digits of precision. We know they are correct because they agree with the  result of Tischer's method, which we prove above to give a correct result to  several more places. Also at the beginning of this report we found a lim inf and lim sup for the power series sum as x->1. That notwithstanding, we do not have a stand-alone proof that the power series itself converges to the limit. An indirect proof is that Tischer's method converges, and one observes that the truncated power series only omits from Tischer's truncated sum terms of the same sign, hence it too must converge. Of course the construction of terms with coefficients of size  O(((log [log[ n]])/n)) forces one to conclude that estimation of the limit by truncation of the power series will be problematic, but nonetheless an independent proof of  convergence would be of interest in its own right.

Acknowledgements

I thank Harold Diamond for several helpful discussions about this problem, and for inviting me to speak about it in the weekly departmental seminar in analytical number theory. I also thank the University of Illinois Math Department for extending its courtesy to me while I was a visitor during the academic year 1998-9.

References

J. Keiper (1992). Numerical computation with Mathematica (tutorial notes). An electronic version may be found at:
 http://mathsource.com/Content22/Enhancements/Numerical/0204-028

A. Knopfmacher (1999a). Usenet news group post to comp.soft-sys.math.mathematica.
http://library.wolfram.com/mathgroup/archive/1999/Jan/msg00023.html

A. Knopfmacher (1999b). Private communication.

A. Knopfmacher and R. Warlimont (1995). Distinct degree factorizations for  polynomials over a finite field. Transactions of the American Mathematical  Society 347(6):2235-2243.

P. Ribenboim (1996). The New Book of Prime Number Records. Springer.

W. Rudin (1976). Principles of mathematical Analysis (third edition). McGraw-Hill.

J. Tischer (1999). Usenet news group post to comp.soft-sys.math.mathematica.
http://library.wolfram.com/mathgroup/archive/1999/Jan/msg00355.html

S. Wolfram (1999). The Mathematica Book (4th edition). Wolfram Media/Cambridge University Press.