21768358452

Committed 06 Feb 2026 10:35PM UTC coverage: 90.064% (-0.003%) from 90.067%

Build # 21768358452

Build Type

Pull #5289

github

Committed by

web-flow

Commit Message

Merge f589db195 into 8ea0ca252

Pull Request Pull Request #5289: Further misc header reductions, forward declarations, etc

Run Details

102238 of 113517 relevant lines covered (90.06%)

11357432.36 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

96.69

/src/lib/math/numbertheory/numthry.cpp

/*
* Number Theory Functions
* (C) 1999-2011,2016,2018,2019 Jack Lloyd
* (C) 2007,2008 Falko Strenzke, FlexSecure GmbH
*
* Botan is released under the Simplified BSD License (see license.txt)
*/

#include <botan/numthry.h>

#include <botan/exceptn.h>
#include <botan/rng.h>
#include <botan/internal/barrett.h>
#include <botan/internal/ct_utils.h>
#include <botan/internal/divide.h>
#include <botan/internal/monty.h>
#include <botan/internal/monty_exp.h>
#include <botan/internal/mp_core.h>
#include <botan/internal/primality.h>
#include <algorithm>

namespace Botan {

/*
* Tonelli-Shanks algorithm
*/
BigInt sqrt_modulo_prime(const BigInt& a, const BigInt& p) {
   BOTAN_ARG_CHECK(p > 1, "invalid prime");
   BOTAN_ARG_CHECK(a < p, "value to solve for must be less than p");
   BOTAN_ARG_CHECK(a >= 0, "value to solve for must not be negative");

   // some very easy cases
   if(p == 2 || a <= 1) {
      return a;
   }

   BOTAN_ARG_CHECK(p.is_odd(), "invalid prime");

   if(jacobi(a, p) != 1) {  // not a quadratic residue
      return BigInt::from_s32(-1);
   }

   auto mod_p = Barrett_Reduction::for_public_modulus(p);
   const Montgomery_Params monty_p(p, mod_p);

   // If p == 3 (mod 4) there is a simple solution
   if(p % 4 == 3) {
      return monty_exp_vartime(monty_p, a, ((p + 1) >> 2)).value();
   }

   // Otherwise we have to use Shanks-Tonelli
   size_t s = low_zero_bits(p - 1);
   BigInt q = p >> s;

   q -= 1;
   q >>= 1;

   BigInt r = monty_exp_vartime(monty_p, a, q).value();
   BigInt n = mod_p.multiply(a, mod_p.square(r));
   r = mod_p.multiply(r, a);

   if(n == 1) {
      return r;
   }

   // find random quadratic nonresidue z
   word z = 2;
   for(;;) {
      if(jacobi(BigInt::from_word(z), p) == -1) {  // found one
         break;
      }

      z += 1;  // try next z

      /*
      * The expected number of tests to find a non-residue modulo a
      * prime is 2. If we have not found one after 256 then almost
      * certainly we have been given a non-prime p.
      */
      if(z >= 256) {
         return BigInt::from_s32(-1);
      }
   }

   BigInt c = monty_exp_vartime(monty_p, BigInt::from_word(z), (q << 1) + 1).value();

   while(n > 1) {
      q = n;

      size_t i = 0;
      while(q != 1) {
         q = mod_p.square(q);
         ++i;

         if(i >= s) {
            return BigInt::from_s32(-1);
         }
      }

      BOTAN_ASSERT_NOMSG(s >= (i + 1));  // No underflow!
      c = monty_exp_vartime(monty_p, c, BigInt::power_of_2(s - i - 1)).value();
      r = mod_p.multiply(r, c);
      c = mod_p.square(c);
      n = mod_p.multiply(n, c);

      // s decreases as the algorithm proceeds
      BOTAN_ASSERT_NOMSG(s >= i);
      s = i;
   }

   return r;
}

/*
* Calculate the Jacobi symbol
*
* See Algorithm 2.149 in Handbook of Applied Cryptography
*/
int32_t jacobi(BigInt a, BigInt n) {
   BOTAN_ARG_CHECK(n.is_odd() && n >= 3, "Argument n must be an odd integer >= 3");

   if(a < 0 || a >= n) {
      a %= n;
   }

   if(a == 0) {
      return 0;
   }
   if(a == 1) {
      return 1;
   }

   int32_t s = 1;

   for(;;) {
      const size_t e = low_zero_bits(a);
      a >>= e;
      const word n_mod_8 = n.word_at(0) % 8;
      const word n_mod_4 = n_mod_8 % 4;

      if(e % 2 == 1 && (n_mod_8 == 3 || n_mod_8 == 5)) {
         s = -s;
      }

      if(n_mod_4 == 3 && a % 4 == 3) {
         s = -s;
      }

      /*
      * The HAC presentation of the algorithm uses recursion, which is not
      * desirable or necessary.
      *
      * Instead we loop accumulating the product of the various jacobi()
      * subcomputations into s, until we reach algorithm termination, which
      * occurs in one of two ways.
      *
      * If a == 1 then the recursion has completed; we can return the value of s.
      *
      * Otherwise, after swapping and reducing, check for a == 0 [this value is
      * called `n1` in HAC's presentation]. This would imply that jacobi(n1,a1)
      * would have the value 0, due to Line 1 in HAC 2.149, in which case the
      * entire product is zero, and we can immediately return that result.
      */

      if(a == 1) {
         return s;
      }

      std::swap(a, n);

      BOTAN_ASSERT_NOMSG(n.is_odd());

      a %= n;

      if(a == 0) {
         return 0;
      }
   }
}

/*
* Square a BigInt
*/
BigInt square(const BigInt& x) {
   BigInt z = x;
   secure_vector<word> ws;
   z.square(ws);
   return z;
}

/*
* Return the number of 0 bits at the end of n
*/
size_t low_zero_bits(const BigInt& n) {
   size_t low_zero = 0;

   auto seen_nonempty_word = CT::Mask<word>::cleared();

   for(size_t i = 0; i != n.size(); ++i) {
      const word x = n.word_at(i);

      // ctz(0) will return sizeof(word)
      const size_t tz_x = ctz(x);

      // if x > 0 we want to count tz_x in total but not any
      // further words, so set the mask after the addition
      low_zero += seen_nonempty_word.if_not_set_return(tz_x);

      seen_nonempty_word |= CT::Mask<word>::expand(x);
   }

   // if we saw no words with x > 0 then n == 0 and the value we have
   // computed is meaningless. Instead return BigInt::zero() in that case.
   return static_cast<size_t>(seen_nonempty_word.if_set_return(low_zero));
}

/*
* Calculate the GCD in constant time
*/
BigInt gcd(const BigInt& a, const BigInt& b) {
   if(a.is_zero()) {
      return abs(b);
   }
   if(b.is_zero()) {
      return abs(a);
   }

   const size_t sz = std::max(a.sig_words(), b.sig_words());
   auto u = BigInt::with_capacity(sz);
   auto v = BigInt::with_capacity(sz);
   u += a;
   v += b;

   CT::poison_all(u, v);

   u.set_sign(BigInt::Positive);
   v.set_sign(BigInt::Positive);

   // In the worst case we have two fully populated big ints. After right
   // shifting so many times, we'll have reached the result for sure.
   const size_t loop_cnt = u.bits() + v.bits();

   // This temporary is big enough to hold all intermediate results of the
   // algorithm. No reallocation will happen during the loop.
   // Note however, that `ct_cond_assign()` will invalidate the 'sig_words'
   // cache, which _does not_ shrink the capacity of the underlying buffer.
   auto tmp = BigInt::with_capacity(sz);
   secure_vector<word> ws(sz * 2);
   size_t factors_of_two = 0;
   for(size_t i = 0; i != loop_cnt; ++i) {
      auto both_odd = CT::Mask<word>::expand_bool(u.is_odd()) & CT::Mask<word>::expand_bool(v.is_odd());

      // Subtract the smaller from the larger if both are odd
      auto u_gt_v = CT::Mask<word>::expand_bool(bigint_cmp(u._data(), u.size(), v._data(), v.size()) > 0);
      bigint_sub_abs(tmp.mutable_data(), u._data(), v._data(), sz, ws.data());
      u.ct_cond_assign((u_gt_v & both_odd).as_bool(), tmp);
      v.ct_cond_assign((~u_gt_v & both_odd).as_bool(), tmp);

      const auto u_is_even = CT::Mask<word>::expand_bool(u.is_even());
      const auto v_is_even = CT::Mask<word>::expand_bool(v.is_even());
      BOTAN_DEBUG_ASSERT((u_is_even | v_is_even).as_bool());

      // When both are even, we're going to eliminate a factor of 2.
      // We have to reapply this factor to the final result.
      factors_of_two += (u_is_even & v_is_even).if_set_return(1);

      // remove one factor of 2, if u is even
      bigint_shr2(tmp.mutable_data(), u._data(), sz, 1);
      u.ct_cond_assign(u_is_even.as_bool(), tmp);

      // remove one factor of 2, if v is even
      bigint_shr2(tmp.mutable_data(), v._data(), sz, 1);
      v.ct_cond_assign(v_is_even.as_bool(), tmp);
   }

   // The GCD (without factors of two) is either in u or v, the other one is
   // zero. The non-zero variable _must_ be odd, because all factors of two were
   // removed in the loop iterations above.
   BOTAN_DEBUG_ASSERT(u.is_zero() || v.is_zero());
   BOTAN_DEBUG_ASSERT(u.is_odd() || v.is_odd());

   // make sure that the GCD (without factors of two) is in u
   u.ct_cond_assign(u.is_even() /* .is_zero() would not be constant time */, v);

   // re-apply the factors of two
   u.ct_shift_left(factors_of_two);

   CT::unpoison_all(u, v);

   return u;
}

/*
* Calculate the LCM
*/
BigInt lcm(const BigInt& a, const BigInt& b) {
   if(a == b) {
      return a;
   }

   auto ab = a * b;
   ab.set_sign(BigInt::Positive);  // ignore the signs of a & b
   const auto g = gcd(a, b);
   return ct_divide(ab, g);
}

/*
* Modular Exponentiation
*/
BigInt power_mod(const BigInt& base, const BigInt& exp, const BigInt& mod) {
   if(mod.is_negative() || mod == 1) {
      return BigInt::zero();
   }

   if(base.is_zero() || mod.is_zero()) {
      if(exp.is_zero()) {
         return BigInt::one();
      }
      return BigInt::zero();
   }

   auto reduce_mod = Barrett_Reduction::for_secret_modulus(mod);

   const size_t exp_bits = exp.bits();

   if(mod.is_odd()) {
      const Montgomery_Params monty_params(mod, reduce_mod);
      return monty_exp(monty_params, ct_modulo(base, mod), exp, exp_bits).value();
   }

   /*
   Support for even modulus is just a convenience and not considered
   cryptographically important, so this implementation is slow ...
   */
   BigInt accum = BigInt::one();
   BigInt g = ct_modulo(base, mod);
   BigInt t;

   for(size_t i = 0; i != exp_bits; ++i) {
      t = reduce_mod.multiply(g, accum);
      g = reduce_mod.square(g);
      accum.ct_cond_assign(exp.get_bit(i), t);
   }
   return accum;
}

BigInt is_perfect_square(const BigInt& C) {
   if(C < 1) {
      throw Invalid_Argument("is_perfect_square requires C >= 1");
   }
   if(C == 1) {
      return BigInt::one();
   }

   const size_t n = C.bits();
   const size_t m = (n + 1) / 2;
   const BigInt B = C + BigInt::power_of_2(m);

   BigInt X = BigInt::power_of_2(m) - 1;
   BigInt X2 = (X * X);

   for(;;) {
      X = (X2 + C) / (2 * X);
      X2 = (X * X);

      if(X2 < B) {
         break;
      }
   }

   if(X2 == C) {
      return X;
   } else {
      return BigInt::zero();
   }
}

/*
* Test for primality using Miller-Rabin
*/
bool is_prime(const BigInt& n, RandomNumberGenerator& rng, size_t prob, bool is_random) {
   if(n == 2) {
      return true;
   }
   if(n <= 1 || n.is_even()) {
      return false;
   }

   const size_t n_bits = n.bits();

   // Fast path testing for small numbers (<= 65521)
   if(n_bits <= 16) {
      const uint16_t num = static_cast<uint16_t>(n.word_at(0));

      return std::binary_search(PRIMES, PRIMES + PRIME_TABLE_SIZE, num);
   }

   auto mod_n = Barrett_Reduction::for_secret_modulus(n);

   if(rng.is_seeded()) {
      const size_t t = miller_rabin_test_iterations(n_bits, prob, is_random);

      if(!is_miller_rabin_probable_prime(n, mod_n, rng, t)) {
         return false;
      }

      if(is_random) {
         return true;
      } else {
         return is_lucas_probable_prime(n, mod_n);
      }
   } else {
      return is_bailie_psw_probable_prime(n, mod_n);
   }
}

}  // namespace Botan

