Line data Source code
1 : #include "./fd_bn254_field_inl.h"
2 : #include "./fd_bn254_glv.h"
3 :
4 : /* G1 */
5 :
6 : /* GLV consts, decl in fd_bn254_glv.h, shared with fd_bn254_g2.c. */
7 :
8 : /* beta in Montgomery form.
9 : 0x30644e72e131a0295e6dd9e7e0acccb0c28f069fbb966e3de4bd44e5607cfd48 */
10 : const fd_bn254_fp_t fd_bn254_const_beta_mont[1] = {{{
11 : 0x3350c88e13e80b9cUL, 0x7dce557cdb5e56b9UL, 0x6001b4b8b615564aUL, 0x2682e617020217e0UL
12 : }}};
13 :
14 : /* Lattice constants, see glv.py */
15 : const ulong na[ 2 ] = { 0x8211bbeb7d4f1128UL, 0x6f4d8248eeb859fcUL };
16 : const ulong nb[ 1 ] = { 0x89d3256894d213e3UL };
17 : const ulong nc[ 2 ] = { 0x0be4e1541221250bUL, 0x6f4d8248eeb859fdUL };
18 :
19 : /* g2 = round(2^256 * N_B / r), 66-bit (2 limbs). Same for G1 and G2. */
20 : const ulong g2_const[ 2 ] = { 0xd91d232ec7e0b3d7UL, 0x0000000000000002UL };
21 :
22 : static inline fd_bn254_g1_t *
23 : fd_bn254_g1_set( fd_bn254_g1_t * r,
24 60213 : fd_bn254_g1_t const * p ) {
25 60213 : fd_bn254_fp_set( &r->X, &p->X );
26 60213 : fd_bn254_fp_set( &r->Y, &p->Y );
27 60213 : fd_bn254_fp_set( &r->Z, &p->Z );
28 60213 : return r;
29 60213 : }
30 :
31 : static inline fd_bn254_g1_t *
32 54 : fd_bn254_g1_set_zero( fd_bn254_g1_t * r ) {
33 : // fd_bn254_fp_set_zero( &r->X );
34 : // fd_bn254_fp_set_zero( &r->Y );
35 54 : fd_bn254_fp_set_zero( &r->Z );
36 54 : return r;
37 54 : }
38 :
39 : static inline fd_bn254_g1_t *
40 : fd_bn254_g1_to_affine( fd_bn254_g1_t * r,
41 60165 : fd_bn254_g1_t const * p ) {
42 60165 : if( FD_UNLIKELY( fd_bn254_fp_is_zero( &p->Z ) || fd_bn254_fp_is_one( &p->Z ) ) ) {
43 30060 : return fd_bn254_g1_set( r, p );
44 30060 : }
45 :
46 30105 : fd_bn254_fp_t iz[1], iz2[1];
47 30105 : fd_bn254_fp_inv( iz, &p->Z );
48 30105 : fd_bn254_fp_sqr( iz2, iz );
49 :
50 : /* X / Z^2, Y / Z^3 */
51 30105 : fd_bn254_fp_mul( &r->X, &p->X, iz2 );
52 30105 : fd_bn254_fp_mul( &r->Y, &p->Y, iz2 );
53 30105 : fd_bn254_fp_mul( &r->Y, &r->Y, iz );
54 30105 : fd_bn254_fp_set_one( &r->Z );
55 30105 : return r;
56 60165 : }
57 :
58 : uchar *
59 : fd_bn254_g1_tobytes( uchar out[64],
60 : fd_bn254_g1_t const * p,
61 60183 : int big_endian ) {
62 60183 : if( FD_UNLIKELY( fd_bn254_g1_is_zero( p ) ) ) {
63 18 : fd_memset( out, 0, 64UL );
64 : /* no flags */
65 18 : return out;
66 18 : }
67 :
68 60165 : fd_bn254_g1_t r[1];
69 60165 : fd_bn254_g1_to_affine( r, p );
70 :
71 60165 : fd_bn254_fp_from_mont( &r->X, &r->X );
72 60165 : fd_bn254_fp_from_mont( &r->Y, &r->Y );
73 :
74 60165 : fd_bn254_fp_tobytes_nm( &out[ 0], &r->X, big_endian );
75 60165 : fd_bn254_fp_tobytes_nm( &out[32], &r->Y, big_endian );
76 : /* no flags */
77 60165 : return out;
78 60183 : }
79 :
80 : /* fd_bn254_g1_affine_add computes r = p + q.
81 : Both p, q are affine, i.e. Z==1. */
82 : fd_bn254_g1_t *
83 : fd_bn254_g1_affine_add( fd_bn254_g1_t * r,
84 : fd_bn254_g1_t const * p,
85 60180 : fd_bn254_g1_t const * q ) {
86 : /* p==0, return q */
87 60180 : if( FD_UNLIKELY( fd_bn254_g1_is_zero( p ) ) ) {
88 21 : return fd_bn254_g1_set( r, q );
89 21 : }
90 : /* q==0, return p */
91 60159 : if( FD_UNLIKELY( fd_bn254_g1_is_zero( q ) ) ) {
92 9 : return fd_bn254_g1_set( r, p );
93 9 : }
94 :
95 60150 : fd_bn254_fp_t lambda[1], x[1], y[1];
96 :
97 : /* same X, either the points are equal or opposite */
98 60150 : if( fd_bn254_fp_eq( &p->X, &q->X ) ) {
99 6 : if( fd_bn254_fp_eq( &p->Y, &q->Y ) ) {
100 : /* p==q => point double: lambda = 3 * x1^2 / (2 * y1) */
101 6 : fd_bn254_fp_sqr( x, &p->X ); /* x = x1^2 */
102 6 : fd_bn254_fp_add( y, x, x ); /* y = 2 x1^2 */
103 6 : fd_bn254_fp_add( x, x, y ); /* x = 3 x1^2 */
104 6 : fd_bn254_fp_add( y, &p->Y, &p->Y );
105 6 : fd_bn254_fp_inv( lambda, y );
106 6 : fd_bn254_fp_mul( lambda, lambda, x );
107 6 : } else {
108 : /* p==-q => r=0 */
109 : /* COV: this may never happen with real data */
110 0 : return fd_bn254_g1_set_zero( r );
111 0 : }
112 60144 : } else {
113 : /* point add: lambda = (y1 - y2) / (x1 - x2) */
114 60144 : fd_bn254_fp_sub( x, &p->X, &q->X );
115 60144 : fd_bn254_fp_sub( y, &p->Y, &q->Y );
116 60144 : fd_bn254_fp_inv( lambda, x );
117 60144 : fd_bn254_fp_mul( lambda, lambda, y );
118 60144 : }
119 :
120 : /* x3 = lambda^2 - x1 - x2 */
121 60150 : fd_bn254_fp_sqr( x, lambda );
122 60150 : fd_bn254_fp_sub( x, x, &p->X );
123 60150 : fd_bn254_fp_sub( x, x, &q->X );
124 :
125 : /* y3 = lambda * (x1 - x3) - y1 */
126 60150 : fd_bn254_fp_sub( y, &p->X, x );
127 60150 : fd_bn254_fp_mul( y, y, lambda );
128 60150 : fd_bn254_fp_sub( y, y, &p->Y );
129 :
130 60150 : fd_bn254_fp_set( &r->X, x );
131 60150 : fd_bn254_fp_set( &r->Y, y );
132 60150 : fd_bn254_fp_set_one( &r->Z );
133 60150 : return r;
134 60150 : }
135 :
136 : /* fd_bn254_g1_dbl computes r = 2p.
137 : https://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#doubling-dbl-2009-l */
138 : fd_bn254_g1_t *
139 : fd_bn254_g1_dbl( fd_bn254_g1_t * r,
140 3790419 : fd_bn254_g1_t const * p ) {
141 : /* p==0, return 0 */
142 3790419 : if( FD_UNLIKELY( fd_bn254_g1_is_zero( p ) ) ) {
143 0 : return fd_bn254_g1_set_zero( r );
144 0 : }
145 :
146 3790419 : fd_bn254_fp_t a[1], b[1], c[1];
147 3790419 : fd_bn254_fp_t d[1], e[1], f[1];
148 :
149 : /* A = X1^2 */
150 3790419 : fd_bn254_fp_sqr( a, &p->X );
151 : /* B = Y1^2 */
152 3790419 : fd_bn254_fp_sqr( b, &p->Y );
153 : /* C = B^2 */
154 3790419 : fd_bn254_fp_sqr( c, b );
155 : /* D = 2*((X1+B)^2-A-C)
156 : (X1+B)^2 = X1^2 + 2*X1*B + B^2
157 : D = 2*(X1^2 + 2*X1*B + B^2 - A - C)
158 : D = 2*(X1^2 + 2*X1*B + B^2 - X1^2 - B^2)
159 : ^ ^ ^ ^
160 : |---------------|-----| |
161 : |------------|
162 : These terms cancel each other out, and we're left with:
163 : D = 2*(2*X1*B) */
164 3790419 : fd_bn254_fp_mul( d, &p->X, b );
165 3790419 : fd_bn254_fp_add( d, d, d );
166 3790419 : fd_bn254_fp_add( d, d, d );
167 : /* E = 3*A */
168 3790419 : fd_bn254_fp_add( e, a, a );
169 3790419 : fd_bn254_fp_add( e, a, e );
170 : /* F = E^2 */
171 3790419 : fd_bn254_fp_sqr( f, e );
172 : /* X3 = F-2*D */
173 3790419 : fd_bn254_fp_add( &r->X, d, d );
174 3790419 : fd_bn254_fp_sub( &r->X, f, &r->X );
175 : /* Z3 = (Y1+Z1)^2-YY-ZZ
176 : note: compute Z3 before Y3 because it depends on p->Y,
177 : that might be overwritten if r==p. */
178 : /* Z3 = 2*Y1*Z1 */
179 3790419 : fd_bn254_fp_mul( &r->Z, &p->Y, &p->Z );
180 3790419 : fd_bn254_fp_add( &r->Z, &r->Z, &r->Z );
181 : /* Y3 = E*(D-X3)-8*C */
182 3790419 : fd_bn254_fp_sub( &r->Y, d, &r->X );
183 3790419 : fd_bn254_fp_mul( &r->Y, e, &r->Y );
184 3790419 : fd_bn254_fp_add( c, c, c ); /* 2*c */
185 3790419 : fd_bn254_fp_add( c, c, c ); /* 4*y */
186 3790419 : fd_bn254_fp_add( c, c, c ); /* 8*y */
187 3790419 : fd_bn254_fp_sub( &r->Y, &r->Y, c );
188 3790419 : return r;
189 3790419 : }
190 :
191 : /* fd_bn254_g1_add_mixed computes r = p + q, when q->Z==1.
192 : http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-madd-2007-bl */
193 : fd_bn254_g1_t *
194 : fd_bn254_g1_add_mixed( fd_bn254_g1_t * r,
195 : fd_bn254_g1_t const * p,
196 2376864 : fd_bn254_g1_t const * q ) {
197 : /* p==0, return q */
198 2376864 : if( FD_UNLIKELY( fd_bn254_g1_is_zero( p ) ) ) {
199 0 : return fd_bn254_g1_set( r, q );
200 0 : }
201 2376864 : fd_bn254_fp_t zz[1], u2[1], s2[1];
202 2376864 : fd_bn254_fp_t h[1], hh[1];
203 2376864 : fd_bn254_fp_t i[1], j[1];
204 2376864 : fd_bn254_fp_t rr[1], v[1];
205 : /* Z1Z1 = Z1^2 */
206 2376864 : fd_bn254_fp_sqr( zz, &p->Z );
207 : /* U2 = X2*Z1Z1 */
208 2376864 : fd_bn254_fp_mul( u2, &q->X, zz );
209 : /* S2 = Y2*Z1*Z1Z1 */
210 2376864 : fd_bn254_fp_mul( s2, &q->Y, &p->Z );
211 2376864 : fd_bn254_fp_mul( s2, s2, zz );
212 :
213 : /* if p==q, call fd_bn254_g1_dbl */
214 2376864 : if( FD_UNLIKELY( fd_bn254_fp_eq( u2, &p->X ) && fd_bn254_fp_eq( s2, &p->Y ) ) ) {
215 : /* COV: this may never happen with real data */
216 0 : return fd_bn254_g1_dbl( r, p );
217 0 : }
218 :
219 : /* H = U2-X1 */
220 2376864 : fd_bn254_fp_sub( h, u2, &p->X );
221 : /* HH = H^2 */
222 2376864 : fd_bn254_fp_sqr( hh, h );
223 : /* I = 4*HH */
224 2376864 : fd_bn254_fp_add( i, hh, hh );
225 2376864 : fd_bn254_fp_add( i, i, i );
226 : /* J = H*I */
227 2376864 : fd_bn254_fp_mul( j, h, i );
228 : /* r = 2*(S2-Y1) */
229 2376864 : fd_bn254_fp_sub( rr, s2, &p->Y );
230 2376864 : fd_bn254_fp_add( rr, rr, rr );
231 : /* V = X1*I */
232 2376864 : fd_bn254_fp_mul( v, &p->X, i );
233 : /* X3 = r^2-J-2*V */
234 2376864 : fd_bn254_fp_sqr( &r->X, rr );
235 2376864 : fd_bn254_fp_sub( &r->X, &r->X, j );
236 2376864 : fd_bn254_fp_sub( &r->X, &r->X, v );
237 2376864 : fd_bn254_fp_sub( &r->X, &r->X, v );
238 : /* Y3 = r*(V-X3)-2*Y1*J
239 : note: i no longer used */
240 2376864 : fd_bn254_fp_mul( i, &p->Y, j ); /* i = Y1*J */
241 2376864 : fd_bn254_fp_add( i, i, i ); /* i = 2*Y1*J */
242 2376864 : fd_bn254_fp_sub( &r->Y, v, &r->X );
243 2376864 : fd_bn254_fp_mul( &r->Y, &r->Y, rr );
244 2376864 : fd_bn254_fp_sub( &r->Y, &r->Y, i );
245 : /* Z3 = (Z1+H)^2-Z1Z1-HH */
246 2376864 : fd_bn254_fp_add( &r->Z, &p->Z, h );
247 2376864 : fd_bn254_fp_sqr( &r->Z, &r->Z );
248 2376864 : fd_bn254_fp_sub( &r->Z, &r->Z, zz );
249 2376864 : fd_bn254_fp_sub( &r->Z, &r->Z, hh );
250 2376864 : return r;
251 2376864 : }
252 :
253 : /* fd_bn254_g1_scalar_mul computes r = [s]P.
254 : p must be in affine form (p->Z == 1).
255 : The result is in projective coordinates. */
256 : fd_bn254_g1_t *
257 : fd_bn254_g1_scalar_mul( fd_bn254_g1_t * r,
258 : fd_bn254_g1_t const * p,
259 30126 : fd_bn254_scalar_t const * s ) {
260 30126 : if( FD_UNLIKELY( fd_uint256_is_zero( s ) || fd_bn254_g1_is_zero( p ) ) ) {
261 3 : return fd_bn254_g1_set_zero( r );
262 3 : }
263 30123 : const ulong g1_const[ 3 ] = { 0x5398fd0300ff6565UL, 0x4ccef014a773d2d2UL, 0x0000000000000002UL };
264 30123 : ulong b1[ 3 ];
265 30123 : ulong b2[ 2 ];
266 30123 : fd_bn254_glv_sxg3( b1, s, g1_const );
267 30123 : fd_bn254_glv_sxg2( b2, s, g2_const );
268 :
269 : /* k1 = s - b1*N_A - b2*N_B (always non-negative for G1) */
270 30123 : fd_uint256_t k1[1];
271 30123 : {
272 30123 : ulong p11[ 4 ];
273 : /* b2*nb will produce at most 3 limbs, but we want the 4th zeroed for the addition. */
274 30123 : ulong p21[ 4 ] = {0};
275 30123 : ulong t[ 4 ];
276 30123 : fd_bn254_glv_mul3x2( p11, b1, na );
277 30123 : fd_bn254_glv_mul2x1( p21, b2, nb );
278 30123 : fd_bn254_glv_add4( t, p11, p21 );
279 30123 : fd_bn254_glv_sub4( k1->limbs, s->limbs, t );
280 30123 : }
281 :
282 : /* k2 = b1*N_B - b2*N_C (may be negative) */
283 30123 : fd_uint256_t k2_abs[1];
284 30123 : int k2_neg = 0;
285 30123 : {
286 30123 : ulong pos[ 4 ], neg[ 4 ];
287 30123 : fd_bn254_glv_mul3x1( pos, b1, nb );
288 30123 : fd_bn254_glv_mul2x2( neg, b2, nc );
289 30123 : ulong borrow = fd_bn254_glv_sub4( k2_abs->limbs, pos, neg );
290 30123 : if( borrow ) {
291 0 : k2_neg = 1;
292 0 : fd_bn254_glv_negate4( k2_abs->limbs );
293 0 : }
294 30123 : }
295 :
296 : /* pt2 = phi(P) = (beta * P.x, P.y). If k2 < 0, negate pt2. */
297 30123 : fd_bn254_g1_t pt2[1];
298 30123 : fd_bn254_fp_mul ( &pt2->X, &p->X, fd_bn254_const_beta_mont );
299 30123 : fd_bn254_fp_set ( &pt2->Y, &p->Y );
300 30123 : fd_bn254_fp_set_one( &pt2->Z );
301 30123 : if( k2_neg ) {
302 0 : fd_bn254_fp_neg( &pt2->Y, &pt2->Y );
303 0 : }
304 :
305 30123 : fd_bn254_g1_t pt12[1];
306 30123 : fd_bn254_g1_affine_add( pt12, p, pt2 );
307 :
308 : /* Shamir's trick: simultaneous double-and-add on k1, k2. */
309 30123 : int i = 255;
310 3921069 : for( ; i>=0; i-- ) {
311 3921069 : int k1b = !!fd_uint256_bit( k1, i );
312 3921069 : int k2b = !!fd_uint256_bit( k2_abs, i );
313 3921069 : if( k1b || k2b ) {
314 30123 : fd_bn254_g1_set( r, ( k1b && k2b ) ? pt12 : ( k1b ? p : pt2 ) );
315 30123 : break;
316 30123 : }
317 3921069 : }
318 30123 : if( FD_UNLIKELY( i<0 ) ) {
319 0 : return fd_bn254_g1_set_zero( r );
320 0 : }
321 3820542 : for( i--; i >= 0; i-- ) {
322 3790419 : fd_bn254_g1_dbl( r, r );
323 3790419 : int k1b = !!fd_uint256_bit( k1, i );
324 3790419 : int k2b = !!fd_uint256_bit( k2_abs, i );
325 3790419 : if( k1b && k2b ) {
326 1262370 : fd_bn254_g1_add_mixed( r, r, pt12 );
327 2528049 : } else if( k1b ) {
328 452739 : fd_bn254_g1_add_mixed( r, r, p );
329 2075310 : } else if( k2b ) {
330 661755 : fd_bn254_g1_add_mixed( r, r, pt2 );
331 661755 : }
332 3790419 : }
333 :
334 30123 : return r;
335 30123 : }
336 :
337 : /* fd_bn254_g1_frombytes_internal extracts (x, y) and performs basic checks.
338 : This is used by fd_bn254_g1_compress() and fd_bn254_g1_frombytes_check_subgroup().
339 : https://github.com/arkworks-rs/algebra/blob/v0.4.2/ec/src/models/short_weierstrass/mod.rs#L173-L178 */
340 : fd_bn254_g1_t *
341 : fd_bn254_g1_frombytes_internal( fd_bn254_g1_t * p,
342 : uchar const in[64],
343 121215 : int big_endian ) {
344 : /* Special case: all zeros => point at infinity */
345 121215 : const uchar zero[64] = { 0 };
346 121215 : if( FD_UNLIKELY( fd_memeq( in, zero, 64 ) ) ) {
347 48 : return fd_bn254_g1_set_zero( p );
348 48 : }
349 :
350 : /* Check x < p */
351 121167 : if( FD_UNLIKELY( !fd_bn254_fp_frombytes_nm( &p->X, &in[0], big_endian, NULL, NULL ) ) ) {
352 0 : return NULL;
353 0 : }
354 :
355 : /* Check flags and y < p */
356 121167 : int is_inf, is_neg;
357 121167 : if( FD_UNLIKELY( !fd_bn254_fp_frombytes_nm( &p->Y, &in[32], big_endian, &is_inf, &is_neg ) ) ) {
358 0 : return NULL;
359 0 : }
360 :
361 121167 : if( FD_UNLIKELY( is_inf ) ) {
362 3 : return fd_bn254_g1_set_zero( p );
363 3 : }
364 :
365 121164 : fd_bn254_fp_set_one( &p->Z );
366 121164 : return p;
367 121167 : }
368 :
369 : /* fd_bn254_g1_frombytes_check_subgroup performs frombytes AND checks subgroup membership. */
370 : fd_bn254_g1_t *
371 : fd_bn254_g1_frombytes_check_subgroup( fd_bn254_g1_t * p,
372 : uchar const in[64],
373 91116 : int big_endian ) {
374 91116 : if( FD_UNLIKELY( !fd_bn254_g1_frombytes_internal( p, in, big_endian ) ) ) {
375 0 : return NULL;
376 0 : }
377 91116 : if( FD_UNLIKELY( fd_bn254_g1_is_zero( p ) ) ) {
378 45 : return p;
379 45 : }
380 :
381 91071 : fd_bn254_fp_to_mont( &p->X, &p->X );
382 91071 : fd_bn254_fp_to_mont( &p->Y, &p->Y );
383 91071 : fd_bn254_fp_set_one( &p->Z );
384 :
385 : /* Check that y^2 = x^3 + b */
386 91071 : fd_bn254_fp_t y2[1], x3b[1];
387 91071 : fd_bn254_fp_sqr( y2, &p->Y );
388 91071 : fd_bn254_fp_sqr( x3b, &p->X );
389 91071 : fd_bn254_fp_mul( x3b, x3b, &p->X );
390 91071 : fd_bn254_fp_add( x3b, x3b, fd_bn254_const_b_mont );
391 91071 : if( FD_UNLIKELY( !fd_bn254_fp_eq( y2, x3b ) ) ) {
392 0 : return NULL;
393 0 : }
394 :
395 : /* G1 has prime order, so we don't need to do any further checks. */
396 :
397 91071 : return p;
398 91071 : }
|