

 This page uses content from Wikipedia. The original article was at Miller–Rabin primality test. The list of authors can be seen in the page history. As with Rosetta Code, the text of Wikipedia is available under the GNU FDL. (See links for details on variance)



The Miller–Rabin primality test or Rabin–Miller primality test is a primality test: an algorithm which determines whether a given number is prime or not.

The algorithm, as modified by Michael O. Rabin to avoid the generalized Riemann hypothesis, is a probabilistic algorithm.

The pseudocode, from Wikipedia is:

Input: n > 2, an odd integer to be tested for primality;
       k, a parameter that determines the accuracy of the test
Output: composite if n is composite, otherwise probably prime
write n − 1 as 2s·d with d odd by factoring powers of 2 from n − 1
LOOP: repeat k times:
   pick a randomly in the range [2, n − 1]
   x ← ad mod n
   if x = 1 or x = n − 1 then do next LOOP
   repeat s − 1 times:
      x ← x2 mod n
      if x = 1 then return composite
      if x = n − 1 then do next LOOP
   return composite
return probably prime

 The nature of the test involves big numbers, so the use of "big numbers" libraries (or similar features of the language of your choice) are suggested, but not mandatory.
 Deterministic variants of the test exist and can be implemented as extra (not mandatory to complete the task)


 
