login
A366189
a(n) is the positive integer k such that the k-th derivative of the n-th Bernoulli polynomial B(n, x) contains only integer coefficients but no lower derivative of B(n, x) has this property.
4
1, 1, 2, 1, 2, 1, 2, 3, 2, 1, 2, 1, 2, 3, 2, 2, 3, 3, 4, 3, 2, 3, 4, 5, 2, 3, 4, 1, 2, 1, 2, 3, 4, 5, 3, 1, 2, 3, 4, 5, 6, 3, 4, 5, 4, 5, 6, 7, 5, 3, 4, 3, 4, 5, 2, 3, 2, 3, 4, 1, 2, 3, 4, 5, 6, 3, 4, 4, 5, 2, 3, 3, 4, 5, 6, 7, 6, 3, 4, 3, 4, 5, 6, 5, 6, 7
OFFSET
1,3
COMMENTS
The 'integral index' k of a rational polynomial p(x) is the smallest integer k such that p^[k](x) is an integer polynomial, where p^[k](x) denotes the k-th derivative of p. (Integer polynomials have integral index 0.) Using this way of speaking, the a(n) are the integral indices of the Bernoulli polynomials.
Conjecture: Every integer appears in this sequence only a finite number of times. (This generalizes the conjectures made in A366186-A366188.)
LINKS
EXAMPLE
B = Bernoulli(8, x).
B = -(1/30) + (2/3)*x^2 - (7/3)*x^4 + (14/3)*x^6 - 4*x^7 + x^8;
B' = (4/3)*x - (28/3)*x^3 + 28*x^5 - 28*x^6 + 8*x^7;
B'' = (4/3) - 28*x^2 + 140*x^4 - 168*x^5 + 56*x^6;
B''' = -56*x + 560*x^3 - 840*x^4 + 336*x^5.
Thus the integral index of B is a(8) = 3.
MAPLE
aList := proc(len) local n, k, d, A;
A := Array([seq(0, n = 0..len-1)]);
for n from 1 to len do
k := 0: d:= 0;
while d <> 1 do
k := k + 1;
d := denom(diff(bernoulli(n, x), `$`(x, k)));
od;
A[n] := k;
od;
convert(A, list) end:
aList(86);
MATHEMATICA
a[n_] := Module[{b, k, x}, b = BernoulliB[n, x]; For[k = 1, True, k++, b = D[b, x]; If[AllTrue[CoefficientList[b, x], IntegerQ], Return[k]]]];
Table[a[n], {n, 1, 100}] (* Jean-François Alcover, Oct 23 2023 *)
PROG
(Python)
from itertools import count
from sympy import Poly, bernoulli, diff
from sympy.abc import x
def A366189(n):
p = Poly(bernoulli(n, x))
for i in count(1):
p = diff(p)
if all(c.is_integer for c in p.coeffs()):
return i # Chai Wah Wu, Oct 03 2023
(SageMath)
def A366189List(len):
A = [0 for _ in range(len)]
P.<x> = ZZ[]
for n in range(len):
ber = bernoulli_polynomial(x, n + 1)
k = 0
while True:
k = k + 1
ber = diff(ber, x)
if ber.denominator() == 1:
A[n] = k; break
return A
print(A366189List(86)) # Peter Luschny, Oct 04 2023
CROSSREFS
Bernoulli polynomials: A196838/A196839 (with rising powers).
Cf. A094960 (m=1), A366169 (m=2), A366186 (m=3), A366187 (m=4), A366188 (m=5).
Sequence in context: A087039 A102096 A322813 * A137866 A361200 A329888
KEYWORD
nonn
AUTHOR
Peter Luschny, Oct 03 2023
STATUS
approved