2020 * implementation for N=30, using 30-bit signed limbs represented as int32_t.
2121 */
2222
23+ #ifdef VERIFY
24+ static const secp256k1_modinv32_signed30 SECP256K1_SIGNED30_ONE = {{1 }};
25+
26+ /* Compute a*factor and put it in r. All but the top limb in r will be in range [0,2^30). */
27+ static void secp256k1_modinv32_mul_30 (secp256k1_modinv32_signed30 * r , const secp256k1_modinv32_signed30 * a , int32_t factor ) {
28+ const int32_t M30 = (int32_t )(UINT32_MAX >> 2 );
29+ int64_t c = 0 ;
30+ int i ;
31+ for (i = 0 ; i < 8 ; ++ i ) {
32+ c += (int64_t )a -> v [i ] * factor ;
33+ r -> v [i ] = (int32_t )c & M30 ; c >>= 30 ;
34+ }
35+ c += (int64_t )a -> v [8 ] * factor ;
36+ VERIFY_CHECK (c == (int32_t )c );
37+ r -> v [8 ] = (int32_t )c ;
38+ }
39+
40+ /* Return -1 for a<b*factor, 0 for a==b*factor, 1 for a>b*factor. */
41+ static int secp256k1_modinv32_mul_cmp_30 (const secp256k1_modinv32_signed30 * a , const secp256k1_modinv32_signed30 * b , int32_t factor ) {
42+ int i ;
43+ secp256k1_modinv32_signed30 am , bm ;
44+ secp256k1_modinv32_mul_30 (& am , a , 1 ); /* Normalize all but the top limb of a. */
45+ secp256k1_modinv32_mul_30 (& bm , b , factor );
46+ for (i = 0 ; i < 8 ; ++ i ) {
47+ /* Verify that all but the top limb of a and b are normalized. */
48+ VERIFY_CHECK (am .v [i ] >> 30 == 0 );
49+ VERIFY_CHECK (bm .v [i ] >> 30 == 0 );
50+ }
51+ for (i = 8 ; i >= 0 ; -- i ) {
52+ if (am .v [i ] < bm .v [i ]) return -1 ;
53+ if (am .v [i ] > bm .v [i ]) return 1 ;
54+ }
55+ return 0 ;
56+ }
57+ #endif
58+
2359/* Take as input a signed30 number in range (-2*modulus,modulus), and add a multiple of the modulus
2460 * to it to bring it to range [0,modulus). If sign < 0, the input will also be negated in the
2561 * process. The input must have limbs in range (-2^30,2^30). The output will have limbs in range
@@ -30,6 +66,17 @@ static void secp256k1_modinv32_normalize_30(secp256k1_modinv32_signed30 *r, int3
3066 r5 = r -> v [5 ], r6 = r -> v [6 ], r7 = r -> v [7 ], r8 = r -> v [8 ];
3167 int32_t cond_add , cond_negate ;
3268
69+ #ifdef VERIFY
70+ /* Verify that all limbs are in range (-2^30,2^30). */
71+ int i ;
72+ for (i = 0 ; i < 9 ; ++ i ) {
73+ VERIFY_CHECK (r -> v [i ] >= - M30 );
74+ VERIFY_CHECK (r -> v [i ] <= M30 );
75+ }
76+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (r , & modinfo -> modulus , -2 ) > 0 ); /* r > -2*modulus */
77+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (r , & modinfo -> modulus , 1 ) < 0 ); /* r < modulus */
78+ #endif
79+
3380 /* In a first step, add the modulus if the input is negative, and then negate if requested.
3481 * This brings r from range (-2*modulus,modulus) to range (-modulus,modulus). As all input
3582 * limbs are in range (-2^30,2^30), this cannot overflow an int32_t. Note that the right
@@ -96,6 +143,20 @@ static void secp256k1_modinv32_normalize_30(secp256k1_modinv32_signed30 *r, int3
96143 r -> v [6 ] = r6 ;
97144 r -> v [7 ] = r7 ;
98145 r -> v [8 ] = r8 ;
146+
147+ #ifdef VERIFY
148+ VERIFY_CHECK (r0 >> 30 == 0 );
149+ VERIFY_CHECK (r1 >> 30 == 0 );
150+ VERIFY_CHECK (r2 >> 30 == 0 );
151+ VERIFY_CHECK (r3 >> 30 == 0 );
152+ VERIFY_CHECK (r4 >> 30 == 0 );
153+ VERIFY_CHECK (r5 >> 30 == 0 );
154+ VERIFY_CHECK (r6 >> 30 == 0 );
155+ VERIFY_CHECK (r7 >> 30 == 0 );
156+ VERIFY_CHECK (r8 >> 30 == 0 );
157+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (r , & modinfo -> modulus , 0 ) >= 0 ); /* r >= 0 */
158+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (r , & modinfo -> modulus , 1 ) < 0 ); /* r < modulus */
159+ #endif
99160}
100161
101162/* Data type for transition matrices (see section 3 of explanation).
@@ -155,12 +216,19 @@ static int32_t secp256k1_modinv32_divsteps_30(int32_t eta, uint32_t f0, uint32_t
155216 g >>= 1 ;
156217 u <<= 1 ;
157218 v <<= 1 ;
219+ /* Bounds on eta that follow from the bounds on iteration count (max 25*30 divsteps). */
220+ VERIFY_CHECK (eta >= -751 && eta <= 751 );
158221 }
159222 /* Return data in t and return value. */
160223 t -> u = (int32_t )u ;
161224 t -> v = (int32_t )v ;
162225 t -> q = (int32_t )q ;
163226 t -> r = (int32_t )r ;
227+ /* The determinant of t must be a power of two. This guarantees that multiplication with t
228+ * does not change the gcd of f and g, apart from adding a power-of-2 factor to it (which
229+ * will be divided out again). As each divstep's individual matrix has determinant 2, the
230+ * aggregate of 30 of them will have determinant 2^30. */
231+ VERIFY_CHECK ((int64_t )t -> u * t -> r - (int64_t )t -> v * t -> q == ((int64_t )1 ) << 30 );
164232 return eta ;
165233}
166234
@@ -211,6 +279,8 @@ static int32_t secp256k1_modinv32_divsteps_30_var(int32_t eta, uint32_t f0, uint
211279 VERIFY_CHECK ((g & 1 ) == 1 );
212280 VERIFY_CHECK ((u * f0 + v * g0 ) == f << (30 - i ));
213281 VERIFY_CHECK ((q * f0 + r * g0 ) == g << (30 - i ));
282+ /* Bounds on eta that follow from the bounds on iteration count (max 25*30 divsteps). */
283+ VERIFY_CHECK (eta >= -751 && eta <= 751 );
214284 /* If eta is negative, negate it and replace f,g with g,-f. */
215285 if (eta < 0 ) {
216286 uint32_t tmp ;
@@ -224,6 +294,7 @@ static int32_t secp256k1_modinv32_divsteps_30_var(int32_t eta, uint32_t f0, uint
224294 * can be done as its sign will flip once that happens. */
225295 limit = ((int )eta + 1 ) > i ? i : ((int )eta + 1 );
226296 /* m is a mask for the bottom min(limit, 8) bits (our table only supports 8 bits). */
297+ VERIFY_CHECK (limit > 0 && limit <= 30 );
227298 m = (UINT32_MAX >> (32 - limit )) & 255U ;
228299 /* Find what multiple of f must be added to g to cancel its bottom min(limit, 8) bits. */
229300 w = (g * inv256 [(f >> 1 ) & 127 ]) & m ;
@@ -238,6 +309,11 @@ static int32_t secp256k1_modinv32_divsteps_30_var(int32_t eta, uint32_t f0, uint
238309 t -> v = (int32_t )v ;
239310 t -> q = (int32_t )q ;
240311 t -> r = (int32_t )r ;
312+ /* The determinant of t must be a power of two. This guarantees that multiplication with t
313+ * does not change the gcd of f and g, apart from adding a power-of-2 factor to it (which
314+ * will be divided out again). As each divstep's individual matrix has determinant 2, the
315+ * aggregate of 30 of them will have determinant 2^30. */
316+ VERIFY_CHECK ((int64_t )t -> u * t -> r - (int64_t )t -> v * t -> q == ((int64_t )1 ) << 30 );
241317 return eta ;
242318}
243319
@@ -254,6 +330,16 @@ static void secp256k1_modinv32_update_de_30(secp256k1_modinv32_signed30 *d, secp
254330 int32_t di , ei , md , me , sd , se ;
255331 int64_t cd , ce ;
256332 int i ;
333+ #ifdef VERIFY
334+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (d , & modinfo -> modulus , -2 ) > 0 ); /* d > -2*modulus */
335+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (d , & modinfo -> modulus , 1 ) < 0 ); /* d < modulus */
336+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (e , & modinfo -> modulus , -2 ) > 0 ); /* e > -2*modulus */
337+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (e , & modinfo -> modulus , 1 ) < 0 ); /* e < modulus */
338+ VERIFY_CHECK ((labs (u ) + labs (v )) >= 0 ); /* |u|+|v| doesn't overflow */
339+ VERIFY_CHECK ((labs (q ) + labs (r )) >= 0 ); /* |q|+|r| doesn't overflow */
340+ VERIFY_CHECK ((labs (u ) + labs (v )) <= M30 + 1 ); /* |u|+|v| <= 2^30 */
341+ VERIFY_CHECK ((labs (q ) + labs (r )) <= M30 + 1 ); /* |q|+|r| <= 2^30 */
342+ #endif
257343 /* [md,me] start as zero; plus [u,q] if d is negative; plus [v,r] if e is negative. */
258344 sd = d -> v [8 ] >> 31 ;
259345 se = e -> v [8 ] >> 31 ;
@@ -288,6 +374,12 @@ static void secp256k1_modinv32_update_de_30(secp256k1_modinv32_signed30 *d, secp
288374 /* What remains is limb 9 of t*[d,e]+modulus*[md,me]; store it as output limb 8. */
289375 d -> v [8 ] = (int32_t )cd ;
290376 e -> v [8 ] = (int32_t )ce ;
377+ #ifdef VERIFY
378+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (d , & modinfo -> modulus , -2 ) > 0 ); /* d > -2*modulus */
379+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (d , & modinfo -> modulus , 1 ) < 0 ); /* d < modulus */
380+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (e , & modinfo -> modulus , -2 ) > 0 ); /* e > -2*modulus */
381+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (e , & modinfo -> modulus , 1 ) < 0 ); /* e < modulus */
382+ #endif
291383}
292384
293385/* Compute (t/2^30) * [f, g], where t is a transition matrix for 30 divsteps.
@@ -341,13 +433,35 @@ static void secp256k1_modinv32(secp256k1_modinv32_signed30 *x, const secp256k1_m
341433 /* Update d,e using that transition matrix. */
342434 secp256k1_modinv32_update_de_30 (& d , & e , & t , modinfo );
343435 /* Update f,g using that transition matrix. */
436+ #ifdef VERIFY
437+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , -1 ) > 0 ); /* f > -modulus */
438+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , 1 ) <= 0 ); /* f <= modulus */
439+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , -1 ) > 0 ); /* g > -modulus */
440+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , 1 ) < 0 ); /* g < modulus */
441+ #endif
344442 secp256k1_modinv32_update_fg_30 (& f , & g , & t );
443+ #ifdef VERIFY
444+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , -1 ) > 0 ); /* f > -modulus */
445+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , 1 ) <= 0 ); /* f <= modulus */
446+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , -1 ) > 0 ); /* g > -modulus */
447+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , 1 ) < 0 ); /* g < modulus */
448+ #endif
345449 }
346450
347451 /* At this point sufficient iterations have been performed that g must have reached 0
348452 * and (if g was not originally 0) f must now equal +/- GCD of the initial f, g
349453 * values i.e. +/- 1, and d now contains +/- the modular inverse. */
350- VERIFY_CHECK ((g .v [0 ] | g .v [1 ] | g .v [2 ] | g .v [3 ] | g .v [4 ] | g .v [5 ] | g .v [6 ] | g .v [7 ] | g .v [8 ]) == 0 );
454+ #ifdef VERIFY
455+ /* g == 0 */
456+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & SECP256K1_SIGNED30_ONE , 0 ) == 0 );
457+ /* |f| == 1, or (x == 0 and d == 0 and |f|=modulus) */
458+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & SECP256K1_SIGNED30_ONE , -1 ) == 0 ||
459+ secp256k1_modinv32_mul_cmp_30 (& f , & SECP256K1_SIGNED30_ONE , 1 ) == 0 ||
460+ (secp256k1_modinv32_mul_cmp_30 (x , & SECP256K1_SIGNED30_ONE , 0 ) == 0 &&
461+ secp256k1_modinv32_mul_cmp_30 (& d , & SECP256K1_SIGNED30_ONE , 0 ) == 0 &&
462+ (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , 1 ) == 0 ||
463+ secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , -1 ) == 0 )));
464+ #endif
351465
352466 /* Optionally negate d, normalize to [0,modulus), and return it. */
353467 secp256k1_modinv32_normalize_30 (& d , f .v [8 ], modinfo );
@@ -361,6 +475,9 @@ static void secp256k1_modinv32_var(secp256k1_modinv32_signed30 *x, const secp256
361475 secp256k1_modinv32_signed30 e = {{1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 }};
362476 secp256k1_modinv32_signed30 f = modinfo -> modulus ;
363477 secp256k1_modinv32_signed30 g = * x ;
478+ #ifdef VERIFY
479+ int i = 0 ;
480+ #endif
364481 int j ;
365482 int32_t eta = -1 ;
366483 int32_t cond ;
@@ -373,6 +490,12 @@ static void secp256k1_modinv32_var(secp256k1_modinv32_signed30 *x, const secp256
373490 /* Update d,e using that transition matrix. */
374491 secp256k1_modinv32_update_de_30 (& d , & e , & t , modinfo );
375492 /* Update f,g using that transition matrix. */
493+ #ifdef VERIFY
494+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , -1 ) > 0 ); /* f > -modulus */
495+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , 1 ) <= 0 ); /* f <= modulus */
496+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , -1 ) > 0 ); /* g > -modulus */
497+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , 1 ) < 0 ); /* g < modulus */
498+ #endif
376499 secp256k1_modinv32_update_fg_30 (& f , & g , & t );
377500 /* If the bottom limb of g is 0, there is a chance g=0. */
378501 if (g .v [0 ] == 0 ) {
@@ -384,10 +507,28 @@ static void secp256k1_modinv32_var(secp256k1_modinv32_signed30 *x, const secp256
384507 /* If so, we're done. */
385508 if (cond == 0 ) break ;
386509 }
510+ #ifdef VERIFY
511+ VERIFY_CHECK (++ i < 25 ); /* We should never need more than 25*30 = 750 divsteps */
512+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , -1 ) > 0 ); /* f > -modulus */
513+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , 1 ) <= 0 ); /* f <= modulus */
514+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , -1 ) > 0 ); /* g > -modulus */
515+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & modinfo -> modulus , 1 ) < 0 ); /* g < modulus */
516+ #endif
387517 }
388518
389519 /* At this point g is 0 and (if g was not originally 0) f must now equal +/- GCD of
390520 * the initial f, g values i.e. +/- 1, and d now contains +/- the modular inverse. */
521+ #ifdef VERIFY
522+ /* g == 0 */
523+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& g , & SECP256K1_SIGNED30_ONE , 0 ) == 0 );
524+ /* |f| == 1, or (x == 0 and d == 0 and |f|=modulus) */
525+ VERIFY_CHECK (secp256k1_modinv32_mul_cmp_30 (& f , & SECP256K1_SIGNED30_ONE , -1 ) == 0 ||
526+ secp256k1_modinv32_mul_cmp_30 (& f , & SECP256K1_SIGNED30_ONE , 1 ) == 0 ||
527+ (secp256k1_modinv32_mul_cmp_30 (x , & SECP256K1_SIGNED30_ONE , 0 ) == 0 &&
528+ secp256k1_modinv32_mul_cmp_30 (& d , & SECP256K1_SIGNED30_ONE , 0 ) == 0 &&
529+ (secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , 1 ) == 0 ||
530+ secp256k1_modinv32_mul_cmp_30 (& f , & modinfo -> modulus , -1 ) == 0 )));
531+ #endif
391532
392533 /* Optionally negate d, normalize to [0,modulus), and return it. */
393534 secp256k1_modinv32_normalize_30 (& d , f .v [8 ], modinfo );
0 commit comments