DIM mersenne(11)
mersenne(1)=7
mersenne(2)=31
mersenne(3)=127
mersenne(4)=8191
mersenne(5)=131071
mersenne(6)=524287
mersenne(7)=2147483647
mersenne(8)=2305843009213693951
mersenne(9)=618970019642690137449562111
mersenne(10)=162259276829213363391578010288127
mersenne(11)=170141183460469231731687303715884105727
 
 
dim SmallPrimes(1000)
data       2,      3,      5,      7,     11,     13,     17,     19,     23,     29
data      31,     37,     41,     43,     47,     53,     59,     61,     67,     71
data      73,     79,     83,     89,     97,    101,    103,    107,    109,    113
data     127,    131,    137,    139,    149,    151,    157,    163,    167,    173
data     179,    181,    191,    193,    197,    199,    211,    223,    227,    229
data     233,    239,    241,    251,    257,    263,    269,    271,    277,    281
data     283,    293,    307,    311,    313,    317,    331,    337,    347,    349
data     353,    359,    367,    373,    379,    383,    389,    397,    401,    409
data     419,    421,    431,    433,    439,    443,    449,    457,    461,    463
data     467,    479,    487,    491,    499,    503,    509,    521,    523,    541
data     547,    557,    563,    569,    571,    577,    587,    593,    599,    601
data     607,    613,    617,    619,    631,    641,    643,    647,    653,    659
data     661,    673,    677,    683,    691,    701,    709,    719,    727,    733
data     739,    743,    751,    757,    761,    769,    773,    787,    797,    809
data     811,    821,    823,    827,    829,    839,    853,    857,    859,    863
data     877,    881,    883,    887,    907,    911,    919,    929,    937,    941
data     947,    953,    967,    971,    977,    983,    991,    997,   1009,   1013
data    1019,   1021,   1031,   1033,   1039,   1049,   1051,   1061,   1063,   1069
data    1087,   1091,   1093,   1097,   1103,   1109,   1117,   1123,   1129,   1151
data    1153,   1163,   1171,   1181,   1187,   1193,   1201,   1213,   1217,   1223
data    1229,   1231,   1237,   1249,   1259,   1277,   1279,   1283,   1289,   1291
data    1297,   1301,   1303,   1307,   1319,   1321,   1327,   1361,   1367,   1373
data    1381,   1399,   1409,   1423,   1427,   1429,   1433,   1439,   1447,   1451
data    1453,   1459,   1471,   1481,   1483,   1487,   1489,   1493,   1499,   1511
data    1523,   1531,   1543,   1549,   1553,   1559,   1567,   1571,   1579,   1583
data    1597,   1601,   1607,   1609,   1613,   1619,   1621,   1627,   1637,   1657
data    1663,   1667,   1669,   1693,   1697,   1699,   1709,   1721,   1723,   1733
data    1741,   1747,   1753,   1759,   1777,   1783,   1787,   1789,   1801,   1811
data    1823,   1831,   1847,   1861,   1867,   1871,   1873,   1877,   1879,   1889
data    1901,   1907,   1913,   1931,   1933,   1949,   1951,   1973,   1979,   1987
data    1993,   1997,   1999,   2003,   2011,   2017,   2027,   2029,   2039,   2053
data    2063,   2069,   2081,   2083,   2087,   2089,   2099,   2111,   2113,   2129
data    2131,   2137,   2141,   2143,   2153,   2161,   2179,   2203,   2207,   2213
data    2221,   2237,   2239,   2243,   2251,   2267,   2269,   2273,   2281,   2287
data    2293,   2297,   2309,   2311,   2333,   2339,   2341,   2347,   2351,   2357
data    2371,   2377,   2381,   2383,   2389,   2393,   2399,   2411,   2417,   2423
data    2437,   2441,   2447,   2459,   2467,   2473,   2477,   2503,   2521,   2531
data    2539,   2543,   2549,   2551,   2557,   2579,   2591,   2593,   2609,   2617
data    2621,   2633,   2647,   2657,   2659,   2663,   2671,   2677,   2683,   2687
data    2689,   2693,   2699,   2707,   2711,   2713,   2719,   2729,   2731,   2741
data    2749,   2753,   2767,   2777,   2789,   2791,   2797,   2801,   2803,   2819
data    2833,   2837,   2843,   2851,   2857,   2861,   2879,   2887,   2897,   2903
data    2909,   2917,   2927,   2939,   2953,   2957,   2963,   2969,   2971,   2999
data    3001,   3011,   3019,   3023,   3037,   3041,   3049,   3061,   3067,   3079
data    3083,   3089,   3109,   3119,   3121,   3137,   3163,   3167,   3169,   3181
data    3187,   3191,   3203,   3209,   3217,   3221,   3229,   3251,   3253,   3257
data    3259,   3271,   3299,   3301,   3307,   3313,   3319,   3323,   3329,   3331
data    3343,   3347,   3359,   3361,   3371,   3373,   3389,   3391,   3407,   3413
data    3433,   3449,   3457,   3461,   3463,   3467,   3469,   3491,   3499,   3511
data    3517,   3527,   3529,   3533,   3539,   3541,   3547,   3557,   3559,   3571
data    3581,   3583,   3593,   3607,   3613,   3617,   3623,   3631,   3637,   3643
data    3659,   3671,   3673,   3677,   3691,   3697,   3701,   3709,   3719,   3727
data    3733,   3739,   3761,   3767,   3769,   3779,   3793,   3797,   3803,   3821
data    3823,   3833,   3847,   3851,   3853,   3863,   3877,   3881,   3889,   3907
data    3911,   3917,   3919,   3923,   3929,   3931,   3943,   3947,   3967,   3989
data    4001,   4003,   4007,   4013,   4019,   4021,   4027,   4049,   4051,   4057
data    4073,   4079,   4091,   4093,   4099,   4111,   4127,   4129,   4133,   4139
data    4153,   4157,   4159,   4177,   4201,   4211,   4217,   4219,   4229,   4231
data    4241,   4243,   4253,   4259,   4261,   4271,   4273,   4283,   4289,   4297
data    4327,   4337,   4339,   4349,   4357,   4363,   4373,   4391,   4397,   4409
data    4421,   4423,   4441,   4447,   4451,   4457,   4463,   4481,   4483,   4493
data    4507,   4513,   4517,   4519,   4523,   4547,   4549,   4561,   4567,   4583
data    4591,   4597,   4603,   4621,   4637,   4639,   4643,   4649,   4651,   4657
data    4663,   4673,   4679,   4691,   4703,   4721,   4723,   4729,   4733,   4751
data    4759,   4783,   4787,   4789,   4793,   4799,   4801,   4813,   4817,   4831
data    4861,   4871,   4877,   4889,   4903,   4909,   4919,   4931,   4933,   4937
data    4943,   4951,   4957,   4967,   4969,   4973,   4987,   4993,   4999,   5003
data    5009,   5011,   5021,   5023,   5039,   5051,   5059,   5077,   5081,   5087
data    5099,   5101,   5107,   5113,   5119,   5147,   5153,   5167,   5171,   5179
data    5189,   5197,   5209,   5227,   5231,   5233,   5237,   5261,   5273,   5279
data    5281,   5297,   5303,   5309,   5323,   5333,   5347,   5351,   5381,   5387
data    5393,   5399,   5407,   5413,   5417,   5419,   5431,   5437,   5441,   5443
data    5449,   5471,   5477,   5479,   5483,   5501,   5503,   5507,   5519,   5521
data    5527,   5531,   5557,   5563,   5569,   5573,   5581,   5591,   5623,   5639
data    5641,   5647,   5651,   5653,   5657,   5659,   5669,   5683,   5689,   5693
data    5701,   5711,   5717,   5737,   5741,   5743,   5749,   5779,   5783,   5791
data    5801,   5807,   5813,   5821,   5827,   5839,   5843,   5849,   5851,   5857
data    5861,   5867,   5869,   5879,   5881,   5897,   5903,   5923,   5927,   5939
data    5953,   5981,   5987,   6007,   6011,   6029,   6037,   6043,   6047,   6053
data    6067,   6073,   6079,   6089,   6091,   6101,   6113,   6121,   6131,   6133
data    6143,   6151,   6163,   6173,   6197,   6199,   6203,   6211,   6217,   6221
data    6229,   6247,   6257,   6263,   6269,   6271,   6277,   6287,   6299,   6301
data    6311,   6317,   6323,   6329,   6337,   6343,   6353,   6359,   6361,   6367
data    6373,   6379,   6389,   6397,   6421,   6427,   6449,   6451,   6469,   6473
data    6481,   6491,   6521,   6529,   6547,   6551,   6553,   6563,   6569,   6571
data    6577,   6581,   6599,   6607,   6619,   6637,   6653,   6659,   6661,   6673
data    6679,   6689,   6691,   6701,   6703,   6709,   6719,   6733,   6737,   6761
data    6763,   6779,   6781,   6791,   6793,   6803,   6823,   6827,   6829,   6833
data    6841,   6857,   6863,   6869,   6871,   6883,   6899,   6907,   6911,   6917
data    6947,   6949,   6959,   6961,   6967,   6971,   6977,   6983,   6991,   6997
data    7001,   7013,   7019,   7027,   7039,   7043,   7057,   7069,   7079,   7103
data    7109,   7121,   7127,   7129,   7151,   7159,   7177,   7187,   7193,   7207
data    7211,   7213,   7219,   7229,   7237,   7243,   7247,   7253,   7283,   7297
data    7307,   7309,   7321,   7331,   7333,   7349,   7351,   7369,   7393,   7411
data    7417,   7433,   7451,   7457,   7459,   7477,   7481,   7487,   7489,   7499
data    7507,   7517,   7523,   7529,   7537,   7541,   7547,   7549,   7559,   7561
data    7573,   7577,   7583,   7589,   7591,   7603,   7607,   7621,   7639,   7643
data    7649,   7669,   7673,   7681,   7687,   7691,   7699,   7703,   7717,   7723
data    7727,   7741,   7753,   7757,   7759,   7789,   7793,   7817,   7823,   7829
data    7841,   7853,   7867,   7873,   7877,   7879,   7883,   7901,   7907,   7919
 
 
print "Liberty Miller Rabin Demonstration"
print "Loading Small Primes"
for i=1 to 1000:   read x : SmallPrimes(i)=x :next :NoOfSmallPrimes=1000
print NoOfSmallPrimes;" Primes Loaded"
 
'Prompt "Enter number to test:";resp$
'x=val(resp$)
'goto [Jump]
 
 
For i=1 to 11
 
 x=mersenne(i)
 
 
 t1=time$("ms")
 [TryAnother]
 print
 
 iterations=0
 [Loop]
    iterations=iterations+1
 
    if MillerRabin(x,7)=1 then
     t2=time$("ms")
     print "Composite, found in ";t2-t1;" milliseconds"
    else
     t2=time$("ms")
     print x;" Probably Prime. Tested in ";t2-t1;" milliseconds"
     playwave "tada.wav", async
 end if
 print
 
next
 
END
 
 
Function GCD( m,n )
' Find greatest common divisor with Extend Euclidian Algorithm
' Knuth Vol 1 P.13 Algorithm E
 
ap =1 :b  =1 :a  =0 :bp =0: c  =m :d  =n
 
[StepE2]
q  = int(c/d) :r  = c-q*d
 
if r<>0 then
    c=d :d=r :t=ap :ap=a :a=t-q*a :t=bp  :bp=b  :b=t-q*b
    'print ap;" ";b;" ";a;" ";bp;" ";c;" ";d;" ";t;" ";q
    goto [StepE2]
end if
 
GCD=a*m+b*n
 
'print ap;" ";b;" ";a;" ";bp;" ";c;" ";d;" ";t;" ";q
 
