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.

Showing 1-5 of 5 results.

A195441 a(n) = denominator(Bernoulli_{n+1}(x) - Bernoulli_{n+1}).

Original entry on oeis.org

1, 1, 2, 1, 6, 2, 6, 3, 10, 2, 6, 2, 210, 30, 6, 3, 30, 10, 210, 42, 330, 30, 30, 30, 546, 42, 14, 2, 30, 2, 462, 231, 3570, 210, 6, 2, 51870, 2730, 210, 42, 2310, 330, 2310, 210, 4830, 210, 210, 210, 6630, 1326, 858, 66, 330, 110, 798, 114, 870, 30, 30, 6
Offset: 0

Views

Author

Peter Luschny, Sep 18 2011

Keywords

Comments

If s(n) is the smallest number such that s(n)*(1^n + 2^n + ... + x^n) is a polynomial in x with integer coefficients then a(n)=s(n)/(n+1) (see A064538).
a(n) is squarefree, by the von Staudt-Clausen theorem on the denominators of Bernoulli numbers. - Kieren MacMillan and Jonathan Sondow, Nov 20 2015
Kellner and Sondow give a detailed analysis of this sequence and provide a simple way to compute the terms without using Bernoulli polynomials and numbers. They prove that a(n) is the product of the primes less than or equal to (n+2)/(2+(n mod 2)) such that the sum of digits of n+1 in base p is at least p. - Peter Luschny, May 14 2017
The equation a(n-1) = denominator(Bernoulli_n(x) - Bernoulli_n) = rad(n+1) has only finitely many solutions, where rad(n) = A007947(n) is the radical of n. It is conjectured that S = {3, 5, 8, 9, 11, 27, 29, 35, 59} is the full set of all such solutions. Note that (S\{8})+1 joined with {1,2} equals A094960. More precisely, the set S implies the finite sequence of A094960. See Kellner 2023. - Bernd C. Kellner, Oct 18 2023
As was observed in the example section of A318256: denominator(B_n(x)) = rad(n+1) if n is in {0, 1, 3, 5, 9, 11, 27, 29, 35, 59} = {A094960(n) - 1: 1 <= n <= 10}. - Peter Luschny, Oct 18 2023

Crossrefs

