A214686 - OEIS

%I #19 Jul 27 2012 12:32:12

%S 1,1,7,1,23,1,47,1,79,1,113,89,23,73,31,1,283,89,113,139,173,67,47,1,

%T 619,131,109,83,113,211,191,1,1087,1,1223,1,1367,1,1511,367,83,1,1847,

%U 1,2017,317,571,241,199,1,2593,367,211,271,223,1,3229,1117,239,1

%N Numerators of a series with denominators n! and sum 1.

%C a(n) is the greatest integer x such that gcd(x,n!) = 1 and x/n! < 1 - sum_{j=2}^{n-1} a(j)/j!.

%C The infinite series sum_{n=2}^infinity a(n)/n! = 1

%C For each n, either a(n) = 1 or a(n) >= n+1.

%H Robert Israel and T. D. Noe, <a href="/A214686/b214686.txt">Table of n, a(n) for n = 2..10000</a> (first 1000 terms from Robert Israel)

%e 1 - a(2)/2! = 1/2 = 3/3!, gcd(2,3!)>1 so a(3) = 1.

%e 1 - a(2)/2! - a(3)/3! = 8/4! so a(4) = 7.

%p N:= 100; a[2]:= 1; R[2]:= 1/2;

%p for j from 3 to N do

%p    T:= R[j-1] *j!;

%p    for x from T-1 by -1 while igcd(x,j!) > 1 do end do;

%p    a[j]:= x;

%p    R[j]:= R[j-1] - x/j!

%p end do:

%p seq(a[j],j=2..N);

%t s = 1; Table[d = n! s; q = If[d <= n, 1, If[d <= n^2, Prime[PrimePi[d]], Print["d > n^2"]; Abort[]]]; s = s - q/n!; q, {n, 2, 100}] (* _T. D. Noe_, Jul 27 2012 *)

%o (Sage)

%o def A214686_list(n) :

%o     a = [1]; R = 1/2

%o     for j in (3..n+1) :

%o         J = factorial(j)

%o         T = R * J

%o         for x in range(T-1, -1, -1) :

%o             if gcd(x, J) == 1 : break

%o         a.append(x)

%o         R = R - x / J

%o     return a

%o A214686_list(51) # _Peter Luschny_, Jul 27 2012

%K nonn

%O 2,3

%A _Robert Israel_, Jul 25 2012