Line data Source code
1 : #ifndef HEADER_fd_src_ballet_bn254_fd_bn254_glv_h
2 : #define HEADER_fd_src_ballet_bn254_fd_bn254_glv_h
3 :
4 : /* Included by fd_bn254_g1.c and fd_bn254_g2.c, should not be used elsewhere. */
5 :
6 : #include "./fd_bn254_internal.h"
7 :
8 : /* Rundown of the BN254 GLV implementation:
9 :
10 : Context: https://bitcointalk.org/index.php?topic=3238.msg45565#msg45565
11 :
12 : BN254 is a "pairing-friendly" curve E: y^2 = x^3 + 3 over a prime field
13 : Fp, where:
14 : p = 21888242871839275222246405745257275088696311157297823662689037894645226208583
15 :
16 : The group G1 = E(Fp) has prime order:
17 : r = 21888242871839275222246405745257275088548364400416034343698204186575808495617
18 :
19 : A naive approach to scalar multiplication, [s]P, requires ~256 doublings
20 : and ~128 additions (double-and-add) in the worst case. The GLV method
21 : takes advantage of an easy to compute endomorphism to cut this roughly
22 : in half.
23 :
24 : For BN254, there exists a cube root of unity beta in Fp such that:
25 : phi: (x, y) -> (beta * x, y)
26 : is an endomorphism of E.
27 : That is, if point P = (x, y) is on the curve, then phi(P) = (beta * x, y)
28 : is also on the curve, since (beta*x)^3 = beta^3 * x^3) = x^3,
29 : so the curve equation is preserved.
30 :
31 : phi satisfies phi(P) = [lambda]P where lambda is a cube root of unity
32 : modulo r (the group order we defined above^).
33 :
34 : Computing phi(P) = (beta * x, y) costs just one FP multiplication, but
35 : is equivalent to a full scalar multiplication by the scalar lambda.
36 :
37 : We want to "decompose" the input scalar in a way where a majority of the
38 : effort is re-used by the lambda computation, leaving the remaining operations
39 : smaller and faster to perform. Given a scalar s, we want to find small
40 : k1, k2 (each ~128-bit) such that:
41 : s = k1 + k2 * lambda (mod r)
42 : Then: [s]P = [k1]P + [k2](lambda * P) = [k1]P + [k2]phi(P).
43 :
44 : Using the "Straus-Shamir trick", https://pmc.ncbi.nlm.nih.gov/articles/PMC9028562/#app1-sensors-22-03083,
45 : this requires only ~128 doublings + ~128 additions instead of ~256 doublings.
46 :
47 : To decompose the scalar, we have a 2-dimensional lattice
48 : L = { (a, b) in Z^2 : a + b*lambda = 0 (mod r) }
49 :
50 : A reduced basis of L is given by three magnitudes:
51 : N_A = 147946756881789319000765030803803410728 (~128 bits, 2 limbs)
52 : N_B = 9931322734385697763 (~64 bits, 1 limb)
53 : N_C = 147946756881789319010696353538189108491 (~128 bits, 2 limbs)
54 :
55 : These are arranged differently depending on the group.
56 : G1:
57 : | +N_A +N_B |
58 : | -N_B +N_C |
59 : G2:
60 : | -N_C -N_B |
61 : | +N_B -N_A |
62 :
63 : We avoid big integer divisions by r using Babai's algorithm with
64 : precomputed fixed-point inverses:
65 : b1 = (s * g1) >> 256, where for group:
66 : G1: g1 = round(2^256 * N_C / r)
67 : G2: g1 = round(2^256 * N_A / r)
68 : b2 = (s * g2) >> 256 g2 = round(2^256 * N_B / r), 66-bit, 2 limbs)
69 :
70 : Then:
71 : G1: k1 = s - b1*N_A - b2*N_B, k2 = b1*N_B - b2*N_C
72 : G2: k1 = s - b1*N_C - b2*N_B, k2 = b2*N_A - b1*N_B
73 :
74 : For G1, k1 >= 0 always, k2 may be negative.
75 : For G2, both k1 and k2 may be negative. */
76 :
77 : /* beta (defined in fd_bn254_g1.c) in Montgomery form.
78 : 0x30644e72e131a0295e6dd9e7e0acccb0c28f069fbb966e3de4bd44e5607cfd48 */
79 : extern const fd_bn254_fp_t fd_bn254_const_beta_mont[1];
80 :
81 : /* Lattice constants, see glv.py. Defined in fd_bn254_field.c. */
82 : extern const ulong na[ 2 ];
83 : extern const ulong nb[ 1 ];
84 : extern const ulong nc[ 2 ];
85 :
86 : /* g1 = round(2^256 * N_C / r) for G1, round(2^256 * N_A / r) for G2
87 : (3 limbs). */
88 : extern const ulong g1_const_g1[ 3 ];
89 : extern const ulong g1_const_g2[ 3 ];
90 :
91 : /* g2 = round(2^256 * N_B / r), 66-bit (2 limbs). Same for G1 and G2. */
92 : extern const ulong g2_const[ 2 ];
93 :
94 : /* Multiply 4-limb scalar s by a 3-limb constant g.
95 : Returns top 3 limbs. */
96 : static inline void
97 : fd_bn254_glv_sxg3( ulong out[ 3 ],
98 : fd_bn254_scalar_t const * s,
99 776 : ulong const g[ 3 ] ) {
100 776 : uint128 s0 = s->limbs[0];
101 776 : uint128 s1 = s->limbs[1];
102 776 : uint128 s2 = s->limbs[2];
103 776 : uint128 s3 = s->limbs[3];
104 776 : uint128 acc;
105 776 : acc = s0 * g[ 0 ];
106 776 : acc = s1 * g[ 0 ] + s0 * g[ 1 ] + (ulong)(acc >> 64);
107 776 : acc = s2 * g[ 0 ] + s1 * g[ 1 ] + s0 * g[ 2 ] + (ulong)(acc >> 64);
108 776 : acc = s3 * g[ 0 ] + s2 * g[ 1 ] + s1 * g[ 2 ] + (ulong)(acc >> 64);
109 776 : acc = s3 * g[ 1 ] + s2 * g[ 2 ] + (ulong)(acc >> 64); out[ 0 ] = (ulong)acc;
110 776 : acc = s3 * g[ 2 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
111 776 : acc = (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
112 776 : }
113 :
114 : /* Same, but for a 2-limb constant g. */
115 : static inline void
116 : fd_bn254_glv_sxg2( ulong out[ 2 ],
117 : fd_bn254_scalar_t const * s,
118 776 : ulong const g[ 2 ] ) {
119 776 : uint128 s0 = s->limbs[0];
120 776 : uint128 s1 = s->limbs[1];
121 776 : uint128 s2 = s->limbs[2];
122 776 : uint128 s3 = s->limbs[3];
123 776 : uint128 acc;
124 776 : acc = s0 * g[ 0 ];
125 776 : acc = s1 * g[ 0 ] + s0 * g[ 1 ] + (ulong)(acc >> 64);
126 776 : acc = s2 * g[ 0 ] + s1 * g[ 1 ] + (ulong)(acc >> 64);
127 776 : acc = s3 * g[ 0 ] + s2 * g[ 1 ] + (ulong)(acc >> 64);
128 776 : acc = s3 * g[ 1 ] + (ulong)(acc >> 64); out[ 0 ] = (ulong)acc;
129 776 : acc = (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
130 776 : }
131 :
132 : /* Multiply 3-limb a by 2-limb n, store low 4 limbs into out. */
133 : static inline void
134 : fd_bn254_glv_mul3x2( ulong out[ 4 ],
135 : ulong const a[ 3 ],
136 775 : ulong const n[ 2 ] ) {
137 775 : uint128 acc;
138 775 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
139 775 : acc = (uint128)a[ 1 ] * n[ 0 ] + (uint128)a[ 0 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
140 775 : acc = (uint128)a[ 2 ] * n[ 0 ] + (uint128)a[ 1 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
141 775 : acc = (uint128)a[ 2 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 3 ] = (ulong)acc;
142 775 : }
143 :
144 : /* Multiply 3-limb by 1-limb, store into 4-limb. */
145 : static inline void
146 : fd_bn254_glv_mul3x1( ulong out[ 4 ],
147 : ulong const a[ 3 ],
148 775 : ulong const n[ 1 ] ) {
149 775 : uint128 acc;
150 775 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
151 775 : acc = (uint128)a[ 1 ] * n[ 0 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
152 775 : acc = (uint128)a[ 2 ] * n[ 0 ] + (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
153 775 : acc = (ulong)(acc >> 64); out[ 3 ] = (ulong)acc;
154 775 : }
155 :
156 : /* Multiply 2-limb by 2-limb, store into 4-limb. */
157 : static inline void
158 : fd_bn254_glv_mul2x2( ulong out[ 4 ],
159 : ulong const a[ 2 ],
160 775 : ulong const n[ 2 ] ) {
161 775 : uint128 acc;
162 775 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
163 775 : acc = (uint128)a[ 1 ] * n[ 0 ] + (uint128)a[ 0 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
164 775 : acc = (uint128)a[ 1 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
165 775 : acc = (ulong)(acc >> 64); out[ 3 ] = (ulong)acc;
166 775 : }
167 :
168 : /* Multiply 2-limb by 1-limb, stores 3-limbs. */
169 : static inline void
170 : fd_bn254_glv_mul2x1( ulong out[ 3 ],
171 : ulong const a[ 2 ],
172 775 : ulong const n[ 1 ] ) {
173 775 : uint128 acc;
174 775 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
175 775 : acc = (uint128)a[ 1 ] * n[ 0 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
176 775 : acc = (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
177 775 : }
178 :
179 : /* 4-limb addition: out = a + b. Returns carry. */
180 : static inline ulong
181 : fd_bn254_glv_add4( ulong out[ 4 ],
182 : ulong const a[ 4 ],
183 775 : ulong const b[ 4 ] ) {
184 775 : ulong carry = 0;
185 3875 : for( int j = 0; j < 4; j++ ) {
186 3100 : uint128 acc = (uint128)a[ j ] + b[ j ] + carry;
187 3100 : out[ j ] = (ulong)acc;
188 3100 : carry = (ulong)(acc >> 64);
189 3100 : }
190 775 : return carry;
191 775 : }
192 :
193 : /* 4-limb subtraction: out = a - b. Returns borrow (1 if a < b). */
194 : static inline ulong
195 : fd_bn254_glv_sub4( ulong out[ 4 ],
196 : ulong const a[ 4 ],
197 1550 : ulong const b[ 4 ] ) {
198 1550 : ulong borrow = 0;
199 7750 : for( int j = 0; j < 4; j++ ) {
200 6200 : ulong av = a[ j ];
201 6200 : ulong bv = b[ j ];
202 6200 : ulong diff = av - bv - borrow;
203 6200 : borrow = ( av < bv || (borrow && av == bv) ) ? 1UL : 0UL;
204 6200 : out[ j ] = diff;
205 6200 : }
206 1550 : return borrow;
207 1550 : }
208 :
209 : /* Two's complement negation of a 4-limb value in-place. */
210 : static inline void
211 28 : fd_bn254_glv_negate4( ulong v[ 4 ] ) {
212 28 : ulong carry = 1;
213 140 : for( int j = 0; j < 4; j++ ) {
214 112 : uint128 sum = (uint128)(~v[ j ]) + carry;
215 112 : v[ j ] = (ulong)sum;
216 112 : carry = (ulong)(sum >> 64);
217 112 : }
218 28 : }
219 :
220 : #endif /* HEADER_fd_src_ballet_bn254_fd_bn254_glv_h */
|
