|
%I #17 Aug 04 2018 03:33:49
%S 25,175,553,949,1541,1975,2701,4537,6931,7501,13825,21349,25273,25477,
%T 29341,38323,42121,42127,47617,49141,52417,63701,80137,83333,85609,
%U 94753,105121,113527,128143,144841,167869,169027,179521,181351,197209,201373
%N Strong pseudoprimes to base 24.
%H R. J. Mathar, <a href="/A020250/b020250.txt">Table of n, a(n) for n = 1..203</a>
%H <a href="/index/Ps#pseudoprimes">Index entries for sequences related to pseudoprimes</a>
%t strongPseudoprimeQ[b_, n_] := Module[{rems = Table[PowerMod[b, (n - 1)/2^expo, n], {expo, 0, IntegerExponent[n - 1, 2]}]}, (rems[[-1]] == 1 || MemberQ[rems, n - 1]) && PowerMod[b, n - 1, n] == 1]; max = 5000; Select[Complement[Range[2, max], Prime[Range[PrimePi[max]]]], strongPseudoprimeQ[24, #] &] (* _Alonso del Arte_, Aug 03 2018 *)
%o (PARI) is(n, a=24)= (bittest(n, 0) && !isprime(n) && n>8) || return; my(s=valuation(n-1, 2)); if(1==a=Mod(a, n)^(n>>s), return(1)); while(a!=-1 && s--, a=a^2); a==-1 \\ _Felix Fröhlich_, Aug 03 2018, adapted from code by _M. F. Hasler_ in A001262
%K nonn
%O 1,1
%A _David W. Wilson_
|
