OFFSET
1,1
COMMENTS
Number of integers b, 1 <= b <= 2n, such that if 2n = 2^k*m with odd m, then the sequence (b^m, b^(2*m), ..., b^(2^k*m)) modulo 2n+1 satisfies the Rabin-Miller test.
Comments from R. J. Mathar, Jul 03 2012 (Start)
The subsequence related to composite 2n+1 is characterized with records in A195328 and associated 2n+1 tabulated in A141768.
Let N = 2n+1 = product_{i=1..s} p_i^r_i be the prime factorization of the odd 2n+1. Related odd parts q and q_i are defined by N-1=2^k*q and p_i-1 = 2^(k_i)*q_i, with sorting such that k_1 <= k_2 <=k_3... Then a(n) = (1+sum_{j=0..k1-1} 2^(j*s)) *product_{i=1..s} gcd(q,qi).
Reduces to A006093 if 2n+1 is prime.
This might be correlated with 2*A195508(n). (End)
REFERENCES
Paulo Ribenboim, The Little Book of Bigger Primes, 2nd ed., Springer-Verlag, New York, 2004, p. 98.
LINKS
Amiram Eldar, Table of n, a(n) for n = 1..10000
F. Arnault, The Rabin-Monier theorem for Lucas pseudoprimes, Mathematics of Computation, vol. 66, no 218, April 1997, pp. 869-881.
Louis Monier, Evaluation and comparison of two efficient primality testing algorithms, Theoretical Computer Science, Vol. 11 (1980), pp. 97-108.
Wikipedia, Strong pseudoprime
FORMULA
For k = 2*n+1, a(k) = k - 1 if k is prime, otherwise, a(k) = (1 + 2^(omega(k)*nu(k)) - 1)/(2^omega(k)-1)) * Product_{p|k} gcd(od(k-1), od(p-1)), where omega(m) is the number of distinct prime factors of m (A001221), od(m) is the largest odd divisor of m (A000265) and nu(m) = min_{p|m} A007814(p-1). - Amiram Eldar, Nov 08 2019
MAPLE
rabinmiller := proc(n, a); k := 0; mu := n-1; while irem(mu, 2)=0 do k := k+1; mu := mu/2 od; G := a&^mu mod(n); h := 0; if G=1 then RETURN(1) else while h<k-1 and G&^2 mod n <>1 do h := h+1; G := G&^2 mod n; od; if h<k and G<> n-1 then RETURN(0) else RETURN(1) fi; if G=1 then RETURN(1); fi; fi; end; compte := proc(n) local l; RETURN(sum('rabinmiller(2*n+1, l)', 'l'=1..2*n)); end;
Maple code from R. J. Mathar, Jul 03 2012 (Start)
A000265 := proc(n)
n/2^padic[ordp](n, 2) ;
end proc:
a := proc(n)
q := A000265(n-1) ;
B := 1;
s := 0 ;
k1 := 10000000000000 ;
for pf in ifactors(n)[2] do
pi := op(1, pf) ;
qi := A000265(pi-1) ;
ki := ilog2((pi-1)/qi) ;
k1 := min(k1, ki) ;
B := B*igcd(q, qi) ;
s := s+1 ;
end do:
1+add(2^(j*s), j=0..k1-1) ;
return B*% ;
end proc:
seq(a(2*n+1), n=1..60) ;
MATHEMATICA
o[n_] := (n-1)/2^IntegerExponent[n-1, 2]; a[n_?PrimeQ] := n-1; a[n_] := Module[{p = FactorInteger[n][[;; , 1]]}, om = Length[p]; Product[GCD[o[n], o[p[[k]]]], {k, 1, om}] * (1 + (2^(om * Min[IntegerExponent[#, 2]& /@ (p - 1)]) - 1)/(2^om - 1))]; Table[a[n], {n, 3, 121, 2}] (* Amiram Eldar, Nov 08 2019 *)
CROSSREFS
KEYWORD
nonn
AUTHOR
J.-F. Guiffes (guiffes.jean-francois(AT)wanadoo.fr), Jun 11 2002
EXTENSIONS
Edited by Max Alekseyev, Sep 20 2018
STATUS
approved