Re Bernoulli numbers - more terms Pari

cino hilliard hillcino368 at hotmail.com
Wed Feb 4 01:32:18 CET 2004


Hi, Neal
GP/PARI CALCULATOR Version 2.2.8 (development CHANGES-1.887)
       i686 running cygwin (ix86 kernel) 32-bit version
      compiled: Jan 13 2004, gcc-3.3.1 (cygming special)
     (readline v4.3 enabled, extended help not available)

system pIV 2.533ghz 2 gig ram xp pro

Here are some more terms in Pari. You will have to allocatemen() for m1 or 
m2 > 7000
I had the set for below m1=1,m2 = 10000 but knew it would crash for more 
memory so
I stopped.  BTW, you don't need arrays for this sequence. Some speed gain is 
realized also
by setting 2*n to n2 befor the bernfrac() call.

\\ c = start index of next term not computed,m1 = start, me=end
bern2(c,m1,m2) =
       {
       for(n=m1,m2,
          n2=n+n;
          a = numerator(bernfrac(n2)/(n2));           \\ A001067
          b = numerator(bernfrac(n2)/(n2*(n2-1)));    \\ A046968
          if(a <> b,print("A("c") = "n","a/b);c++)
          )
}

1,574,37
2,1185,103
3,1240,37
4,1269,59
5,1376,131
6,1906,37
7,1910,67
8,2572,37
9,2689,283
10,2980,59
11,3238,37
12,3384,101
13,3801,691
14,3904,37
15,4121,67
16,4570,37
17,4691,59
18,4789,157
19,5236,37

system pIII 1 ghz 256k ram xp pro
This broke down at 4570 giving a huge number as the result. must be memory 
pr p3 /pari bug.
2689,4691,4789  are primes.

Looks like you have 1 more yet one more sequence.

Cino




>From: "N. J. A. Sloane" <njas at research.att.com>
>Reply-To: njas at research.att.com
>To: eclark at math.usf.edu, somos at grail.cba.csuohio.edu, seqfan at ext.jussieu.fr
>Subject: Re Bernoulli numbers
>Date: Tue, 3 Feb 2004 15:39:10 -0500 (EST)
>
>Thanks, Edwin!
>
>%I A090495
>%S A090495 574,1185,1240,1269,1376,1906,1910
>%N A090495 Numbers n such that numerator(Bernoulli(2*n)/(2*n)) is different 
>from numerator(Bernoulli(2*n)/(2*n*(2*n-1))).
>%C A090495 n such that A001067 is different from A046968.
>%O A090495 1,1
>%K A090495 nonn,more
>%Y A090495 Cf. A090496, A001067, A046968.
>%p A090495 a:=n->numer(bernoulli(2*n)/(2*n)): 
>b:=n->numer(bernoulli(2*n)/(2*n*(2*n-1))): for n from 1 to 2000 do if 
>a(n)<>b(n) then print(n,a(n)/b(n));  fi; od:
>%t A090495 a[n_] := Numerator[BernoulliB[2n]/(2n)] (* A001067 *); b[n_] := 
>Numerator[BernoulliB[2n]/(2n(2n-1))] (* A046968 *); For[n=1, n <= 580, n++, 
>If[ a[n] != b[n], Print[n, " ", a[n]/b[n]] ] ]
>%A A090495 njas, Feb 03 2004
>%E A090495 a(1) discovered by Michael Somos, Feb 01 2004. a(2)-a(7) from 
>Edwin Clark, Feb 03, 2004.
>
>%I A090496
>%S A090496 37,103,37,59,131,37,67
>%N A090496 Ratio of numerator(Bernoulli(2*n)/(2*n)) to 
>numerator(Bernoulli(2*n)/(2*n*(2*n-1))) for n's for which they are 
>different.
>%C A090496 A001067(n) / A046968(n) when they are different.
>%O A090496 1,1
>%K A090496 nonn,more
>%Y A090496 Cf. A090496, A001067, A046968.
>%A A090496 njas, Feb 03 2004
>%E A090496 a(1) discovered by Michael Somos, Feb 01 2004. a(2)-a(7) from 
>Edwin Clark, Feb 03, 2004.
>
>NJAS

_________________________________________________________________
Scope out the new MSN Plus Internet Software — optimizes dial-up to the max! 
   http://join.msn.com/?pgmarket=en-us&page=byoa/plus&ST=1






More information about the SeqFan mailing list