legacy-libs/grpc-cloned/deps/grpc/third_party/boringssl/crypto/fipsmodule/bn/prime.c

   1 /* Copyright (C) 1995-1998 Eric Young (eay@cryptsoft.com)
   2  * All rights reserved.
   3  *
   4  * This package is an SSL implementation written
   5  * by Eric Young (eay@cryptsoft.com).
   6  * The implementation was written so as to conform with Netscapes SSL.
   7  *
   8  * This library is free for commercial and non-commercial use as long as
   9  * the following conditions are aheared to.  The following conditions
  10  * apply to all code found in this distribution, be it the RC4, RSA,
  11  * lhash, DES, etc., code; not just the SSL code.  The SSL documentation
  12  * included with this distribution is covered by the same copyright terms
  13  * except that the holder is Tim Hudson (tjh@cryptsoft.com).
  14  *
  15  * Copyright remains Eric Young's, and as such any Copyright notices in
  16  * the code are not to be removed.
  17  * If this package is used in a product, Eric Young should be given attribution
  18  * as the author of the parts of the library used.
  19  * This can be in the form of a textual message at program startup or
  20  * in documentation (online or textual) provided with the package.
  21  *
  22  * Redistribution and use in source and binary forms, with or without
  23  * modification, are permitted provided that the following conditions
  24  * are met:
  25  * 1. Redistributions of source code must retain the copyright
  26  *    notice, this list of conditions and the following disclaimer.
  27  * 2. Redistributions in binary form must reproduce the above copyright
  28  *    notice, this list of conditions and the following disclaimer in the
  29  *    documentation and/or other materials provided with the distribution.
  30  * 3. All advertising materials mentioning features or use of this software
  31  *    must display the following acknowledgement:
  32  *    "This product includes cryptographic software written by
  33  *     Eric Young (eay@cryptsoft.com)"
  34  *    The word 'cryptographic' can be left out if the rouines from the library
  35  *    being used are not cryptographic related :-).
  36  * 4. If you include any Windows specific code (or a derivative thereof) from
  37  *    the apps directory (application code) you must include an acknowledgement:
  38  *    "This product includes software written by Tim Hudson (tjh@cryptsoft.com)"
  39  *
  40  * THIS SOFTWARE IS PROVIDED BY ERIC YOUNG ``AS IS'' AND
  41  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  42  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  43  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
  44  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
  45  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
  46  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
  47  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
  48  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
  49  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
  50  * SUCH DAMAGE.
  51  *
  52  * The licence and distribution terms for any publically available version or
  53  * derivative of this code cannot be changed.  i.e. this code cannot simply be
  54  * copied and put under another distribution licence
  55  * [including the GNU Public Licence.]
  56  */
  57 /* ====================================================================
  58  * Copyright (c) 1998-2001 The OpenSSL Project.  All rights reserved.
  59  *
  60  * Redistribution and use in source and binary forms, with or without
  61  * modification, are permitted provided that the following conditions
  62  * are met:
  63  *
  64  * 1. Redistributions of source code must retain the above copyright
  65  *    notice, this list of conditions and the following disclaimer.
  66  *
  67  * 2. Redistributions in binary form must reproduce the above copyright
  68  *    notice, this list of conditions and the following disclaimer in
  69  *    the documentation and/or other materials provided with the
  70  *    distribution.
  71  *
  72  * 3. All advertising materials mentioning features or use of this
  73  *    software must display the following acknowledgment:
  74  *    "This product includes software developed by the OpenSSL Project
  75  *    for use in the OpenSSL Toolkit. (http://www.openssl.org/)"
  76  *
  77  * 4. The names "OpenSSL Toolkit" and "OpenSSL Project" must not be used to
  78  *    endorse or promote products derived from this software without
  79  *    prior written permission. For written permission, please contact
  80  *    openssl-core@openssl.org.
  81  *
  82  * 5. Products derived from this software may not be called "OpenSSL"
  83  *    nor may "OpenSSL" appear in their names without prior written
  84  *    permission of the OpenSSL Project.
  85  *
  86  * 6. Redistributions of any form whatsoever must retain the following
  87  *    acknowledgment:
  88  *    "This product includes software developed by the OpenSSL Project
  89  *    for use in the OpenSSL Toolkit (http://www.openssl.org/)"
  90  *
  91  * THIS SOFTWARE IS PROVIDED BY THE OpenSSL PROJECT ``AS IS'' AND ANY
  92  * EXPRESSED OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  93  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
  94  * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE OpenSSL PROJECT OR
  95  * ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
  96  * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
  97  * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
  98  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
  99  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
 100  * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
 101  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
 102  * OF THE POSSIBILITY OF SUCH DAMAGE.
 103  * ====================================================================
 104  *
 105  * This product includes cryptographic software written by Eric Young
 106  * (eay@cryptsoft.com).  This product includes software written by Tim
 107  * Hudson (tjh@cryptsoft.com). */
 108
 109 #include <openssl/bn.h>
 110
 111 #include <openssl/err.h>
 112 #include <openssl/mem.h>
 113
 114 #include "internal.h"
 115 #include "../../internal.h"
 116
 117
 118 // The quick sieve algorithm approach to weeding out primes is Philip
 119 // Zimmermann's, as implemented in PGP.  I have had a read of his comments and
 120 // implemented my own version.
 121
 122 #define NUMPRIMES 2048
 123
 124 // primes contains all the primes that fit into a uint16_t.
 125 static const uint16_t primes[NUMPRIMES] = {
 126     2,     3,     5,     7,     11,    13,    17,    19,    23,    29,    31,
 127     37,    41,    43,    47,    53,    59,    61,    67,    71,    73,    79,
 128     83,    89,    97,    101,   103,   107,   109,   113,   127,   131,   137,
 129     139,   149,   151,   157,   163,   167,   173,   179,   181,   191,   193,
 130     197,   199,   211,   223,   227,   229,   233,   239,   241,   251,   257,
 131     263,   269,   271,   277,   281,   283,   293,   307,   311,   313,   317,
 132     331,   337,   347,   349,   353,   359,   367,   373,   379,   383,   389,
 133     397,   401,   409,   419,   421,   431,   433,   439,   443,   449,   457,
 134     461,   463,   467,   479,   487,   491,   499,   503,   509,   521,   523,
 135     541,   547,   557,   563,   569,   571,   577,   587,   593,   599,   601,
 136     607,   613,   617,   619,   631,   641,   643,   647,   653,   659,   661,
 137     673,   677,   683,   691,   701,   709,   719,   727,   733,   739,   743,
 138     751,   757,   761,   769,   773,   787,   797,   809,   811,   821,   823,
 139     827,   829,   839,   853,   857,   859,   863,   877,   881,   883,   887,
 140     907,   911,   919,   929,   937,   941,   947,   953,   967,   971,   977,
 141     983,   991,   997,   1009,  1013,  1019,  1021,  1031,  1033,  1039,  1049,
 142     1051,  1061,  1063,  1069,  1087,  1091,  1093,  1097,  1103,  1109,  1117,
 143     1123,  1129,  1151,  1153,  1163,  1171,  1181,  1187,  1193,  1201,  1213,
 144     1217,  1223,  1229,  1231,  1237,  1249,  1259,  1277,  1279,  1283,  1289,
 145     1291,  1297,  1301,  1303,  1307,  1319,  1321,  1327,  1361,  1367,  1373,
 146     1381,  1399,  1409,  1423,  1427,  1429,  1433,  1439,  1447,  1451,  1453,
 147     1459,  1471,  1481,  1483,  1487,  1489,  1493,  1499,  1511,  1523,  1531,
 148     1543,  1549,  1553,  1559,  1567,  1571,  1579,  1583,  1597,  1601,  1607,
 149     1609,  1613,  1619,  1621,  1627,  1637,  1657,  1663,  1667,  1669,  1693,
 150     1697,  1699,  1709,  1721,  1723,  1733,  1741,  1747,  1753,  1759,  1777,
 151     1783,  1787,  1789,  1801,  1811,  1823,  1831,  1847,  1861,  1867,  1871,
 152     1873,  1877,  1879,  1889,  1901,  1907,  1913,  1931,  1933,  1949,  1951,
 153     1973,  1979,  1987,  1993,  1997,  1999,  2003,  2011,  2017,  2027,  2029,
 154     2039,  2053,  2063,  2069,  2081,  2083,  2087,  2089,  2099,  2111,  2113,
 155     2129,  2131,  2137,  2141,  2143,  2153,  2161,  2179,  2203,  2207,  2213,
 156     2221,  2237,  2239,  2243,  2251,  2267,  2269,  2273,  2281,  2287,  2293,
 157     2297,  2309,  2311,  2333,  2339,  2341,  2347,  2351,  2357,  2371,  2377,
 158     2381,  2383,  2389,  2393,  2399,  2411,  2417,  2423,  2437,  2441,  2447,
 159     2459,  2467,  2473,  2477,  2503,  2521,  2531,  2539,  2543,  2549,  2551,
 160     2557,  2579,  2591,  2593,  2609,  2617,  2621,  2633,  2647,  2657,  2659,
 161     2663,  2671,  2677,  2683,  2687,  2689,  2693,  2699,  2707,  2711,  2713,
 162     2719,  2729,  2731,  2741,  2749,  2753,  2767,  2777,  2789,  2791,  2797,
 163     2801,  2803,  2819,  2833,  2837,  2843,  2851,  2857,  2861,  2879,  2887,
 164     2897,  2903,  2909,  2917,  2927,  2939,  2953,  2957,  2963,  2969,  2971,
 165     2999,  3001,  3011,  3019,  3023,  3037,  3041,  3049,  3061,  3067,  3079,
 166     3083,  3089,  3109,  3119,  3121,  3137,  3163,  3167,  3169,  3181,  3187,
 167     3191,  3203,  3209,  3217,  3221,  3229,  3251,  3253,  3257,  3259,  3271,
 168     3299,  3301,  3307,  3313,  3319,  3323,  3329,  3331,  3343,  3347,  3359,
 169     3361,  3371,  3373,  3389,  3391,  3407,  3413,  3433,  3449,  3457,  3461,
 170     3463,  3467,  3469,  3491,  3499,  3511,  3517,  3527,  3529,  3533,  3539,
 171     3541,  3547,  3557,  3559,  3571,  3581,  3583,  3593,  3607,  3613,  3617,
 172     3623,  3631,  3637,  3643,  3659,  3671,  3673,  3677,  3691,  3697,  3701,
 173     3709,  3719,  3727,  3733,  3739,  3761,  3767,  3769,  3779,  3793,  3797,
 174     3803,  3821,  3823,  3833,  3847,  3851,  3853,  3863,  3877,  3881,  3889,
 175     3907,  3911,  3917,  3919,  3923,  3929,  3931,  3943,  3947,  3967,  3989,
 176     4001,  4003,  4007,  4013,  4019,  4021,  4027,  4049,  4051,  4057,  4073,
 177     4079,  4091,  4093,  4099,  4111,  4127,  4129,  4133,  4139,  4153,  4157,
 178     4159,  4177,  4201,  4211,  4217,  4219,  4229,  4231,  4241,  4243,  4253,
 179     4259,  4261,  4271,  4273,  4283,  4289,  4297,  4327,  4337,  4339,  4349,
 180     4357,  4363,  4373,  4391,  4397,  4409,  4421,  4423,  4441,  4447,  4451,
 181     4457,  4463,  4481,  4483,  4493,  4507,  4513,  4517,  4519,  4523,  4547,
 182     4549,  4561,  4567,  4583,  4591,  4597,  4603,  4621,  4637,  4639,  4643,
 183     4649,  4651,  4657,  4663,  4673,  4679,  4691,  4703,  4721,  4723,  4729,
 184     4733,  4751,  4759,  4783,  4787,  4789,  4793,  4799,  4801,  4813,  4817,
 185     4831,  4861,  4871,  4877,  4889,  4903,  4909,  4919,  4931,  4933,  4937,
 186     4943,  4951,  4957,  4967,  4969,  4973,  4987,  4993,  4999,  5003,  5009,
 187     5011,  5021,  5023,  5039,  5051,  5059,  5077,  5081,  5087,  5099,  5101,
 188     5107,  5113,  5119,  5147,  5153,  5167,  5171,  5179,  5189,  5197,  5209,
 189     5227,  5231,  5233,  5237,  5261,  5273,  5279,  5281,  5297,  5303,  5309,
 190     5323,  5333,  5347,  5351,  5381,  5387,  5393,  5399,  5407,  5413,  5417,
 191     5419,  5431,  5437,  5441,  5443,  5449,  5471,  5477,  5479,  5483,  5501,
 192     5503,  5507,  5519,  5521,  5527,  5531,  5557,  5563,  5569,  5573,  5581,
 193     5591,  5623,  5639,  5641,  5647,  5651,  5653,  5657,  5659,  5669,  5683,
 194     5689,  5693,  5701,  5711,  5717,  5737,  5741,  5743,  5749,  5779,  5783,
 195     5791,  5801,  5807,  5813,  5821,  5827,  5839,  5843,  5849,  5851,  5857,
 196     5861,  5867,  5869,  5879,  5881,  5897,  5903,  5923,  5927,  5939,  5953,
 197     5981,  5987,  6007,  6011,  6029,  6037,  6043,  6047,  6053,  6067,  6073,
 198     6079,  6089,  6091,  6101,  6113,  6121,  6131,  6133,  6143,  6151,  6163,
 199     6173,  6197,  6199,  6203,  6211,  6217,  6221,  6229,  6247,  6257,  6263,
 200     6269,  6271,  6277,  6287,  6299,  6301,  6311,  6317,  6323,  6329,  6337,
 201     6343,  6353,  6359,  6361,  6367,  6373,  6379,  6389,  6397,  6421,  6427,
 202     6449,  6451,  6469,  6473,  6481,  6491,  6521,  6529,  6547,  6551,  6553,
 203     6563,  6569,  6571,  6577,  6581,  6599,  6607,  6619,  6637,  6653,  6659,
 204     6661,  6673,  6679,  6689,  6691,  6701,  6703,  6709,  6719,  6733,  6737,
 205     6761,  6763,  6779,  6781,  6791,  6793,  6803,  6823,  6827,  6829,  6833,
 206     6841,  6857,  6863,  6869,  6871,  6883,  6899,  6907,  6911,  6917,  6947,
 207     6949,  6959,  6961,  6967,  6971,  6977,  6983,  6991,  6997,  7001,  7013,
 208     7019,  7027,  7039,  7043,  7057,  7069,  7079,  7103,  7109,  7121,  7127,
 209     7129,  7151,  7159,  7177,  7187,  7193,  7207,  7211,  7213,  7219,  7229,
 210     7237,  7243,  7247,  7253,  7283,  7297,  7307,  7309,  7321,  7331,  7333,
 211     7349,  7351,  7369,  7393,  7411,  7417,  7433,  7451,  7457,  7459,  7477,
 212     7481,  7487,  7489,  7499,  7507,  7517,  7523,  7529,  7537,  7541,  7547,
 213     7549,  7559,  7561,  7573,  7577,  7583,  7589,  7591,  7603,  7607,  7621,
 214     7639,  7643,  7649,  7669,  7673,  7681,  7687,  7691,  7699,  7703,  7717,
 215     7723,  7727,  7741,  7753,  7757,  7759,  7789,  7793,  7817,  7823,  7829,
 216     7841,  7853,  7867,  7873,  7877,  7879,  7883,  7901,  7907,  7919,  7927,
 217     7933,  7937,  7949,  7951,  7963,  7993,  8009,  8011,  8017,  8039,  8053,
 218     8059,  8069,  8081,  8087,  8089,  8093,  8101,  8111,  8117,  8123,  8147,
 219     8161,  8167,  8171,  8179,  8191,  8209,  8219,  8221,  8231,  8233,  8237,
 220     8243,  8263,  8269,  8273,  8287,  8291,  8293,  8297,  8311,  8317,  8329,
 221     8353,  8363,  8369,  8377,  8387,  8389,  8419,  8423,  8429,  8431,  8443,
 222     8447,  8461,  8467,  8501,  8513,  8521,  8527,  8537,  8539,  8543,  8563,
 223     8573,  8581,  8597,  8599,  8609,  8623,  8627,  8629,  8641,  8647,  8663,
 224     8669,  8677,  8681,  8689,  8693,  8699,  8707,  8713,  8719,  8731,  8737,
 225     8741,  8747,  8753,  8761,  8779,  8783,  8803,  8807,  8819,  8821,  8831,
 226     8837,  8839,  8849,  8861,  8863,  8867,  8887,  8893,  8923,  8929,  8933,
 227     8941,  8951,  8963,  8969,  8971,  8999,  9001,  9007,  9011,  9013,  9029,
 228     9041,  9043,  9049,  9059,  9067,  9091,  9103,  9109,  9127,  9133,  9137,
 229     9151,  9157,  9161,  9173,  9181,  9187,  9199,  9203,  9209,  9221,  9227,
 230     9239,  9241,  9257,  9277,  9281,  9283,  9293,  9311,  9319,  9323,  9337,
 231     9341,  9343,  9349,  9371,  9377,  9391,  9397,  9403,  9413,  9419,  9421,
 232     9431,  9433,  9437,  9439,  9461,  9463,  9467,  9473,  9479,  9491,  9497,
 233     9511,  9521,  9533,  9539,  9547,  9551,  9587,  9601,  9613,  9619,  9623,
 234     9629,  9631,  9643,  9649,  9661,  9677,  9679,  9689,  9697,  9719,  9721,
 235     9733,  9739,  9743,  9749,  9767,  9769,  9781,  9787,  9791,  9803,  9811,
 236     9817,  9829,  9833,  9839,  9851,  9857,  9859,  9871,  9883,  9887,  9901,
 237     9907,  9923,  9929,  9931,  9941,  9949,  9967,  9973,  10007, 10009, 10037,
 238     10039, 10061, 10067, 10069, 10079, 10091, 10093, 10099, 10103, 10111, 10133,
 239     10139, 10141, 10151, 10159, 10163, 10169, 10177, 10181, 10193, 10211, 10223,
 240     10243, 10247, 10253, 10259, 10267, 10271, 10273, 10289, 10301, 10303, 10313,
 241     10321, 10331, 10333, 10337, 10343, 10357, 10369, 10391, 10399, 10427, 10429,
 242     10433, 10453, 10457, 10459, 10463, 10477, 10487, 10499, 10501, 10513, 10529,
 243     10531, 10559, 10567, 10589, 10597, 10601, 10607, 10613, 10627, 10631, 10639,
 244     10651, 10657, 10663, 10667, 10687, 10691, 10709, 10711, 10723, 10729, 10733,
 245     10739, 10753, 10771, 10781, 10789, 10799, 10831, 10837, 10847, 10853, 10859,
 246     10861, 10867, 10883, 10889, 10891, 10903, 10909, 10937, 10939, 10949, 10957,
 247     10973, 10979, 10987, 10993, 11003, 11027, 11047, 11057, 11059, 11069, 11071,
 248     11083, 11087, 11093, 11113, 11117, 11119, 11131, 11149, 11159, 11161, 11171,
 249     11173, 11177, 11197, 11213, 11239, 11243, 11251, 11257, 11261, 11273, 11279,
 250     11287, 11299, 11311, 11317, 11321, 11329, 11351, 11353, 11369, 11383, 11393,
 251     11399, 11411, 11423, 11437, 11443, 11447, 11467, 11471, 11483, 11489, 11491,
 252     11497, 11503, 11519, 11527, 11549, 11551, 11579, 11587, 11593, 11597, 11617,
 253     11621, 11633, 11657, 11677, 11681, 11689, 11699, 11701, 11717, 11719, 11731,
 254     11743, 11777, 11779, 11783, 11789, 11801, 11807, 11813, 11821, 11827, 11831,
 255     11833, 11839, 11863, 11867, 11887, 11897, 11903, 11909, 11923, 11927, 11933,
 256     11939, 11941, 11953, 11959, 11969, 11971, 11981, 11987, 12007, 12011, 12037,
 257     12041, 12043, 12049, 12071, 12073, 12097, 12101, 12107, 12109, 12113, 12119,
 258     12143, 12149, 12157, 12161, 12163, 12197, 12203, 12211, 12227, 12239, 12241,
 259     12251, 12253, 12263, 12269, 12277, 12281, 12289, 12301, 12323, 12329, 12343,
 260     12347, 12373, 12377, 12379, 12391, 12401, 12409, 12413, 12421, 12433, 12437,
 261     12451, 12457, 12473, 12479, 12487, 12491, 12497, 12503, 12511, 12517, 12527,
 262     12539, 12541, 12547, 12553, 12569, 12577, 12583, 12589, 12601, 12611, 12613,
 263     12619, 12637, 12641, 12647, 12653, 12659, 12671, 12689, 12697, 12703, 12713,
 264     12721, 12739, 12743, 12757, 12763, 12781, 12791, 12799, 12809, 12821, 12823,
 265     12829, 12841, 12853, 12889, 12893, 12899, 12907, 12911, 12917, 12919, 12923,
 266     12941, 12953, 12959, 12967, 12973, 12979, 12983, 13001, 13003, 13007, 13009,
 267     13033, 13037, 13043, 13049, 13063, 13093, 13099, 13103, 13109, 13121, 13127,
 268     13147, 13151, 13159, 13163, 13171, 13177, 13183, 13187, 13217, 13219, 13229,
 269     13241, 13249, 13259, 13267, 13291, 13297, 13309, 13313, 13327, 13331, 13337,
 270     13339, 13367, 13381, 13397, 13399, 13411, 13417, 13421, 13441, 13451, 13457,
 271     13463, 13469, 13477, 13487, 13499, 13513, 13523, 13537, 13553, 13567, 13577,
 272     13591, 13597, 13613, 13619, 13627, 13633, 13649, 13669, 13679, 13681, 13687,
 273     13691, 13693, 13697, 13709, 13711, 13721, 13723, 13729, 13751, 13757, 13759,
 274     13763, 13781, 13789, 13799, 13807, 13829, 13831, 13841, 13859, 13873, 13877,
 275     13879, 13883, 13901, 13903, 13907, 13913, 13921, 13931, 13933, 13963, 13967,
 276     13997, 13999, 14009, 14011, 14029, 14033, 14051, 14057, 14071, 14081, 14083,
 277     14087, 14107, 14143, 14149, 14153, 14159, 14173, 14177, 14197, 14207, 14221,
 278     14243, 14249, 14251, 14281, 14293, 14303, 14321, 14323, 14327, 14341, 14347,
 279     14369, 14387, 14389, 14401, 14407, 14411, 14419, 14423, 14431, 14437, 14447,
 280     14449, 14461, 14479, 14489, 14503, 14519, 14533, 14537, 14543, 14549, 14551,
 281     14557, 14561, 14563, 14591, 14593, 14621, 14627, 14629, 14633, 14639, 14653,
 282     14657, 14669, 14683, 14699, 14713, 14717, 14723, 14731, 14737, 14741, 14747,
 283     14753, 14759, 14767, 14771, 14779, 14783, 14797, 14813, 14821, 14827, 14831,
 284     14843, 14851, 14867, 14869, 14879, 14887, 14891, 14897, 14923, 14929, 14939,
 285     14947, 14951, 14957, 14969, 14983, 15013, 15017, 15031, 15053, 15061, 15073,
 286     15077, 15083, 15091, 15101, 15107, 15121, 15131, 15137, 15139, 15149, 15161,
 287     15173, 15187, 15193, 15199, 15217, 15227, 15233, 15241, 15259, 15263, 15269,
 288     15271, 15277, 15287, 15289, 15299, 15307, 15313, 15319, 15329, 15331, 15349,
 289     15359, 15361, 15373, 15377, 15383, 15391, 15401, 15413, 15427, 15439, 15443,
 290     15451, 15461, 15467, 15473, 15493, 15497, 15511, 15527, 15541, 15551, 15559,
 291     15569, 15581, 15583, 15601, 15607, 15619, 15629, 15641, 15643, 15647, 15649,
 292     15661, 15667, 15671, 15679, 15683, 15727, 15731, 15733, 15737, 15739, 15749,
 293     15761, 15767, 15773, 15787, 15791, 15797, 15803, 15809, 15817, 15823, 15859,
 294     15877, 15881, 15887, 15889, 15901, 15907, 15913, 15919, 15923, 15937, 15959,
 295     15971, 15973, 15991, 16001, 16007, 16033, 16057, 16061, 16063, 16067, 16069,
 296     16073, 16087, 16091, 16097, 16103, 16111, 16127, 16139, 16141, 16183, 16187,
 297     16189, 16193, 16217, 16223, 16229, 16231, 16249, 16253, 16267, 16273, 16301,
 298     16319, 16333, 16339, 16349, 16361, 16363, 16369, 16381, 16411, 16417, 16421,
 299     16427, 16433, 16447, 16451, 16453, 16477, 16481, 16487, 16493, 16519, 16529,
 300     16547, 16553, 16561, 16567, 16573, 16603, 16607, 16619, 16631, 16633, 16649,
 301     16651, 16657, 16661, 16673, 16691, 16693, 16699, 16703, 16729, 16741, 16747,
 302     16759, 16763, 16787, 16811, 16823, 16829, 16831, 16843, 16871, 16879, 16883,
 303     16889, 16901, 16903, 16921, 16927, 16931, 16937, 16943, 16963, 16979, 16981,
 304     16987, 16993, 17011, 17021, 17027, 17029, 17033, 17041, 17047, 17053, 17077,
 305     17093, 17099, 17107, 17117, 17123, 17137, 17159, 17167, 17183, 17189, 17191,
 306     17203, 17207, 17209, 17231, 17239, 17257, 17291, 17293, 17299, 17317, 17321,
 307     17327, 17333, 17341, 17351, 17359, 17377, 17383, 17387, 17389, 17393, 17401,
 308     17417, 17419, 17431, 17443, 17449, 17467, 17471, 17477, 17483, 17489, 17491,
 309     17497, 17509, 17519, 17539, 17551, 17569, 17573, 17579, 17581, 17597, 17599,
 310     17609, 17623, 17627, 17657, 17659, 17669, 17681, 17683, 17707, 17713, 17729,
 311     17737, 17747, 17749, 17761, 17783, 17789, 17791, 17807, 17827, 17837, 17839,
 312     17851, 17863,
 313 };
 314
 315 // BN_prime_checks_for_size returns the number of Miller-Rabin iterations
 316 // necessary for a 'bits'-bit prime, in order to maintain an error rate greater
 317 // than the security level for an RSA prime of that many bits (calculated using
 318 // the FIPS SP 800-57 security level and 186-4 Section F.1; original paper:
 319 // Damgaard, Landrock, Pomerance: Average case error estimates for the strong
 320 // probable prime test. -- Math. Comp. 61 (1993) 177-194)
 321 static int BN_prime_checks_for_size(int bits) {
 322   if (bits >= 3747) {
 323     return 3;
 324   }
 325   if (bits >= 1345) {
 326     return 4;
 327   }
 328   if (bits >= 476) {
 329     return 5;
 330   }
 331   if (bits >= 400) {
 332     return 6;
 333   }
 334   if (bits >= 308) {
 335     return 8;
 336   }
 337   if (bits >= 205) {
 338     return 13;
 339   }
 340   if (bits >= 155) {
 341     return 19;
 342   }
 343   return 28;
 344 }
 345
 346 // BN_PRIME_CHECKS_BLINDED is the iteration count for blinding the constant-time
 347 // primality test. See |BN_primality_test| for details. This number is selected
 348 // so that, for a candidate N-bit RSA prime, picking |BN_PRIME_CHECKS_BLINDED|
 349 // random N-bit numbers will have at least |BN_prime_checks_for_size(N)| values
 350 // in range with high probability.
 351 //
 352 // The following Python script computes the blinding factor needed for the
 353 // corresponding iteration count.
 354 /*
 355 import math
 356
 357 # We choose candidate RSA primes between sqrt(2)/2 * 2^N and 2^N and select
 358 # witnesses by generating random N-bit numbers. Thus the probability of
 359 # selecting one in range is at least sqrt(2)/2.
 360 p = math.sqrt(2) / 2
 361
 362 # Target around 2^-8 probability of the blinding being insufficient given that
 363 # key generation is a one-time, noisy operation.
 364 epsilon = 2**-8
 365
 366 def choose(a, b):
 367   r = 1
 368   for i in xrange(b):
 369     r *= a - i
 370     r /= (i + 1)
 371   return r
 372
 373 def failure_rate(min_uniform, iterations):
 374   """ Returns the probability that, for |iterations| candidate witnesses, fewer
 375       than |min_uniform| of them will be uniform. """
 376   prob = 0.0
 377   for i in xrange(min_uniform):
 378     prob += (choose(iterations, i) *
 379              p**i * (1-p)**(iterations - i))
 380   return prob
 381
 382 for min_uniform in (3, 4, 5, 6, 8, 13, 19, 28):
 383   # Find the smallest number of iterations under the target failure rate.
 384   iterations = min_uniform
 385   while True:
 386     prob = failure_rate(min_uniform, iterations)
 387     if prob < epsilon:
 388       print min_uniform, iterations, prob
 389       break
 390     iterations += 1
 391
 392 Output:
 393   3 9 0.00368894873911
 394   4 11 0.00363319494662
 395   5 13 0.00336215573898
 396   6 15 0.00300145783158
 397   8 19 0.00225214119331
 398   13 27 0.00385610026955
 399   19 38 0.0021410539126
 400   28 52 0.00325405801769
 401
 402 16 iterations suffices for 400-bit primes and larger (6 uniform samples needed),
 403 which is already well below the minimum acceptable key size for RSA.
 404 */
 405 #define BN_PRIME_CHECKS_BLINDED 16
 406
 407 static int probable_prime(BIGNUM *rnd, int bits);
 408 static int probable_prime_dh(BIGNUM *rnd, int bits, const BIGNUM *add,
 409                              const BIGNUM *rem, BN_CTX *ctx);
 410 static int probable_prime_dh_safe(BIGNUM *rnd, int bits, const BIGNUM *add,
 411                                   const BIGNUM *rem, BN_CTX *ctx);
 412
 413 void BN_GENCB_set(BN_GENCB *callback,
 414                   int (*f)(int event, int n, struct bn_gencb_st *),
 415                   void *arg) {
 416   callback->callback = f;
 417   callback->arg = arg;
 418 }
 419
 420 int BN_GENCB_call(BN_GENCB *callback, int event, int n) {
 421   if (!callback) {
 422     return 1;
 423   }
 424
 425   return callback->callback(event, n, callback);
 426 }
 427
 428 int BN_generate_prime_ex(BIGNUM *ret, int bits, int safe, const BIGNUM *add,
 429                          const BIGNUM *rem, BN_GENCB *cb) {
 430   BIGNUM *t;
 431   int found = 0;
 432   int i, j, c1 = 0;
 433   BN_CTX *ctx;
 434   int checks = BN_prime_checks_for_size(bits);
 435
 436   if (bits < 2) {
 437     // There are no prime numbers this small.
 438     OPENSSL_PUT_ERROR(BN, BN_R_BITS_TOO_SMALL);
 439     return 0;
 440   } else if (bits == 2 && safe) {
 441     // The smallest safe prime (7) is three bits.
 442     OPENSSL_PUT_ERROR(BN, BN_R_BITS_TOO_SMALL);
 443     return 0;
 444   }
 445
 446   ctx = BN_CTX_new();
 447   if (ctx == NULL) {
 448     goto err;
 449   }
 450   BN_CTX_start(ctx);
 451   t = BN_CTX_get(ctx);
 452   if (!t) {
 453     goto err;
 454   }
 455
 456 loop:
 457   // make a random number and set the top and bottom bits
 458   if (add == NULL) {
 459     if (!probable_prime(ret, bits)) {
 460       goto err;
 461     }
 462   } else {
 463     if (safe) {
 464       if (!probable_prime_dh_safe(ret, bits, add, rem, ctx)) {
 465         goto err;
 466       }
 467     } else {
 468       if (!probable_prime_dh(ret, bits, add, rem, ctx)) {
 469         goto err;
 470       }
 471     }
 472   }
 473
 474   if (!BN_GENCB_call(cb, BN_GENCB_GENERATED, c1++)) {
 475     // aborted
 476     goto err;
 477   }
 478
 479   if (!safe) {
 480     i = BN_is_prime_fasttest_ex(ret, checks, ctx, 0, cb);
 481     if (i == -1) {
 482       goto err;
 483     } else if (i == 0) {
 484       goto loop;
 485     }
 486   } else {
 487     // for "safe prime" generation, check that (p-1)/2 is prime. Since a prime
 488     // is odd, We just need to divide by 2
 489     if (!BN_rshift1(t, ret)) {
 490       goto err;
 491     }
 492
 493     for (i = 0; i < checks; i++) {
 494       j = BN_is_prime_fasttest_ex(ret, 1, ctx, 0, NULL);
 495       if (j == -1) {
 496         goto err;
 497       } else if (j == 0) {
 498         goto loop;
 499       }
 500
 501       j = BN_is_prime_fasttest_ex(t, 1, ctx, 0, NULL);
 502       if (j == -1) {
 503         goto err;
 504       } else if (j == 0) {
 505         goto loop;
 506       }
 507
 508       if (!BN_GENCB_call(cb, i, c1 - 1)) {
 509         goto err;
 510       }
 511       // We have a safe prime test pass
 512     }
 513   }
 514
 515   // we have a prime :-)
 516   found = 1;
 517
 518 err:
 519   if (ctx != NULL) {
 520     BN_CTX_end(ctx);
 521     BN_CTX_free(ctx);
 522   }
 523
 524   return found;
 525 }
 526
 527 // The following functions use a Barrett reduction variant to avoid leaking the
 528 // numerator. See http://ridiculousfish.com/blog/posts/labor-of-division-episode-i.html
 529 //
 530 // We use 32-bit numerator and 16-bit divisor for simplicity. This allows
 531 // computing |m| and |q| without architecture-specific code.
 532
 533 // mod_u16 returns |n| mod |d|. |p| and |m| are the "magic numbers" for |d| (see
 534 // reference). For proof of correctness in Coq, see
 535 // https://github.com/davidben/fiat-crypto/blob/barrett/src/Arithmetic/BarrettReduction/RidiculousFish.v
 536 // Note the Coq version of |mod_u16| additionally includes the computation of
 537 // |p| and |m| from |bn_mod_u16_consttime| below.
 538 static uint16_t mod_u16(uint32_t n, uint16_t d, uint32_t p, uint32_t m) {
 539   // Compute floor(n/d) per steps 3 through 5.
 540   uint32_t q = ((uint64_t)m * n) >> 32;
 541   // Note there is a typo in the reference. We right-shift by one, not two.
 542   uint32_t t = ((n - q) >> 1) + q;
 543   t = t >> (p - 1);
 544
 545   // Multiply and subtract to get the remainder.
 546   n -= d * t;
 547   assert(n < d);
 548   return n;
 549 }
 550
 551 // shift_and_add_mod_u16 returns |r| * 2^32 + |a| mod |d|. |p| and |m| are the
 552 // "magic numbers" for |d| (see reference).
 553 static uint16_t shift_and_add_mod_u16(uint16_t r, uint32_t a, uint16_t d,
 554                                       uint32_t p, uint32_t m) {
 555   // Incorporate |a| in two 16-bit chunks.
 556   uint32_t t = r;
 557   t <<= 16;
 558   t |= a >> 16;
 559   t = mod_u16(t, d, p, m);
 560
 561   t <<= 16;
 562   t |= a & 0xffff;
 563   t = mod_u16(t, d, p, m);
 564   return t;
 565 }
 566
 567 uint16_t bn_mod_u16_consttime(const BIGNUM *bn, uint16_t d) {
 568   if (d <= 1) {
 569     return 0;
 570   }
 571
 572   // Compute the "magic numbers" for |d|. See steps 1 and 2.
 573   // This computes p = ceil(log_2(d)).
 574   uint32_t p = BN_num_bits_word(d - 1);
 575   // This operation is not constant-time, but |p| and |d| are public values.
 576   // Note that |p| is at most 16, so the computation fits in |uint64_t|.
 577   assert(p <= 16);
 578   uint32_t m = ((UINT64_C(1) << (32 + p)) + d - 1) / d;
 579
 580   uint16_t ret = 0;
 581   for (int i = bn->width - 1; i >= 0; i--) {
 582 #if BN_BITS2 == 32
 583     ret = shift_and_add_mod_u16(ret, bn->d[i], d, p, m);
 584 #elif BN_BITS2 == 64
 585     ret = shift_and_add_mod_u16(ret, bn->d[i] >> 32, d, p, m);
 586     ret = shift_and_add_mod_u16(ret, bn->d[i] & 0xffffffff, d, p, m);
 587 #else
 588 #error "Unknown BN_ULONG size"
 589 #endif
 590   }
 591   return ret;
 592 }
 593
 594 static int bn_trial_division(uint16_t *out, const BIGNUM *bn) {
 595   for (int i = 1; i < NUMPRIMES; i++) {
 596     if (bn_mod_u16_consttime(bn, primes[i]) == 0) {
 597       *out = primes[i];
 598       return 1;
 599     }
 600   }
 601   return 0;
 602 }
 603
 604 int bn_odd_number_is_obviously_composite(const BIGNUM *bn) {
 605   uint16_t prime;
 606   return bn_trial_division(&prime, bn) && !BN_is_word(bn, prime);
 607 }
 608
 609 int BN_primality_test(int *is_probably_prime, const BIGNUM *w,
 610                       int iterations, BN_CTX *ctx, int do_trial_division,
 611                       BN_GENCB *cb) {
 612   *is_probably_prime = 0;
 613
 614   // To support RSA key generation, this function should treat |w| as secret if
 615   // it is a large prime. Composite numbers are discarded, so they may return
 616   // early.
 617
 618   if (BN_cmp(w, BN_value_one()) <= 0) {
 619     return 1;
 620   }
 621
 622   if (!BN_is_odd(w)) {
 623     // The only even prime is two.
 624     *is_probably_prime = BN_is_word(w, 2);
 625     return 1;
 626   }
 627
 628   // Miller-Rabin does not work for three.
 629   if (BN_is_word(w, 3)) {
 630     *is_probably_prime = 1;
 631     return 1;
 632   }
 633
 634   if (do_trial_division) {
 635     // Perform additional trial division checks to discard small primes.
 636     uint16_t prime;
 637     if (bn_trial_division(&prime, w)) {
 638       *is_probably_prime = BN_is_word(w, prime);
 639       return 1;
 640     }
 641     if (!BN_GENCB_call(cb, 1, -1)) {
 642       return 0;
 643     }
 644   }
 645
 646   if (iterations == BN_prime_checks) {
 647     iterations = BN_prime_checks_for_size(BN_num_bits(w));
 648   }
 649
 650   // See C.3.1 from FIPS 186-4.
 651   int ret = 0;
 652   BN_MONT_CTX *mont = NULL;
 653   BN_CTX_start(ctx);
 654   BIGNUM *w1 = BN_CTX_get(ctx);
 655   if (w1 == NULL ||
 656       !bn_usub_consttime(w1, w, BN_value_one())) {
 657     goto err;
 658   }
 659
 660   // Write w1 as m * 2^a (Steps 1 and 2).
 661   int w_len = BN_num_bits(w);
 662   int a = BN_count_low_zero_bits(w1);
 663   BIGNUM *m = BN_CTX_get(ctx);
 664   if (m == NULL ||
 665       !bn_rshift_secret_shift(m, w1, a, ctx)) {
 666     goto err;
 667   }
 668
 669   // Montgomery setup for computations mod w. Additionally, compute 1 and w - 1
 670   // in the Montgomery domain for later comparisons.
 671   BIGNUM *b = BN_CTX_get(ctx);
 672   BIGNUM *z = BN_CTX_get(ctx);
 673   BIGNUM *one_mont = BN_CTX_get(ctx);
 674   BIGNUM *w1_mont = BN_CTX_get(ctx);
 675   mont = BN_MONT_CTX_new_for_modulus(w, ctx);
 676   if (b == NULL || z == NULL || one_mont == NULL || w1_mont == NULL ||
 677       mont == NULL ||
 678       !bn_one_to_montgomery(one_mont, mont, ctx) ||
 679       // w - 1 is -1 mod w, so we can compute it in the Montgomery domain, -R,
 680       // with a subtraction. (|one_mont| cannot be zero.)
 681       !bn_usub_consttime(w1_mont, w, one_mont)) {
 682     goto err;
 683   }
 684
 685   // The following loop performs in inner iteration of the Miller-Rabin
 686   // Primality test (Step 4).
 687   //
 688   // The algorithm as specified in FIPS 186-4 leaks information on |w|, the RSA
 689   // private key. Instead, we run through each iteration unconditionally,
 690   // performing modular multiplications, masking off any effects to behave
 691   // equivalently to the specified algorithm.
 692   //
 693   // We also blind the number of values of |b| we try. Steps 4.1–4.2 say to
 694   // discard out-of-range values. To avoid leaking information on |w|, we use
 695   // |bn_rand_secret_range| which, rather than discarding bad values, adjusts
 696   // them to be in range. Though not uniformly selected, these adjusted values
 697   // are still usable as Rabin-Miller checks.
 698   //
 699   // Rabin-Miller is already probabilistic, so we could reach the desired
 700   // confidence levels by just suitably increasing the iteration count. However,
 701   // to align with FIPS 186-4, we use a more pessimal analysis: we do not count
 702   // the non-uniform values towards the iteration count. As a result, this
 703   // function is more complex and has more timing risk than necessary.
 704   //
 705   // We count both total iterations and uniform ones and iterate until we've
 706   // reached at least |BN_PRIME_CHECKS_BLINDED| and |iterations|, respectively.
 707   // If the latter is large enough, it will be the limiting factor with high
 708   // probability and we won't leak information.
 709   //
 710   // Note this blinding does not impact most calls when picking primes because
 711   // composites are rejected early. Only the two secret primes see extra work.
 712
 713   crypto_word_t uniform_iterations = 0;
 714   // Using |constant_time_lt_w| seems to prevent the compiler from optimizing
 715   // this into two jumps.
 716   for (int i = 1; (i <= BN_PRIME_CHECKS_BLINDED) |
 717                   constant_time_lt_w(uniform_iterations, iterations);
 718        i++) {
 719     int is_uniform;
 720     if (// Step 4.1-4.2
 721         !bn_rand_secret_range(b, &is_uniform, 2, w1) ||
 722         // Step 4.3
 723         !BN_mod_exp_mont_consttime(z, b, m, w, ctx, mont)) {
 724       goto err;
 725     }
 726     uniform_iterations += is_uniform;
 727
 728     // loop_done is all ones if the loop has completed and all zeros otherwise.
 729     crypto_word_t loop_done = 0;
 730     // next_iteration is all ones if we should continue to the next iteration
 731     // (|b| is not a composite witness for |w|). This is equivalent to going to
 732     // step 4.7 in the original algorithm.
 733     crypto_word_t next_iteration = 0;
 734
 735     // Step 4.4. If z = 1 or z = w-1, mask off the loop and continue to the next
 736     // iteration (go to step 4.7).
 737     loop_done = BN_equal_consttime(z, BN_value_one()) |
 738                 BN_equal_consttime(z, w1);
 739     loop_done = 0 - loop_done;   // Make it all zeros or all ones.
 740     next_iteration = loop_done;  // Go to step 4.7 if |loop_done|.
 741
 742     // Step 4.5. We use Montgomery-encoding for better performance and to avoid
 743     // timing leaks.
 744     if (!BN_to_montgomery(z, z, mont, ctx)) {
 745       goto err;
 746     }
 747
 748     // To avoid leaking |a|, we run the loop to |w_len| and mask off all
 749     // iterations once |j| = |a|.
 750     for (int j = 1; j < w_len; j++) {
 751       loop_done |= constant_time_eq_int(j, a);
 752
 753       // Step 4.5.1.
 754       if (!BN_mod_mul_montgomery(z, z, z, mont, ctx)) {
 755         goto err;
 756       }
 757
 758       // Step 4.5.2. If z = w-1 and the loop is not done, run through the next
 759       // iteration.
 760       crypto_word_t z_is_w1_mont = BN_equal_consttime(z, w1_mont) & ~loop_done;
 761       z_is_w1_mont = 0 - z_is_w1_mont;  // Make it all zeros or all ones.
 762       loop_done |= z_is_w1_mont;
 763       next_iteration |= z_is_w1_mont;  // Go to step 4.7 if |z_is_w1_mont|.
 764
 765       // Step 4.5.3. If z = 1 and the loop is not done, w is composite and we
 766       // may exit in variable time.
 767       if (BN_equal_consttime(z, one_mont) & ~loop_done) {
 768         assert(!next_iteration);
 769         break;
 770       }
 771     }
 772
 773     if (!next_iteration) {
 774       // Step 4.6. We did not see z = w-1 before z = 1, so w must be composite.
 775       // (For any prime, the value of z immediately preceding 1 must be -1.
 776       // There are no non-trivial square roots of 1 modulo a prime.)
 777       *is_probably_prime = 0;
 778       ret = 1;
 779       goto err;
 780     }
 781
 782     // Step 4.7
 783     if (!BN_GENCB_call(cb, 1, i)) {
 784       goto err;
 785     }
 786   }
 787
 788   assert(uniform_iterations >= (crypto_word_t)iterations);
 789   *is_probably_prime = 1;
 790   ret = 1;
 791
 792 err:
 793   BN_MONT_CTX_free(mont);
 794   BN_CTX_end(ctx);
 795   return ret;
 796 }
 797
 798 int BN_is_prime_ex(const BIGNUM *candidate, int checks, BN_CTX *ctx, BN_GENCB *cb) {
 799   return BN_is_prime_fasttest_ex(candidate, checks, ctx, 0, cb);
 800 }
 801
 802 int BN_is_prime_fasttest_ex(const BIGNUM *a, int checks, BN_CTX *ctx,
 803                             int do_trial_division, BN_GENCB *cb) {
 804   int is_probably_prime;
 805   if (!BN_primality_test(&is_probably_prime, a, checks, ctx, do_trial_division,
 806                          cb)) {
 807     return -1;
 808   }
 809   return is_probably_prime;
 810 }
 811
 812 int BN_enhanced_miller_rabin_primality_test(
 813     enum bn_primality_result_t *out_result, const BIGNUM *w, int iterations,
 814     BN_CTX *ctx, BN_GENCB *cb) {
 815   // Enhanced Miller-Rabin is only valid on odd integers greater than 3.
 816   if (!BN_is_odd(w) || BN_cmp_word(w, 3) <= 0) {
 817     OPENSSL_PUT_ERROR(BN, BN_R_INVALID_INPUT);
 818     return 0;
 819   }
 820
 821   if (iterations == BN_prime_checks) {
 822     iterations = BN_prime_checks_for_size(BN_num_bits(w));
 823   }
 824
 825   int ret = 0;
 826   BN_MONT_CTX *mont = NULL;
 827
 828   BN_CTX_start(ctx);
 829
 830   BIGNUM *w1 = BN_CTX_get(ctx);
 831   if (w1 == NULL ||
 832       !BN_copy(w1, w) ||
 833       !BN_sub_word(w1, 1)) {
 834     goto err;
 835   }
 836
 837   // Write w1 as m*2^a (Steps 1 and 2).
 838   int a = 0;
 839   while (!BN_is_bit_set(w1, a)) {
 840     a++;
 841   }
 842   BIGNUM *m = BN_CTX_get(ctx);
 843   if (m == NULL ||
 844       !BN_rshift(m, w1, a)) {
 845     goto err;
 846   }
 847
 848   BIGNUM *b = BN_CTX_get(ctx);
 849   BIGNUM *g = BN_CTX_get(ctx);
 850   BIGNUM *z = BN_CTX_get(ctx);
 851   BIGNUM *x = BN_CTX_get(ctx);
 852   BIGNUM *x1 = BN_CTX_get(ctx);
 853   if (b == NULL ||
 854       g == NULL ||
 855       z == NULL ||
 856       x == NULL ||
 857       x1 == NULL) {
 858     goto err;
 859   }
 860
 861   // Montgomery setup for computations mod w
 862   mont = BN_MONT_CTX_new_for_modulus(w, ctx);
 863   if (mont == NULL) {
 864     goto err;
 865   }
 866
 867   // The following loop performs in inner iteration of the Enhanced Miller-Rabin
 868   // Primality test (Step 4).
 869   for (int i = 1; i <= iterations; i++) {
 870     // Step 4.1-4.2
 871     if (!BN_rand_range_ex(b, 2, w1)) {
 872       goto err;
 873     }
 874
 875     // Step 4.3-4.4
 876     if (!BN_gcd(g, b, w, ctx)) {
 877       goto err;
 878     }
 879     if (BN_cmp_word(g, 1) > 0) {
 880       *out_result = bn_composite;
 881       ret = 1;
 882       goto err;
 883     }
 884
 885     // Step 4.5
 886     if (!BN_mod_exp_mont(z, b, m, w, ctx, mont)) {
 887       goto err;
 888     }
 889
 890     // Step 4.6
 891     if (BN_is_one(z) || BN_cmp(z, w1) == 0) {
 892       goto loop;
 893     }
 894
 895     // Step 4.7
 896     for (int j = 1; j < a; j++) {
 897       if (!BN_copy(x, z) || !BN_mod_mul(z, x, x, w, ctx)) {
 898         goto err;
 899       }
 900       if (BN_cmp(z, w1) == 0) {
 901         goto loop;
 902       }
 903       if (BN_is_one(z)) {
 904         goto composite;
 905       }
 906     }
 907
 908     // Step 4.8-4.9
 909     if (!BN_copy(x, z) || !BN_mod_mul(z, x, x, w, ctx)) {
 910       goto err;
 911     }
 912
 913     // Step 4.10-4.11
 914     if (!BN_is_one(z) && !BN_copy(x, z)) {
 915       goto err;
 916     }
 917
 918  composite:
 919     // Step 4.12-4.14
 920     if (!BN_copy(x1, x) ||
 921         !BN_sub_word(x1, 1) ||
 922         !BN_gcd(g, x1, w, ctx)) {
 923       goto err;
 924     }
 925     if (BN_cmp_word(g, 1) > 0) {
 926       *out_result = bn_composite;
 927     } else {
 928       *out_result = bn_non_prime_power_composite;
 929     }
 930
 931     ret = 1;
 932     goto err;
 933
 934  loop:
 935     // Step 4.15
 936     if (!BN_GENCB_call(cb, 1, i)) {
 937       goto err;
 938     }
 939   }
 940
 941   *out_result = bn_probably_prime;
 942   ret = 1;
 943
 944 err:
 945   BN_MONT_CTX_free(mont);
 946   BN_CTX_end(ctx);
 947
 948   return ret;
 949 }
 950
 951 static int probable_prime(BIGNUM *rnd, int bits) {
 952   int i;
 953   uint16_t mods[NUMPRIMES];
 954   BN_ULONG delta;
 955   BN_ULONG maxdelta = BN_MASK2 - primes[NUMPRIMES - 1];
 956   char is_single_word = bits <= BN_BITS2;
 957
 958 again:
 959   if (!BN_rand(rnd, bits, BN_RAND_TOP_TWO, BN_RAND_BOTTOM_ODD)) {
 960     return 0;
 961   }
 962
 963   // we now have a random number 'rnd' to test.
 964   for (i = 1; i < NUMPRIMES; i++) {
 965     mods[i] = bn_mod_u16_consttime(rnd, primes[i]);
 966   }
 967   // If bits is so small that it fits into a single word then we
 968   // additionally don't want to exceed that many bits.
 969   if (is_single_word) {
 970     BN_ULONG size_limit;
 971     if (bits == BN_BITS2) {
 972       // Avoid undefined behavior.
 973       size_limit = ~((BN_ULONG)0) - BN_get_word(rnd);
 974     } else {
 975       size_limit = (((BN_ULONG)1) << bits) - BN_get_word(rnd) - 1;
 976     }
 977     if (size_limit < maxdelta) {
 978       maxdelta = size_limit;
 979     }
 980   }
 981   delta = 0;
 982
 983 loop:
 984   if (is_single_word) {
 985     BN_ULONG rnd_word = BN_get_word(rnd);
 986
 987     // In the case that the candidate prime is a single word then
 988     // we check that:
 989     //   1) It's greater than primes[i] because we shouldn't reject
 990     //      3 as being a prime number because it's a multiple of
 991     //      three.
 992     //   2) That it's not a multiple of a known prime. We don't
 993     //      check that rnd-1 is also coprime to all the known
 994     //      primes because there aren't many small primes where
 995     //      that's true.
 996     for (i = 1; i < NUMPRIMES && primes[i] < rnd_word; i++) {
 997       if ((mods[i] + delta) % primes[i] == 0) {
 998         delta += 2;
 999         if (delta > maxdelta) {
1000           goto again;
1001         }
1002         goto loop;
1003       }
1004     }
1005   } else {
1006     for (i = 1; i < NUMPRIMES; i++) {
1007       // check that rnd is not a prime and also
1008       // that gcd(rnd-1,primes) == 1 (except for 2)
1009       if (((mods[i] + delta) % primes[i]) <= 1) {
1010         delta += 2;
1011         if (delta > maxdelta) {
1012           goto again;
1013         }
1014         goto loop;
1015       }
1016     }
1017   }
1018
1019   if (!BN_add_word(rnd, delta)) {
1020     return 0;
1021   }
1022   if (BN_num_bits(rnd) != (unsigned)bits) {
1023     goto again;
1024   }
1025
1026   return 1;
1027 }
1028
1029 static int probable_prime_dh(BIGNUM *rnd, int bits, const BIGNUM *add,
1030                              const BIGNUM *rem, BN_CTX *ctx) {
1031   int i, ret = 0;
1032   BIGNUM *t1;
1033
1034   BN_CTX_start(ctx);
1035   if ((t1 = BN_CTX_get(ctx)) == NULL) {
1036     goto err;
1037   }
1038
1039   if (!BN_rand(rnd, bits, BN_RAND_TOP_ONE, BN_RAND_BOTTOM_ODD)) {
1040     goto err;
1041   }
1042
1043   // we need ((rnd-rem) % add) == 0
1044
1045   if (!BN_mod(t1, rnd, add, ctx)) {
1046     goto err;
1047   }
1048   if (!BN_sub(rnd, rnd, t1)) {
1049     goto err;
1050   }
1051   if (rem == NULL) {
1052     if (!BN_add_word(rnd, 1)) {
1053       goto err;
1054     }
1055   } else {
1056     if (!BN_add(rnd, rnd, rem)) {
1057       goto err;
1058     }
1059   }
1060   // we now have a random number 'rand' to test.
1061
1062 loop:
1063   for (i = 1; i < NUMPRIMES; i++) {
1064     // check that rnd is a prime
1065     if (bn_mod_u16_consttime(rnd, primes[i]) <= 1) {
1066       if (!BN_add(rnd, rnd, add)) {
1067         goto err;
1068       }
1069       goto loop;
1070     }
1071   }
1072
1073   ret = 1;
1074
1075 err:
1076   BN_CTX_end(ctx);
1077   return ret;
1078 }
1079
1080 static int probable_prime_dh_safe(BIGNUM *p, int bits, const BIGNUM *padd,
1081                                   const BIGNUM *rem, BN_CTX *ctx) {
1082   int i, ret = 0;
1083   BIGNUM *t1, *qadd, *q;
1084
1085   bits--;
1086   BN_CTX_start(ctx);
1087   t1 = BN_CTX_get(ctx);
1088   q = BN_CTX_get(ctx);
1089   qadd = BN_CTX_get(ctx);
1090   if (qadd == NULL) {
1091     goto err;
1092   }
1093
1094   if (!BN_rshift1(qadd, padd)) {
1095     goto err;
1096   }
1097
1098   if (!BN_rand(q, bits, BN_RAND_TOP_ONE, BN_RAND_BOTTOM_ODD)) {
1099     goto err;
1100   }
1101
1102   // we need ((rnd-rem) % add) == 0
1103   if (!BN_mod(t1, q, qadd, ctx)) {
1104     goto err;
1105   }
1106
1107   if (!BN_sub(q, q, t1)) {
1108     goto err;
1109   }
1110
1111   if (rem == NULL) {
1112     if (!BN_add_word(q, 1)) {
1113       goto err;
1114     }
1115   } else {
1116     if (!BN_rshift1(t1, rem)) {
1117       goto err;
1118     }
1119     if (!BN_add(q, q, t1)) {
1120       goto err;
1121     }
1122   }
1123
1124   // we now have a random number 'rand' to test.
1125   if (!BN_lshift1(p, q)) {
1126     goto err;
1127   }
1128   if (!BN_add_word(p, 1)) {
1129     goto err;
1130   }
1131
1132 loop:
1133   for (i = 1; i < NUMPRIMES; i++) {
1134     // check that p and q are prime
1135     // check that for p and q
1136     // gcd(p-1,primes) == 1 (except for 2)
1137     if (bn_mod_u16_consttime(p, primes[i]) == 0 ||
1138         bn_mod_u16_consttime(q, primes[i]) == 0) {
1139       if (!BN_add(p, p, padd)) {
1140         goto err;
1141       }
1142       if (!BN_add(q, q, qadd)) {
1143         goto err;
1144       }
1145       goto loop;
1146     }
1147   }
1148
1149   ret = 1;
1150
1151 err:
1152   BN_CTX_end(ctx);
1153   return ret;
1154 }