cp's OEIS Frontend

This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.

A337748 T(m,n) is the least k such that the partial sum of the series Sum_{j=0..k} 1/(m*j+1) is > n, read by ascending antidiagonals.

Original entry on oeis.org

1, 1, 3, 1, 7, 10, 1, 17, 56, 30, 1, 43, 353, 418, 82, 1, 111, 2374, 7100, 3091, 226, 1, 289, 16497, 129644, 142624, 22845, 615, 1, 762, 116759, 2448436, 7078315, 2864677, 168803, 1673, 1, 2021, 835695, 47104189, 363380155, 386462913, 57538579, 1247297, 4549, 1, 5389, 6025478, 916451548, 19003186159, 53930396770, 21100160132, 1155693253, 9216353, 12366
Offset: 1

Views

Author

Gerhard Kirchner, Sep 18 2020

Keywords

Comments

H(m,k) = Sum_{j=0..k} 1/(m*j+1) may be thought of as a generalized harmonic sequence: H(1,k) = H(k+1), where H(k) = Sum_{j=1..k} 1/j is the harmonic and H(2,k) the odd harmonic sequence. As a consequence, T(1,n) = A002387(n)-1 and T(2,n) = A092315(n). For m > 2, the sequence T(m,n) is not currently in the OEIS. Therefore, I give a detailed derivation of the formulas below, see link "Generalized harmonic series". The parameter c(m) in the formula for T(m,n): c(1) = gamma (Euler's constant), c(2) = gamma+1-log(2). For m > 2, c(m) has to be determined individually. The algorithm is described in the formula section though it is not only one formula.
In the Maxima code, two tests are implemented, see Formula. The first test will probably be passed without exception, the second one limits the terms. For greater T(m,n), the precision must be increased.

Examples

			T(m,n) begins:
m=1: 1   3    10      30        82 ...
m=2: 1   7    56     418      3091 ...
m=3: 1  17   353    7100    142624 ...
m=4: 1  43  2374  129644   7078315 ...
m=5: 1 111 16497 2448436 363380155 ...
     .............................
Evaluation of c(3):
1. S(3,99) = 0.1472791427382, see formulas (P1) and (P2).
2. F(3,100) = 10^(-2)/9 + 10^(-4)/27 + 10^(-6)/243 - 10^(-8)/243 ..., see (P5) = 0.001114818889, see formulas (P3) and (P4).
The next summand is 10^(-10)/729 ~ 10^(-13), the precision is approximately  p = 13 digits.
3. beta(3) = 0.1483939616270.
4. c(3) = 3.1320337800204.
Evaluation of T(3,5):
5. Summing up directly with k = 142624: H(3,k-1) = 4.9999999, H(3,k) = 5.0000022. => T(3,5) = 142624.
6. Using the formula with k1 = exp(3*5-c(3))-5/6 = 142623.0446139:
   T(3,5) = floor(k1+1) = 142624.
7,1. floor(k1) = floor(k2) is satisfied with k2 = k1-1/(24*k1) = 142623.0446136.
7,2. With d = 0.0446 and p = 13, k1 < d*10^(-13) is satisfied, too.
		

Crossrefs

Programs

  • Maxima
    block(p:100, km: 100, eps:10^(-p), lim:10^(p), cc:[],
    /*"harm-bfile.txt" with terms < lim is saved*/
    gam: bfloat(%gamma),
    fpprec:p, load(newton1),
    fl: openw("harm-bfile.txt"),
    c(m):= block(x:0, delta: 1, r:0,
    if m=1 then return(gam) else
      (for k from 1 thru km-1 do x: x+bfloat(1/(m*k*(m*k+1))),
       while delta=0 or delta >eps do
          (r: r+1, su:0,
          for j from 1 thru r-1 do su: su+a[r-j]*(-1)^j*binomial(r,j+1),
          ar: bfloat(-((-1)^r/m^(r+1)+su)/r),  ad: ar/km^r, a: append(a,[ar]),
          x: x+bfloat(ad), delta: abs(ad)), return(gam+m*(1-x)))),
    mn:0, ok: true,  nr:0,
    while ok do (mn: mn+1, a:[], cc: append(cc, [c(mn)]), m:mn+1,
       while m>1 and ok do  (m: m-1, n: mn+1-m,
       k1: exp(m*n-cc[m])-(m+2)/(2*m), k2: k1-1/(24*k1),
       ka: floor(k1+1), kb: floor(k2+1),
       if ka>1 then
         (d: kb-k2, if d>1-d then d: 1-d, if ka#kb or kb> d*lim then  ok: false),
      if ok then (nr: nr+1,  printf(fl, concat(nr, " ", ka)), newline(fl) ) ) ),
    close(fl));

Formula

T(m,n) = floor(k1+1) with k1 = exp(m*n-c(m))-(m+2)/(2*m).
For avoiding bad terms, k1 has to pass two tests:
Test 1: floor(k1) = floor(k2) with k2 = k1-1/(24*k1).
Test 2: k1 < d*10^p, where p is the precision of c(m) (number of digits) and d the distance between k1 and the next integer floor(k1) or floor(k1)+1.
The parameter c(m):
c(m) = gamma+m*(1-beta(m)) with (P1) beta(m) = lim_{k->oo} S(m,k) and S(m,k) = Sum_{j=1..k}(1/(m*j)-1/(m*j+1)).
For faster convergence, the infinite sum must be split up: (P2) beta(m) = S(m,k-1) + F(m,k).
F(m,k) can be written as:
(P3) F(m,k) = Sum_{r=1..s-1} b(r)/k^r + R(s) with R(s) ~ b(s)/k^s.
Example: For the precision p = 100 digits, k = 100 is large enough to find s with b(s)/k^s < 10^(-p) where s >= 73.
The coefficients b(r) are defined by the recurrence b(1) = 1/m^2 and, for r>1: (P4) b(r) = ((-1/m)^(r+1)-Sum_{j=1..r-1} (-1)^j*b(r-j)*binomial(r,j+1))/r.
1 m-1 m^2-3*m+2 m^2-2*m+1
(P5) b(1) = -----, b(2) = -----, b(3) = -----------, b(4) = -----------,
m^2 2*m^3 6*m^4 4*m^5