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 : lambda = 4407920970296243842393367215006156084916469457145843978461
34 : (a root of X^2 + X + 1 = 0 mod r)
35 :
36 : Computing phi(P) = (beta * x, y) costs just one FP multiplication, but
37 : is equivalent to a full scalar multiplication by the ~254-bit scalar lambda.
38 :
39 : We want to "decompose" the input scalar in a way where a majority of the
40 : effort is re-used by the lambda computation, leaving the remaining operations
41 : smaller and faster to perform. Given a scalar s, we want to find small
42 : k1, k2 (each ~128-bit) such that:
43 : s = k1 + k2 * lambda (mod r)
44 : Then: [s]P = [k1]P + [k2](lambda * P) = [k1]P + [k2]phi(P).
45 :
46 : Using the "Straus-Shamir trick", https://pmc.ncbi.nlm.nih.gov/articles/PMC9028562/#app1-sensors-22-03083,
47 : this requires only ~128 doublings + ~128 additions instead of ~256 doublings.
48 :
49 : To decompose the scalar, we have a 2-dimensional lattice
50 : L = { (a, b) in Z^2 : a + b*lambda = 0 (mod r) }
51 :
52 : A reduced basis of L is given by three magnitudes:
53 : N_A = 147946756881789319000765030803803410728 (~128 bits, 2 limbs)
54 : N_B = 9931322734385697763 (~64 bits, 1 limb)
55 : N_C = 147946756881789319010696353538189108491 (~128 bits, 2 limbs)
56 :
57 : These are arranged differently depending on the group.
58 : G1:
59 : | +N_A +N_B |
60 : | -N_B +N_C |
61 : G2:
62 : | -N_C -N_B |
63 : | +N_B -N_A |
64 :
65 : We avoid big integer divisions by r using Babai's algorithm with
66 : precomputed fixed-point inverses:
67 : b1 = (s * g1) >> 256, where for group:
68 : G1: g1 = round(2^256 * N_C / r)
69 : G2: g1 = round(2^256 * N_A / r)
70 : b2 = (s * g2) >> 256 g2 = round(2^256 * N_B / r), 66-bit, 2 limbs)
71 :
72 : Then:
73 : G1: k1 = s - b1*N_A - b2*N_B, k2 = b1*N_B - b2*N_C
74 : G2: k1 = s - b1*N_C - b2*N_B, k2 = b2*N_A - b1*N_B
75 :
76 : For G1, k1 >= 0 always, k2 may be negative.
77 : For G2, both k1 and k2 may be negative. */
78 :
79 : /* beta (defined in fd_bn254_g1.c) in Montgomery form.
80 : 0x30644e72e131a0295e6dd9e7e0acccb0c28f069fbb966e3de4bd44e5607cfd48 */
81 : extern const fd_bn254_fp_t fd_bn254_const_beta_mont[1];
82 :
83 : /* Lattice constants, see glv.py. Defined in fd_bn254_field.c. */
84 : extern const ulong na[ 2 ];
85 : extern const ulong nb[ 1 ];
86 : extern const ulong nc[ 2 ];
87 :
88 : /* g1 = round(2^256 * N_C / r) for G1, round(2^256 * N_A / r) for G2
89 : (3 limbs). */
90 : extern const ulong g1_const_g1[ 3 ];
91 : extern const ulong g1_const_g2[ 3 ];
92 :
93 : /* g2 = round(2^256 * N_B / r), 66-bit (2 limbs). Same for G1 and G2. */
94 : extern const ulong g2_const[ 2 ];
95 :
96 : /* Multiply 4-limb scalar s by a 3-limb constant g.
97 : Returns top 3 limbs. */
98 : static inline void
99 : fd_bn254_glv_sxg3( ulong out[ 3 ],
100 : fd_bn254_scalar_t const * s,
101 776 : ulong const g[ 3 ] ) {
102 776 : uint128 s0 = s->limbs[0];
103 776 : uint128 s1 = s->limbs[1];
104 776 : uint128 s2 = s->limbs[2];
105 776 : uint128 s3 = s->limbs[3];
106 776 : uint128 acc;
107 776 : acc = s0 * g[ 0 ];
108 776 : acc = s1 * g[ 0 ] + s0 * g[ 1 ] + (ulong)(acc >> 64);
109 776 : acc = s2 * g[ 0 ] + s1 * g[ 1 ] + s0 * g[ 2 ] + (ulong)(acc >> 64);
110 776 : acc = s3 * g[ 0 ] + s2 * g[ 1 ] + s1 * g[ 2 ] + (ulong)(acc >> 64);
111 776 : acc = s3 * g[ 1 ] + s2 * g[ 2 ] + (ulong)(acc >> 64); out[ 0 ] = (ulong)acc;
112 776 : acc = s3 * g[ 2 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
113 776 : acc = (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
114 776 : }
115 :
116 : /* Same, but for a 2-limb constant g. */
117 : static inline void
118 : fd_bn254_glv_sxg2( ulong out[ 2 ],
119 : fd_bn254_scalar_t const * s,
120 776 : ulong const g[ 2 ] ) {
121 776 : uint128 s0 = s->limbs[0];
122 776 : uint128 s1 = s->limbs[1];
123 776 : uint128 s2 = s->limbs[2];
124 776 : uint128 s3 = s->limbs[3];
125 776 : uint128 acc;
126 776 : acc = s0 * g[ 0 ];
127 776 : acc = s1 * g[ 0 ] + s0 * g[ 1 ] + (ulong)(acc >> 64);
128 776 : acc = s2 * g[ 0 ] + s1 * g[ 1 ] + (ulong)(acc >> 64);
129 776 : acc = s3 * g[ 0 ] + s2 * g[ 1 ] + (ulong)(acc >> 64);
130 776 : acc = s3 * g[ 1 ] + (ulong)(acc >> 64); out[ 0 ] = (ulong)acc;
131 776 : acc = (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
132 776 : }
133 :
134 : /* Multiply 3-limb a by 2-limb n, store low 4 limbs into out. */
135 : static inline void
136 : fd_bn254_glv_mul3x2( ulong out[ 4 ],
137 : ulong const a[ 3 ],
138 776 : ulong const n[ 2 ] ) {
139 776 : uint128 acc;
140 776 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
141 776 : acc = (uint128)a[ 1 ] * n[ 0 ] + (uint128)a[ 0 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
142 776 : acc = (uint128)a[ 2 ] * n[ 0 ] + (uint128)a[ 1 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
143 776 : acc = (uint128)a[ 2 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 3 ] = (ulong)acc;
144 776 : }
145 :
146 : /* Multiply 3-limb by 1-limb, store into 4-limb. */
147 : static inline void
148 : fd_bn254_glv_mul3x1( ulong out[ 4 ],
149 : ulong const a[ 3 ],
150 776 : ulong const n[ 1 ] ) {
151 776 : uint128 acc;
152 776 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
153 776 : acc = (uint128)a[ 1 ] * n[ 0 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
154 776 : acc = (uint128)a[ 2 ] * n[ 0 ] + (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
155 776 : acc = (ulong)(acc >> 64); out[ 3 ] = (ulong)acc;
156 776 : }
157 :
158 : /* Multiply 2-limb by 2-limb, store into 4-limb. */
159 : static inline void
160 : fd_bn254_glv_mul2x2( ulong out[ 4 ],
161 : ulong const a[ 2 ],
162 776 : ulong const n[ 2 ] ) {
163 776 : uint128 acc;
164 776 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
165 776 : acc = (uint128)a[ 1 ] * n[ 0 ] + (uint128)a[ 0 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
166 776 : acc = (uint128)a[ 1 ] * n[ 1 ] + (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
167 776 : acc = (ulong)(acc >> 64); out[ 3 ] = (ulong)acc;
168 776 : }
169 :
170 : /* Multiply 2-limb by 1-limb, stores 3-limbs. */
171 : static inline void
172 : fd_bn254_glv_mul2x1( ulong out[ 3 ],
173 : ulong const a[ 2 ],
174 776 : ulong const n[ 1 ] ) {
175 776 : uint128 acc;
176 776 : acc = (uint128)a[ 0 ] * n[ 0 ]; out[ 0 ] = (ulong)acc;
177 776 : acc = (uint128)a[ 1 ] * n[ 0 ] + (ulong)(acc >> 64); out[ 1 ] = (ulong)acc;
178 776 : acc = (ulong)(acc >> 64); out[ 2 ] = (ulong)acc;
179 776 : }
180 :
181 : /* 4-limb addition: out = a + b. Returns carry. */
182 : static inline ulong
183 : fd_bn254_glv_add4( ulong out[ 4 ],
184 : ulong const a[ 4 ],
185 776 : ulong const b[ 4 ] ) {
186 776 : ulong carry = 0;
187 3880 : for( int j = 0; j < 4; j++ ) {
188 3104 : uint128 acc = (uint128)a[ j ] + b[ j ] + carry;
189 3104 : out[ j ] = (ulong)acc;
190 3104 : carry = (ulong)(acc >> 64);
191 3104 : }
192 776 : return carry;
193 776 : }
194 :
195 : /* 4-limb subtraction: out = a - b. Returns borrow (1 if a < b). */
196 : static inline ulong
197 : fd_bn254_glv_sub4( ulong out[ 4 ],
198 : ulong const a[ 4 ],
199 1551 : ulong const b[ 4 ] ) {
200 1551 : ulong borrow = 0;
201 7755 : for( int j = 0; j < 4; j++ ) {
202 6204 : ulong av = a[ j ];
203 6204 : ulong bv = b[ j ];
204 6204 : ulong diff = av - bv - borrow;
205 6204 : borrow = ( av < bv || (borrow && av == bv) ) ? 1UL : 0UL;
206 6204 : out[ j ] = diff;
207 6204 : }
208 1551 : return borrow;
209 1551 : }
210 :
211 : /* Two's complement negation of a 4-limb value in-place. */
212 : static inline void
213 28 : fd_bn254_glv_negate4( ulong v[ 4 ] ) {
214 28 : ulong carry = 1;
215 140 : for( int j = 0; j < 4; j++ ) {
216 112 : uint128 sum = (uint128)(~v[ j ]) + carry;
217 112 : v[ j ] = (ulong)sum;
218 112 : carry = (ulong)(sum >> 64);
219 112 : }
220 28 : }
221 :
222 : #endif /* HEADER_fd_src_ballet_bn254_fd_bn254_glv_h */
|