1	/*
2	* Number Theory Functions
3	* (C) 1999-2011,2016,2018,2019 Jack Lloyd
4	* (C) 2007,2008 Falko Strenzke, FlexSecure GmbH
5	*
6	* Botan is released under the Simplified BSD License (see license.txt)
7	*/
8
9	#include <botan/numthry.h>
10
11	#include <botan/exceptn.h>
12	#include <botan/rng.h>
13	#include <botan/internal/barrett.h>
14	#include <botan/internal/ct_utils.h>
15	#include <botan/internal/divide.h>
16	#include <botan/internal/monty.h>
17	#include <botan/internal/monty_exp.h>
18	#include <botan/internal/mp_core.h>
19	#include <botan/internal/primality.h>
20	#include <algorithm>
21
22	namespace Botan {
23
24	/*
25	* Tonelli-Shanks algorithm
26	*/
27	BigInt sqrt_modulo_prime(const BigInt& a, const BigInt& p) {	1,473✔
28	BOTAN_ARG_CHECK(p > 1, "invalid prime");	1,473✔
29	BOTAN_ARG_CHECK(a < p, "value to solve for must be less than p");	1,473✔
30	BOTAN_ARG_CHECK(a >= 0, "value to solve for must not be negative");	1,472✔
31
32	// some very easy cases
33	if(p == 2 \|\| a <= 1) {	2,942✔
34	return a;	3✔
35	}
36
37	BOTAN_ARG_CHECK(p.is_odd(), "invalid prime");	2,938✔
38
39	if(jacobi(a, p) != 1) { // not a quadratic residue	1,469✔
40	return BigInt::from_s32(-1);	127✔
41	}
42
43	auto mod_p = Barrett_Reduction::for_public_modulus(p);	1,342✔
44	const Montgomery_Params monty_p(p, mod_p);	1,342✔
45
46	// If p == 3 (mod 4) there is a simple solution
47	if(p % 4 == 3) {	1,342✔
48	return monty_exp_vartime(monty_p, a, ((p + 1) >> 2)).value();	984✔
49	}
50
51	// Otherwise we have to use Shanks-Tonelli
52	size_t s = low_zero_bits(p - 1);	358✔
53	BigInt q = p >> s;	358✔
54
55	q -= 1;	358✔
56	q >>= 1;	358✔
57
58	BigInt r = monty_exp_vartime(monty_p, a, q).value();	358✔
59	BigInt n = mod_p.multiply(a, mod_p.square(r));	358✔
60	r = mod_p.multiply(r, a);	358✔
61
62	if(n == 1) {	358✔
63	return r;	171✔
64	}
65
66	// find random quadratic nonresidue z
67	word z = 2;
68	for(;;) {	503✔
69	if(jacobi(BigInt::from_word(z), p) == -1) { // found one	503✔
70	break;
71	}
72
73	z += 1; // try next z	317✔
74
75	/*
76	* The expected number of tests to find a non-residue modulo a
77	* prime is 2. If we have not found one after 256 then almost
78	* certainly we have been given a non-prime p.
79	*/
80	if(z >= 256) {	317✔
81	return BigInt::from_s32(-1);	1✔
82	}
83	}
84
85	BigInt c = monty_exp_vartime(monty_p, BigInt::from_word(z), (q << 1) + 1).value();	372✔
86
87	while(n > 1) {	589✔
88	q = n;	409✔
89
90	size_t i = 0;	409✔
91	while(q != 1) {	11,771✔
92	q = mod_p.square(q);	11,368✔
93	++i;	11,368✔
94
95	if(i >= s) {	11,368✔
96	return BigInt::from_s32(-1);	6✔
97	}
98	}
99
100	BOTAN_ASSERT_NOMSG(s >= (i + 1)); // No underflow!	403✔
101	c = monty_exp_vartime(monty_p, c, BigInt::power_of_2(s - i - 1)).value();	403✔
102	r = mod_p.multiply(r, c);	403✔
103	c = mod_p.square(c);	403✔
104	n = mod_p.multiply(n, c);	403✔
105
106	// s decreases as the algorithm proceeds
107	BOTAN_ASSERT_NOMSG(s >= i);	403✔
108	s = i;
109	}
110
111	return r;	180✔
112	}	3,042✔
113
114	/*
115	* Calculate the Jacobi symbol
116	*
117	* See Algorithm 2.149 in Handbook of Applied Cryptography
118	*/
119	int32_t jacobi(BigInt a, BigInt n) {	100,593✔
120	BOTAN_ARG_CHECK(n.is_odd() && n >= 3, "Argument n must be an odd integer >= 3");	301,779✔
121
122	if(a < 0 \|\| a >= n) {	170,569✔
123	a %= n;	30,889✔
124	}
125
126	if(a == 0) {	100,593✔
127	return 0;
128	}
129	if(a == 1) {	100,570✔
130	return 1;
131	}
132
133	int32_t s = 1;
134
135	for(;;) {	385,734✔
136	const size_t e = low_zero_bits(a);	385,734✔
137	a >>= e;	385,734✔
138	const word n_mod_8 = n.word_at(0) % 8;	385,734✔
139	const word n_mod_4 = n_mod_8 % 4;	385,734✔
140
141	if(e % 2 == 1 && (n_mod_8 == 3 \|\| n_mod_8 == 5)) {	385,734✔
142	s = -s;	66,452✔
143	}
144
145	if(n_mod_4 == 3 && a % 4 == 3) {	385,734✔
146	s = -s;	67,457✔
147	}
148
149	/*
150	* The HAC presentation of the algorithm uses recursion, which is not
151	* desirable or necessary.
152	*
153	* Instead we loop accumulating the product of the various jacobi()
154	* subcomputations into s, until we reach algorithm termination, which
155	* occurs in one of two ways.
156	*
157	* If a == 1 then the recursion has completed; we can return the value of s.
158	*
159	* Otherwise, after swapping and reducing, check for a == 0 [this value is
160	* called `n1` in HAC's presentation]. This would imply that jacobi(n1,a1)
161	* would have the value 0, due to Line 1 in HAC 2.149, in which case the
162	* entire product is zero, and we can immediately return that result.
163	*/
164
165	if(a == 1) {	385,734✔
166	return s;	83,944✔
167	}
168
169	std::swap(a, n);	301,790✔
170
171	BOTAN_ASSERT_NOMSG(n.is_odd());	603,580✔
172
173	a %= n;	301,790✔
174
175	if(a == 0) {	301,790✔
176	return 0;
177	}
178	}
179	}
180
181	/*
182	* Square a BigInt
183	*/
184	BigInt square(const BigInt& x) {	1,021✔
185	BigInt z = x;	1,021✔
186	secure_vector<word> ws;	1,021✔
187	z.square(ws);	1,021✔
188	return z;	1,021✔
189	}	1,021✔
190
191	/*
192	* Return the number of 0 bits at the end of n
193	*/
194	size_t low_zero_bits(const BigInt& n) {	601,369✔
195	size_t low_zero = 0;	601,369✔
196
197	auto seen_nonempty_word = CT::Mask<word>::cleared();	601,369✔
198
199	for(size_t i = 0; i != n.size(); ++i) {	6,091,794✔
200	const word x = n.word_at(i);	5,490,425✔
201
202	// ctz(0) will return sizeof(word)
203	const size_t tz_x = ctz(x);	5,490,425✔
204
205	// if x > 0 we want to count tz_x in total but not any
206	// further words, so set the mask after the addition
207	low_zero += seen_nonempty_word.if_not_set_return(tz_x);	5,490,425✔
208
209	seen_nonempty_word \|= CT::Mask<word>::expand(x);	5,490,425✔
210	}
211
212	// if we saw no words with x > 0 then n == 0 and the value we have
213	// computed is meaningless. Instead return BigInt::zero() in that case.
214	return static_cast<size_t>(seen_nonempty_word.if_set_return(low_zero));	601,369✔
215	}
216
217	/*
218	* Calculate the GCD in constant time
219	*/
220	BigInt gcd(const BigInt& a, const BigInt& b) {	26,569✔
221	if(a.is_zero()) {	28,696✔
222	return abs(b);	3✔
223	}
224	if(b.is_zero()) {	27,006✔
225	return abs(a);	26,569✔
226	}
227
228	const size_t sz = std::max(a.sig_words(), b.sig_words());	26,564✔
229	auto u = BigInt::with_capacity(sz);	26,564✔
230	auto v = BigInt::with_capacity(sz);	26,564✔
231	u += a;	26,564✔
232	v += b;	26,564✔
233
234	CT::poison_all(u, v);	26,564✔
235
236	u.set_sign(BigInt::Positive);	26,564✔
237	v.set_sign(BigInt::Positive);	26,564✔
238
239	// In the worst case we have two fully populated big ints. After right
240	// shifting so many times, we'll have reached the result for sure.
241	const size_t loop_cnt = u.bits() + v.bits();	26,564✔
242
243	// This temporary is big enough to hold all intermediate results of the
244	// algorithm. No reallocation will happen during the loop.
245	// Note however, that `ct_cond_assign()` will invalidate the 'sig_words'
246	// cache, which _does not_ shrink the capacity of the underlying buffer.
247	auto tmp = BigInt::with_capacity(sz);	26,564✔
248	secure_vector<word> ws(sz * 2);	26,564✔
249	size_t factors_of_two = 0;	26,564✔
250	for(size_t i = 0; i != loop_cnt; ++i) {	5,653,958✔
251	auto both_odd = CT::Mask<word>::expand_bool(u.is_odd()) & CT::Mask<word>::expand_bool(v.is_odd());	16,882,182✔
252
253	// Subtract the smaller from the larger if both are odd
254	auto u_gt_v = CT::Mask<word>::expand_bool(bigint_cmp(u._data(), u.size(), v._data(), v.size()) > 0);	5,627,394✔
255	bigint_sub_abs(tmp.mutable_data(), u._data(), v._data(), sz, ws.data());	5,627,394✔
256	u.ct_cond_assign((u_gt_v & both_odd).as_bool(), tmp);	5,627,394✔
257	v.ct_cond_assign((~u_gt_v & both_odd).as_bool(), tmp);	5,627,394✔
258
259	const auto u_is_even = CT::Mask<word>::expand_bool(u.is_even());	11,254,788✔
260	const auto v_is_even = CT::Mask<word>::expand_bool(v.is_even());	11,254,788✔
261	BOTAN_DEBUG_ASSERT((u_is_even \| v_is_even).as_bool());	5,627,394✔
262
263	// When both are even, we're going to eliminate a factor of 2.
264	// We have to reapply this factor to the final result.
265	factors_of_two += (u_is_even & v_is_even).if_set_return(1);	5,627,394✔
266
267	// remove one factor of 2, if u is even
268	bigint_shr2(tmp.mutable_data(), u._data(), sz, 1);	5,627,394✔
269	u.ct_cond_assign(u_is_even.as_bool(), tmp);	5,627,394✔
270
271	// remove one factor of 2, if v is even
272	bigint_shr2(tmp.mutable_data(), v._data(), sz, 1);	5,627,394✔
273	v.ct_cond_assign(v_is_even.as_bool(), tmp);	5,627,394✔
274	}
275
276	// The GCD (without factors of two) is either in u or v, the other one is
277	// zero. The non-zero variable _must_ be odd, because all factors of two were
278	// removed in the loop iterations above.
279	BOTAN_DEBUG_ASSERT(u.is_zero() \|\| v.is_zero());	26,564✔
280	BOTAN_DEBUG_ASSERT(u.is_odd() \|\| v.is_odd());	26,564✔
281
282	// make sure that the GCD (without factors of two) is in u
283	u.ct_cond_assign(u.is_even() /* .is_zero() would not be constant time */, v);	53,128✔
284
285	// re-apply the factors of two
286	u.ct_shift_left(factors_of_two);	26,564✔
287
288	CT::unpoison_all(u, v);	26,564✔
289
290	return u;	26,564✔
291	}	26,564✔
292
293	/*
294	* Calculate the LCM
295	*/
296	BigInt lcm(const BigInt& a, const BigInt& b) {	4,609✔
297	if(a == b) {	4,609✔
298	return a;	×
299	}
300
301	auto ab = a * b;	4,609✔
302	ab.set_sign(BigInt::Positive); // ignore the signs of a & b	4,609✔
303	const auto g = gcd(a, b);	4,609✔
304	return ct_divide(ab, g);	4,609✔
305	}	4,609✔
306
307	/*
308	* Modular Exponentiation
309	*/
310	BigInt power_mod(const BigInt& base, const BigInt& exp, const BigInt& mod) {	65✔
311	if(mod.is_negative() \|\| mod == 1) {	130✔
312	return BigInt::zero();	×
313	}
314
315	if(base.is_zero() \|\| mod.is_zero()) {	188✔
316	if(exp.is_zero()) {	2✔
317	return BigInt::one();	1✔
318	}
319	return BigInt::zero();	×
320	}
321
322	auto reduce_mod = Barrett_Reduction::for_secret_modulus(mod);	64✔
323
324	const size_t exp_bits = exp.bits();	64✔
325
326	if(mod.is_odd()) {	64✔
327	const Montgomery_Params monty_params(mod, reduce_mod);	50✔
328	return monty_exp(monty_params, ct_modulo(base, mod), exp, exp_bits).value();	50✔
329	}	50✔
330
331	/*
332	Support for even modulus is just a convenience and not considered
333	cryptographically important, so this implementation is slow ...
334	*/
335	BigInt accum = BigInt::one();	14✔
336	BigInt g = ct_modulo(base, mod);	14✔
337	BigInt t;	14✔
338
339	for(size_t i = 0; i != exp_bits; ++i) {	3,071✔
340	t = reduce_mod.multiply(g, accum);	3,057✔
341	g = reduce_mod.square(g);	3,057✔
342	accum.ct_cond_assign(exp.get_bit(i), t);	6,114✔
343	}
344	return accum;	14✔
345	}	78✔
346
347	BigInt is_perfect_square(const BigInt& C) {	703✔
348	if(C < 1) {	703✔
349	throw Invalid_Argument("is_perfect_square requires C >= 1");	×
350	}
351	if(C == 1) {	703✔
352	return BigInt::one();	×
353	}
354
355	const size_t n = C.bits();	703✔
356	const size_t m = (n + 1) / 2;	703✔
357	const BigInt B = C + BigInt::power_of_2(m);	703✔
358
359	BigInt X = BigInt::power_of_2(m) - 1;	1,406✔
360	BigInt X2 = (X * X);	703✔
361
362	for(;;) {	1,657✔
363	X = (X2 + C) / (2 * X);	1,657✔
364	X2 = (X * X);	1,657✔
365
366	if(X2 < B) {	1,657✔
367	break;
368	}
369	}
370
371	if(X2 == C) {	703✔
372	return X;	63✔
373	} else {
374	return BigInt::zero();	640✔
375	}
376	}	703✔
377
378	/*
379	* Test for primality using Miller-Rabin
380	*/
381	bool is_prime(const BigInt& n, RandomNumberGenerator& rng, size_t prob, bool is_random) {	52,402✔
382	if(n == 2) {	52,402✔
383	return true;
384	}
385	if(n <= 1 \|\| n.is_even()) {	104,798✔
386	return false;
387	}
388
389	const size_t n_bits = n.bits();	52,331✔
390
391	// Fast path testing for small numbers (<= 65521)
392	if(n_bits <= 16) {	52,331✔
393	const uint16_t num = static_cast<uint16_t>(n.word_at(0));	32,802✔
394
395	return std::binary_search(PRIMES, PRIMES + PRIME_TABLE_SIZE, num);	32,802✔
396	}
397
398	auto mod_n = Barrett_Reduction::for_secret_modulus(n);	19,529✔
399
400	if(rng.is_seeded()) {	19,529✔
401	const size_t t = miller_rabin_test_iterations(n_bits, prob, is_random);	19,529✔
402
403	if(!is_miller_rabin_probable_prime(n, mod_n, rng, t)) {	19,529✔
404	return false;
405	}
406
407	if(is_random) {	3,328✔
408	return true;
409	} else {
410	return is_lucas_probable_prime(n, mod_n);	3,154✔
411	}
412	} else {
413	return is_bailie_psw_probable_prime(n, mod_n);	×
414	}
415	}	19,529✔
416
417	} // namespace Botan

randombit / botan / 21768358452

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous