mcs/class/Mono.Security/Mono.Math.Prime/PrimalityTests.cs

   1 //
   2 // Mono.Math.Prime.PrimalityTests.cs - Test for primality
   3 //
   4 // Authors:
   5 //      Ben Maurer
   6 //
   7 // Copyright (c) 2003 Ben Maurer. All rights reserved
   8 //
   9
  10 //
  11 // Permission is hereby granted, free of charge, to any person obtaining
  12 // a copy of this software and associated documentation files (the
  13 // "Software"), to deal in the Software without restriction, including
  14 // without limitation the rights to use, copy, modify, merge, publish,
  15 // distribute, sublicense, and/or sell copies of the Software, and to
  16 // permit persons to whom the Software is furnished to do so, subject to
  17 // the following conditions:
  18 //
  19 // The above copyright notice and this permission notice shall be
  20 // included in all copies or substantial portions of the Software.
  21 //
  22 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
  23 // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
  24 // MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
  25 // NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
  26 // LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
  27 // OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
  28 // WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
  29 //
  30
  31 using System;
  32
  33 namespace Mono.Math.Prime {
  34
  35 #if INSIDE_CORLIB
  36         internal
  37 #else
  38         public
  39 #endif
  40         delegate bool PrimalityTest (BigInteger bi, ConfidenceFactor confidence);
  41
  42 #if INSIDE_CORLIB
  43         internal
  44 #else
  45         public
  46 #endif
  47         sealed class PrimalityTests {
  48
  49                 private PrimalityTests ()
  50                 {
  51                 }
  52
  53                 #region SPP Test
  54
  55                 private static int GetSPPRounds (BigInteger bi, ConfidenceFactor confidence)
  56                 {
  57                         int bc = bi.BitCount();
  58
  59                         int Rounds;
  60
  61                         // Data from HAC, 4.49
  62                         if      (bc <= 100 ) Rounds = 27;
  63                         else if (bc <= 150 ) Rounds = 18;
  64                         else if (bc <= 200 ) Rounds = 15;
  65                         else if (bc <= 250 ) Rounds = 12;
  66                         else if (bc <= 300 ) Rounds =  9;
  67                         else if (bc <= 350 ) Rounds =  8;
  68                         else if (bc <= 400 ) Rounds =  7;
  69                         else if (bc <= 500 ) Rounds =  6;
  70                         else if (bc <= 600 ) Rounds =  5;
  71                         else if (bc <= 800 ) Rounds =  4;
  72                         else if (bc <= 1250) Rounds =  3;
  73                         else                 Rounds =  2;
  74
  75                         switch (confidence) {
  76                                 case ConfidenceFactor.ExtraLow:
  77                                         Rounds >>= 2;
  78                                         return Rounds != 0 ? Rounds : 1;
  79                                 case ConfidenceFactor.Low:
  80                                         Rounds >>= 1;
  81                                         return Rounds != 0 ? Rounds : 1;
  82                                 case ConfidenceFactor.Medium:
  83                                         return Rounds;
  84                                 case ConfidenceFactor.High:
  85                                         return Rounds << 1;
  86                                 case ConfidenceFactor.ExtraHigh:
  87                                         return Rounds << 2;
  88                                 case ConfidenceFactor.Provable:
  89                                         throw new Exception ("The Rabin-Miller test can not be executed in a way such that its results are provable");
  90                                 default:
  91                                         throw new ArgumentOutOfRangeException ("confidence");
  92                         }
  93                 }
  94
  95                 /// <summary>
  96                 ///     Probabilistic prime test based on Rabin-Miller's test
  97                 /// </summary>
  98                 /// <param name="bi" type="BigInteger.BigInteger">
  99                 ///     <para>
 100                 ///         The number to test.
 101                 ///     </para>
 102                 /// </param>
 103                 /// <param name="confidence" type="int">
 104                 ///     <para>
 105                 ///     The number of chosen bases. The test has at least a
 106                 ///     1/4^confidence chance of falsely returning True.
 107                 ///     </para>
 108                 /// </param>
 109                 /// <returns>
 110                 ///     <para>
 111                 ///             True if "this" is a strong pseudoprime to randomly chosen bases.
 112                 ///     </para>
 113                 ///     <para>
 114                 ///             False if "this" is definitely NOT prime.
 115                 ///     </para>
 116                 /// </returns>
 117                 public static bool RabinMillerTest (BigInteger bi, ConfidenceFactor confidence)
 118                 {
 119                         int Rounds = GetSPPRounds (bi, confidence);
 120
 121                         // calculate values of s and t
 122                         BigInteger p_sub1 = bi - 1;
 123                         int s = p_sub1.LowestSetBit ();
 124
 125                         BigInteger t = p_sub1 >> s;
 126
 127                         int bits = bi.BitCount ();
 128                         BigInteger a = null;
 129                         BigInteger.ModulusRing mr = new BigInteger.ModulusRing (bi);
 130
 131                         // Applying optimization from HAC section 4.50 (base == 2)
 132                         // not a really random base but an interesting (and speedy) one
 133                         BigInteger b = mr.Pow (2, t);
 134                         if (b != 1) {
 135                                 bool result = false;
 136                                 for (int j=0; j < s; j++) {
 137                                         if (b == p_sub1) {         // a^((2^j)*t) mod p = p-1 for some 0 <= j <= s-1
 138                                                 result = true;
 139                                                 break;
 140                                         }
 141
 142                                         b = (b * b) % bi;
 143                                 }
 144                                 if (!result)
 145                                         return false;
 146                         }
 147
 148                         // still here ? start at round 1 (round 0 was a == 2)
 149                         for (int round = 1; round < Rounds; round++) {
 150                                 while (true) {                     // generate a < n
 151                                         a = BigInteger.GenerateRandom (bits);
 152
 153                                         // make sure "a" is not 0 (and not 2 as we have already tested that)
 154                                         if (a > 2 && a < bi)
 155                                                 break;
 156                                 }
 157
 158                                 if (a.GCD (bi) != 1)
 159                                         return false;
 160
 161                                 b = mr.Pow (a, t);
 162
 163                                 if (b == 1)
 164                                         continue;              // a^t mod p = 1
 165
 166                                 bool result = false;
 167                                 for (int j = 0; j < s; j++) {
 168
 169                                         if (b == p_sub1) {         // a^((2^j)*t) mod p = p-1 for some 0 <= j <= s-1
 170                                                 result = true;
 171                                                 break;
 172                                         }
 173
 174                                         b = (b * b) % bi;
 175                                 }
 176
 177                                 if (!result)
 178                                         return false;
 179                         }
 180                         return true;
 181                 }
 182
 183                 public static bool SmallPrimeSppTest (BigInteger bi, ConfidenceFactor confidence)
 184                 {
 185                         int Rounds = GetSPPRounds (bi, confidence);
 186
 187                         // calculate values of s and t
 188                         BigInteger p_sub1 = bi - 1;
 189                         int s = p_sub1.LowestSetBit ();
 190
 191                         BigInteger t = p_sub1 >> s;
 192
 193
 194                         BigInteger.ModulusRing mr = new BigInteger.ModulusRing (bi);
 195
 196                         for (int round = 0; round < Rounds; round++) {
 197
 198                                 BigInteger b = mr.Pow (BigInteger.smallPrimes [round], t);
 199
 200                                 if (b == 1) continue;              // a^t mod p = 1
 201
 202                                 bool result = false;
 203                                 for (int j = 0; j < s; j++) {
 204
 205                                         if (b == p_sub1) {         // a^((2^j)*t) mod p = p-1 for some 0 <= j <= s-1
 206                                                 result = true;
 207                                                 break;
 208                                         }
 209
 210                                         b = (b * b) % bi;
 211                                 }
 212
 213                                 if (result == false)
 214                                         return false;
 215                         }
 216                         return true;
 217                 }
 218
 219                 #endregion
 220
 221                 // TODO: Implement the Lucus test
 222                 // TODO: Implement other new primality tests
 223                 // TODO: Implement primality proving
 224         }
 225 }