Strong Pseudoprime

Author

Eric W. Weisstein
April 6, 2006

This notebook downloaded from http://mathworld.wolfram.com/notebooks/PrimeNumbers/StrongPseudoprime.nb.

For more information, see Eric's MathWorld entry http://mathworld.wolfram.com/StrongPseudoprime.html.

©2006 Wolfram Research, Inc. except for portions noted otherwise

Code

StrongPseudoprimeQ[a_,n_Integer]:=Module[{s,d},
n>2&&Mod[n,2]==1&&!PrimeQ[n]&&(PowerMod[a,d=(n-1)/2^(s=IntegerExponent[n-1,2]),n]==1||If[s==1,PowerMod[a,d,n]==n-1,(PowerMod[a,d #,n]==n-1)&/@(Or@@(2^Range[0,s-1]))
]
)
]

Enumeration

Table[{Prime[n],Select[Range[10^5],StrongPseudoprimeQ[Prime[n],#]&]},{n,10}]//ColumnForm

Table[{n,Select[Range[10^4],StrongPseudoprimeQ[n,#]&]},{n,10}]//ColumnForm

Smallest to small bases

First[spsp2=Select[Range[1*^8],StrongPseudoprimeQ[2,#]&]]//Timing

First[spsp3=Select[spsp2,StrongPseudoprimeQ[3,#]&]]//Timing

First[spsp5=Select[spsp3,StrongPseudoprimeQ[5,#]&]]//Timing

spsp7=Select[spsp5,StrongPseudoprimeQ[7,#]&]//Timing

Smallest strong pseudoprime to all of the first k primes taken as bases

list={2047,1373653,25326001,3215031751,2152302898747,3474749660383,341550071728321,341550071728321};

MapIndexed[{Prime[#2[[1]]],Function[f,StrongPseudoprimeQ[f,#1]]/@Prime[Range[#2[[1]]]]}&,list]//ColumnForm