End Function 'Extended Euclidian GCD
 
 function IsEven( x )
    if ( x MOD 2 )=0 then
        IsEven=1
    else
         IsEven=0
    end if
end function
 
 
function IsOdd( x )
    if ( x MOD 2 )=0 then
        IsOdd=0
    else
        IsOdd=1
    end if
end function
 
 
Function FastExp(x, y, N)
 
  if (y=1) then                  'MOD(x,N)
      FastExp=x-int(x/N)*N
      goto [ExitFunction]
  end if
 
 
  if ( y and 1) = 0  then
 
     dum1=y/2
     dum2=y-int(y/2)*2              'MOD(y,2)
 
     temp=FastExp(x,dum1,N)
     z=temp*temp
     FastExp=z-int(z/N)*N            'MOD(temp*temp,N)
     goto [ExitFunction]
  else
 
     dum1=y-1
     dum1=dum1/2
     temp=FastExp(x,dum1,N)
     dum2=temp*temp
     temp=dum2-int(dum2/N)*N            'MOD(dum2,N)
 
     z=temp*x
     FastExp=z-int(z/N)*N             'MOD(temp*x,N)
     goto [ExitFunction]
  end if
  [ExitFunction]
 
end function
 
 
Function MillerRabin(n,b)
 
'print "Miller Rabin"
't1=time$("ms")
 
  if IsEven(n) then
    MillerRabin=1
    goto [ExtFn]
  end if
 
  i=0
  [Loop]
    i=i+1
    if i>1000 then goto [Continue]
    if ( n MOD SmallPrimes(i) )=0 then
      MillerRabin=0
      goto [ExtFn]
    end if
  goto [Loop]
  [Continue]
 
  if GCD(n,b)>1 then
    MillerRabin=1
    goto [ExtFn]
  end if
 
  q=n-1
 
  t=0
 
  while  (int(q) AND 1 )=0
   t=t+1
   q=int(q/2)
  wend
 
 
  r=FastExp(b, q, n)
 
  if ( r <> 1 ) then
    e=0
    while ( e < (t-1) )
      if ( r <> (n-1) ) then
        r=FastExp(r, r, n)
        else
        Exit While
      end if
 
      e=e+1
    wend
    [ExitLoop]
  end if
 
 
  if ( (r=1) OR (r=(n-1)) ) then
      MillerRabin=0
    else
      MillerRabin=1
  end if
 
[ExtFn]
 
End Function
 