Programs

  • Julia
    using Nemo, Primes
    function A195441(n::Int)
        n < 4 && return ZZ([1,1,2,1][n+1])
        P = primes(2, div(n+2, 2+n%2))
        prod([ZZ(p) for p in P if p <= sum(digits(n+1, base=p))])
    end
    println([A195441(n) for n in 0:59]) # Peter Luschny, May 14 2017
    
  • Maple
    A195441 := n -> denom(bernoulli(n+1, x)-bernoulli(n+1)):
    seq(A195441(i),i=0..59);
    # Formula of Kellner and Sondow:
    a := proc(n) local s; s := (p,n) -> add(i,i=convert(n,base,p));
    select(isprime,[$2..(n+2)/(2+irem(n,2))]); mul(i,i=select(p->s(p,n+1)>=p,%)) end: seq(a(n), n=0..59); # Peter Luschny, May 14 2017
  • Mathematica
    a[n_] := Denominator[Together[(BernoulliB[n + 1, x] - BernoulliB[n + 1])]]; Table[a[n], {n, 0, 59}] (* Jonathan Sondow, Nov 20 2015 *)
    SD[n_, p_] := If[n < 1 || p < 2, 0, Plus @@ IntegerDigits[n, p]]; DD[n_] := Times @@ Select[Prime[Range[PrimePi[(n+2)/(2+Mod[n, 2])]]], SD[n+1, #] >= # &]; Table[DD[n], {n, 0, 59}] (* Bernd C. Kellner, Oct 18 2023 *)
  • PARI
    a(n) = {my(vp = Vec(bernpol(n+1, x)-bernfrac(n+1))); lcm(vector(#vp, k, denominator(vp[k])));} \\ Michel Marcus, Feb 08 2016
    
  • Python
    from math import prod
    from sympy.ntheory.factor_ import primerange, digits
    def A195441(n): return prod(p for p in primerange((n+2)//(2|n&1)+1) if sum(digits(n+1,p)[1:])>=p) # Chai Wah Wu, Oct 04 2023
  • Sage
    A195441 = lambda n: mul([p for p in (2..(n+2)//(2+n%2)) if is_prime(p) and sum((n+1).digits(base=p))>=p])
    print([A195441(n) for n in (0..59)]) # Peter Luschny, May 14 2017
    

Formula

a(n) = A064538(n)/(n+1). - Jonathan Sondow, Nov 12 2015
A001221(a(n)) = A001222(a(n)). - Kieren MacMillan and Jonathan Sondow, Nov 20 2015
a(2*n)/a(2*n+1) = A286516(n+1). - Bernd C. Kellner and Jonathan Sondow, May 24 2017
a(n) = A007947(A338025(n+1)). - Harald Hofstätter, Oct 10 2020
From Bernd C. Kellner, Oct 18 2023: (Start)
Note that the formulas here are shifted in index by 1 due to the definition of a(n) using index n+1!
a(n) = A324369(n+1) * A324370(n+1).
a(n) = A144845(n) / A324371(n+1).
a(n-1) = lcm(a(n), rad(n+1)), if n >= 3 is odd.
If n+1 is composite, then rad(n+1) divides a(n-1).
If m is a Carmichael number (A002997), then m divides both a(m-1) and a(m-2).
See papers of Kellner and Kellner & Sondow. (End)

Extensions

Definition simplified by Jonathan Sondow, Nov 20 2015

A064538 a(n) is the smallest positive integer such that a(n)*(1^n + 2^n + ... + x^n) is a polynomial in x with integer coefficients.

Original entry on oeis.org

1, 2, 6, 4, 30, 12, 42, 24, 90, 20, 66, 24, 2730, 420, 90, 48, 510, 180, 3990, 840, 6930, 660, 690, 720, 13650, 1092, 378, 56, 870, 60, 14322, 7392, 117810, 7140, 210, 72, 1919190, 103740, 8190, 1680, 94710, 13860, 99330, 9240, 217350, 9660, 9870, 10080, 324870
Offset: 0

Views

Author

Floor van Lamoen, Oct 08 2001

Keywords

Comments

a(n) is a multiple of n+1. - Vladimir Shevelev, Dec 20 2011
Let P_n(x) = 1^n + 2^n + ... + x^n = Sum_{i=1..n+1}c_i*x^i. Let P^*n(x) = Sum{i=1..n+1}(c_i/(i+1))*(x^(i+1)-x). Then b(n) = (n+1)*a(n+1)is the smallest positive integer such that b(n)*P^*n(x) is a polynomial with integer coefficients. Proof follows from the recursion P(n+1)(x) = x + (n+1)*P^*n(x). As a corollary, note that, if p is the maximal prime divisor of a(n), then p<=n+1. - _Vladimir Shevelev, Dec 21 2011
The recursion P_(n+1)(x) = x + (n+1)*P^*n(x) is due to Abramovich (1973); see also Shevelev (2007). - _Jonathan Sondow, Nov 16 2015
The sum S_m(n) = Sum_{k=0..n} k^m can be written as S_m(n) = n(n+1)(2n+1)P_m(n)/a(m) for even m>1, or S_m(n) = n^2*(n+1)^2*P_m(n)/a(m) for odd m>1, where a(m) is the LCM of the denominators of the coefficients of the polynomial P_m/a(m), i.e., the smallest integer such that P_m defined in this way has integer coefficients. (Cf. Michon link.) - M. F. Hasler, Mar 10 2013
a(n)/(n+1) is squarefree, by Faulhaber's formula and the von Staudt-Clausen theorem on the denominators of Bernoulli numbers. - Kieren MacMillan and Jonathan Sondow, Nov 20 2015
a(n) equals n+1 times the product of the primes p <= (n+2)/(2+(n mod 2)) such that the sum of the base-p digits of n+1 is at least p. - Bernd C. Kellner and Jonathan Sondow, May 24 2017

Examples

			1^3 + 2^3 + ... + x^3 = (x(x+1))^2/4 so a(3)=4.
1^4 + 2^4 + ... + x^4 = x(x+1)(2x+1)(3x^2+3x-1)/30, therefore a(4)=30.
		

References

  • M. Abramowitz and I. A. Stegun, eds., Handbook of Mathematical Functions, National Bureau of Standards Applied Math. Series 55, 1964 (and various reprints), p. 804, Eq. 23.1.4.

Crossrefs

Programs

  • Maple
    A064538 := n -> denom((bernoulli(n+1,x)-bernoulli(n+1))/(n+1)): # Peter Luschny, Aug 19 2011
    # Formula of Kellner and Sondow (2017):
    a := proc(n) local s; s := (p,n) -> add(i,i=convert(n,base,p));
    select(isprime,[$2..(n+2)/(2+irem(n,2))]);
    (n+1)*mul(i,i=select(p->s(p,n+1)>=p,%)) end: seq(a(n), n=0..48); # Peter Luschny, May 14 2017
  • Mathematica
    A064538[n_] := Denominator[ Together[ (BernoulliB[n+1, x] - BernoulliB[n+1])/(n+1)]];
    Table[A064538[n], {n, 0, 44}] (* Jean-François Alcover, Feb 21 2012, after Maple *)
  • PARI
    a(n) = {my(vp = Vec(bernpol(n+1, x)-bernfrac(n+1))/(n+1)); lcm(vector(#vp, k, denominator(vp[k])));} \\ Michel Marcus, Feb 07 2016
    
  • Python
    from _future_ import division
    from sympy.ntheory.factor_ import digits, nextprime
    def A064538(n):
        p, m = 2, n+1
        while p <= (n+2)//(2+ (n% 2)):
            if sum(d for d in digits(n+1,p)[1:]) >= p:
                m *= p
            p = nextprime(p)
        return m # Chai Wah Wu, Mar 07 2018
  • Sage
    A064538 = lambda n: (n+1)*mul([p for p in (2..(n+2)//(2+n%2)) if is_prime(p) and sum((n+1).digits(base=p)) >= p])
    print([A064538(n) for n in (0..48)]) # Peter Luschny, May 14 2017
    

Formula

a(n) = (n+1)*A195441(n). - Jonathan Sondow, Nov 12 2015
A001221(a(n)/(n+1)) = A001222(a(n)/(n+1)). - Kieren MacMillan and Jonathan Sondow, Nov 20 2015
rad(a(n)) = A007947(a(n)) = A144845(n) = A324369(n+1) * A324370(n+1) * A324371(n+1). - Bernd C. Kellner, Oct 12 2023

A286516 a(n) = b(2*n-1)/b(2*n) where b(n) = A195441(n-1) = denominator(Bernoulli_{n}(x) - Bernoulli_{n}).

Original entry on oeis.org

1, 2, 3, 2, 5, 3, 7, 2, 3, 5, 11, 1, 13, 7, 15, 2, 17, 3, 19, 5, 7, 11, 23, 1, 5, 13, 3, 7, 29, 5, 31, 2, 11, 17, 7, 1, 37, 19, 13, 5, 41, 21, 43, 11, 3, 23, 47, 1, 7, 5, 17, 13, 53, 3, 11, 7, 19, 29, 59, 1, 61, 31, 7, 2, 65, 11, 67, 17, 23, 5, 71, 1, 73, 37
Offset: 1

Views

Author

Keywords

Comments

a(n) is an integer for all n, a(n) is odd if n is not a power of 2, a(2^k)=2 for all k>=1, a(n)=1 infinitely often, and a(n)=p infinitely often for every prime p. See Cor. 2 and Cor. 3 in "The denominators of power sums of arithmetic progressions". See also "Power-sum denominators".

Crossrefs

Programs

  • Mathematica
    b[n_] := Denominator[ Together[ BernoulliB[n, x] - BernoulliB[n]]]; Table[
    b[2 n - 1]/b[2 n], {n, 1, 74}]

Formula

a(n) = A195441(2*n-2) / A195441(2*n-1).

A318256 a(n) = (denominator of B(n,x)) / (the squarefree kernel of n+1), where B(n,x) is the n-th Bernoulli polynomial.

Original entry on oeis.org

1, 1, 2, 1, 6, 1, 6, 3, 10, 1, 6, 1, 210, 15, 2, 3, 30, 5, 210, 21, 110, 15, 30, 5, 546, 21, 14, 1, 30, 1, 462, 231, 1190, 105, 6, 1, 51870, 1365, 70, 21, 2310, 55, 2310, 105, 322, 105, 210, 35, 6630, 663, 286, 33, 330, 55, 798, 57, 290, 15, 30, 1, 930930, 15015
Offset: 0

Views

Author

Peter Luschny, Sep 12 2018

Keywords

Examples

			a(59) = 1 because there exist no number which satisfies the definition (and the product of an empty set is 1).
a(60) = 930930 because {2, 3, 5, 7, 11, 13, 31} are the only primes which satisfy the definition.
The denominator of the Bernoulli polynomial B_n(x) equals the squarefree kernel of n+1 if n is in {0, 1, 3, 5, 9, 11, 27, 29, 35, 59}. These might be the only numbers with this property.
		

Crossrefs

a(n) = A144845(n) / A007947(n+1).
Cf. A324370 (same sequence with offset 1).

Programs

  • Maple
    a := n -> denom(bernoulli(n, x)) / mul(p, p in numtheory:-factorset(n+1)):
    seq(a(n), n=0..61);
  • Mathematica
    sfk[n_] := Times @@ FactorInteger[n][[All, 1]];
    a[n_] := (BernoulliB[n, x] // Together // Denominator)/sfk[n+1];
    Table[a[n], {n, 0, 61}] (* Jean-François Alcover, Feb 14 2019 *)
  • Sage
    def A318256(n): return mul([p for p in (2..(n+2)//(2+n%2))
                    if is_prime(p)
                    and not p.divides(n+1)
                    and sum((n+1).digits(base=p)) >= p])
    print([A318256(n) for n in (0..61)])

Formula

Let Q(n) = {p <= floor((n + 2)/(2 + n mod 2)) and p is prime and p does not divide n + 1 and the sum of the digits in base p of n+1 is at least p} then a(n) = Product_{p in Q(n)} p. (See the Kellner & Sondow links.)
a(n) = denominator(Bernoulli'(n+1, x)), where ' denotes d/dx. - Peter Luschny, Oct 15 2023

A286763 Numbers that appear in A195441 at least once for two consecutive indices.

Original entry on oeis.org

1, 30, 210, 330, 2310, 3990, 6090, 14790, 43890, 66990, 82110, 125970, 144210, 181830, 881790, 1009470, 1067430, 1217370, 2284590, 2381190, 17687670, 18888870, 26265030, 35068110, 39544890, 47763870, 115223790, 127652070, 406816410, 497668710, 741110370, 1024748670
Offset: 1

Views

Author

Peter Luschny, May 14 2017

Keywords

Comments

The sequence is infinite; see Cor. 3 in "The denominators of power sums of arithmetic progressions". - Bernd C. Kellner and Jonathan Sondow, May 24 2017

Examples

			A195441(21) = A195441(22) = 30, so 30 is in the sequence. - _Jonathan Sondow_, Dec 11 2018
		

Crossrefs

Programs

  • Julia
    function A286763_search()
        A = fmpz[]; a = fmpz(0)
        for n in 0:10000
            u = A195441(n)
            a == u && push!(A, a)
            a = u
        end
        S = sort([a for a in Set(A)])
    S[1:32] end
    println(A286763_search())
  • Mathematica
    Take[#, 32] &@ Union@ SequenceCases[ Table[ Denominator[ Together[ (BernoulliB[n + 1, x] - BernoulliB[n + 1])]], {n, 0, 2000}], w_ /; And[SameQ @@ w, Length@ w >= 2]][[All, 1]] (* Michael De Vlieger, Sep 22 2017, after Jonathan Sondow at A195441 *)
Showing 1-5 of 5 results.