field.go raw
1 package gnarl
2
3 // Montgomery field arithmetic over Z_P for the Gnarl signature scheme.
4 //
5 // Representation: 4×uint64 little-endian limbs in Montgomery form.
6 // A value a is stored as aR mod P, where R = 2^256.
7 // Multiplication: montMul(aR, bR) = abR mod P (one R cancels).
8 //
9 // P is a 216-bit prime (27 bytes), but we use 256-bit Montgomery form
10 // internally for uniformity with the CIOS algorithm. The top ~40 bits
11 // of limb[3] are always zero.
12 //
13 // P ≡ 2 mod 3, P mod 5 = 2, 5 is QNR, N = P+1 = 6Q where Q is prime.
14
15 import "math/bits"
16
17 // fe is a field element in Montgomery form: 4 uint64 limbs, little-endian.
18 type fe [4]uint64
19
20 // Prime P limbs (little-endian).
21 // P hex = 0x9563a6b6d81bb9b02e5e5d121e79ccd682cc9931f9791e0f9ee4f5
22 // P bits = 216
23 var pLimbs = fe{
24 0x31f9791e0f9ee4f5,
25 0x121e79ccd682cc99,
26 0xb6d81bb9b02e5e5d,
27 0x00000000009563a6,
28 }
29
30 // R mod P (the Montgomery form of 1).
31 var feOne = fe{
32 0x761ae156c8e8e77c,
33 0x023540a7764c229f,
34 0x315e327188dada36,
35 0x00000000001012bc,
36 }
37
38 // R^2 mod P (used to convert normal form → Montgomery form).
39 var feR2 = fe{
40 0xf56927577991be00,
41 0xf6eb880a87ed3681,
42 0x18d8b12f949488c0,
43 0x0000000000684a62,
44 }
45
46 // feZero is the zero element (Montgomery form of 0 is 0).
47 var feZero fe
48
49 // pPrime = -P^{-1} mod 2^64, the Montgomery reduction constant.
50 const pPrime uint64 = 0xf07d39ef3ea058a3
51
52 // Tonelli-Shanks constants for square root.
53 // P - 1 = 2^s * tsOddQ, where s = 2.
54 const tsS = 2
55
56 // tsOddQ = (P-1) / 4, the odd part of P-1.
57 var tsOddQ = fe{
58 0x4c7e5e4783e7b93d,
59 0x44879e7335a0b326,
60 0xadb606ee6c0b9797,
61 0x00000000002558e9,
62 }
63
64 // tsQm1h = (tsOddQ - 1) / 2, used in Tonelli-Shanks.
65 var tsQm1h = fe{
66 0x263f2f23c1f3dc9e,
67 0xa243cf399ad05993,
68 0xd6db03773605cbcb,
69 0x000000000012ac74,
70 }
71
72 // thirdPNorm = (P-1)/3 in normal (non-Montgomery) form, for y-coordinate
73 // 3-way canonicalization. With P ≡ 2 mod 3, the torus has index 6 = 2×3,
74 // and we canonicalize y to [0, P/3) instead of [0, P/2).
75 var thirdPNorm = fe{
76 0xbb53285f5a8a4c51,
77 0xb0b4d3444780eedd,
78 0x3cf2b3e8900f74c9,
79 0x000000000031cbe2,
80 }
81
82 // halfPNorm = (P-1)/2 in normal (non-Montgomery) form, for alternate checks.
83 var halfPNorm = fe{
84 0x98fcbc8f07cf727a,
85 0x890f3ce66b41664c,
86 0x5b6c0ddcd8172f2e,
87 0x00000000004ab1d3,
88 }
89
90 // pMinus2 = P - 2, used for Fermat inversion.
91 var pMinus2 = fe{
92 0x31f9791e0f9ee4f3,
93 0x121e79ccd682cc99,
94 0xb6d81bb9b02e5e5d,
95 0x00000000009563a6,
96 }
97
98 // feSet copies src into dst.
99 func feSet(dst, src *fe) {
100 *dst = *src
101 }
102
103 // feIsZero returns 1 if a == 0 (constant-time).
104 func feIsZero(a *fe) int {
105 z := a[0] | a[1] | a[2] | a[3]
106 z = (z | (0 - z)) >> 63
107 return 1 - int(z)
108 }
109
110 // feEqual returns 1 if a == b (constant-time).
111 func feEqual(a, b *fe) int {
112 z := (a[0] ^ b[0]) | (a[1] ^ b[1]) | (a[2] ^ b[2]) | (a[3] ^ b[3])
113 z = (z | (0 - z)) >> 63
114 return 1 - int(z)
115 }
116
117 // feReduce reduces a to [0, P) if a >= P. Constant-time.
118 func feReduce(a *fe) {
119 d0, borrow := bits.Sub64(a[0], pLimbs[0], 0)
120 d1, borrow := bits.Sub64(a[1], pLimbs[1], borrow)
121 d2, borrow := bits.Sub64(a[2], pLimbs[2], borrow)
122 d3, borrow := bits.Sub64(a[3], pLimbs[3], borrow)
123
124 mask := uint64(0) - borrow
125 a[0] = (a[0] & mask) | (d0 & ^mask)
126 a[1] = (a[1] & mask) | (d1 & ^mask)
127 a[2] = (a[2] & mask) | (d2 & ^mask)
128 a[3] = (a[3] & mask) | (d3 & ^mask)
129 }
130
131 // feAdd computes r = a + b mod P.
132 func feAdd(r, a, b *fe) {
133 var carry uint64
134 r[0], carry = bits.Add64(a[0], b[0], 0)
135 r[1], carry = bits.Add64(a[1], b[1], carry)
136 r[2], carry = bits.Add64(a[2], b[2], carry)
137 r[3], carry = bits.Add64(a[3], b[3], carry)
138
139 d0, borrow := bits.Sub64(r[0], pLimbs[0], 0)
140 d1, borrow := bits.Sub64(r[1], pLimbs[1], borrow)
141 d2, borrow := bits.Sub64(r[2], pLimbs[2], borrow)
142 d3, borrow := bits.Sub64(r[3], pLimbs[3], borrow)
143
144 underflow := borrow &^ carry
145 mask := uint64(0) - underflow
146 r[0] = (r[0] & mask) | (d0 & ^mask)
147 r[1] = (r[1] & mask) | (d1 & ^mask)
148 r[2] = (r[2] & mask) | (d2 & ^mask)
149 r[3] = (r[3] & mask) | (d3 & ^mask)
150 }
151
152 // feSub computes r = a - b mod P.
153 func feSub(r, a, b *fe) {
154 var borrow uint64
155 r[0], borrow = bits.Sub64(a[0], b[0], 0)
156 r[1], borrow = bits.Sub64(a[1], b[1], borrow)
157 r[2], borrow = bits.Sub64(a[2], b[2], borrow)
158 r[3], borrow = bits.Sub64(a[3], b[3], borrow)
159
160 mask := uint64(0) - borrow
161 var carry uint64
162 r[0], carry = bits.Add64(r[0], pLimbs[0]&mask, 0)
163 r[1], carry = bits.Add64(r[1], pLimbs[1]&mask, carry)
164 r[2], carry = bits.Add64(r[2], pLimbs[2]&mask, carry)
165 r[3], _ = bits.Add64(r[3], pLimbs[3]&mask, carry)
166 }
167
168 // feNeg computes r = -a mod P.
169 func feNeg(r, a *fe) {
170 z := a[0] | a[1] | a[2] | a[3]
171 nonzero := z | (0 - z)
172 mask := 0 - (nonzero >> 63)
173
174 var borrow uint64
175 r[0], borrow = bits.Sub64(pLimbs[0], a[0], 0)
176 r[1], borrow = bits.Sub64(pLimbs[1], a[1], borrow)
177 r[2], borrow = bits.Sub64(pLimbs[2], a[2], borrow)
178 r[3], _ = bits.Sub64(pLimbs[3], a[3], borrow)
179 r[0] &= mask
180 r[1] &= mask
181 r[2] &= mask
182 r[3] &= mask
183 }
184
185 // montMulGeneric computes r = a * b * R^{-1} mod P using fully-unrolled CIOS.
186 //
187 // All carry propagation uses bits.Add64 to avoid silent overflow on hi += c.
188 func montMulGeneric(r, a, b *fe) {
189 var t0, t1, t2, t3, t4 uint64
190 var hi, lo, c0, c1, m, carry uint64
191
192 // ---- i = 0 ----
193 // Part A: t += a[0] * b
194 hi, lo = bits.Mul64(a[0], b[0])
195 t0, c0 = bits.Add64(lo, t0, 0)
196 carry, _ = bits.Add64(hi, 0, c0)
197
198 hi, lo = bits.Mul64(a[0], b[1])
199 lo, c0 = bits.Add64(lo, t1, 0)
200 hi, c1 = bits.Add64(hi, 0, c0)
201 t1, c0 = bits.Add64(lo, carry, 0)
202 carry, _ = bits.Add64(hi, c1, c0)
203
204 hi, lo = bits.Mul64(a[0], b[2])
205 lo, c0 = bits.Add64(lo, t2, 0)
206 hi, c1 = bits.Add64(hi, 0, c0)
207 t2, c0 = bits.Add64(lo, carry, 0)
208 carry, _ = bits.Add64(hi, c1, c0)
209
210 hi, lo = bits.Mul64(a[0], b[3])
211 lo, c0 = bits.Add64(lo, t3, 0)
212 hi, c1 = bits.Add64(hi, 0, c0)
213 t3, c0 = bits.Add64(lo, carry, 0)
214 t4, _ = bits.Add64(hi, c1, c0)
215
216 // Part B: Montgomery reduction
217 m = t0 * pPrime
218 hi, lo = bits.Mul64(m, pLimbs[0])
219 _, c0 = bits.Add64(t0, lo, 0)
220 carry, _ = bits.Add64(hi, 0, c0)
221
222 hi, lo = bits.Mul64(m, pLimbs[1])
223 lo, c0 = bits.Add64(lo, t1, 0)
224 hi, c1 = bits.Add64(hi, 0, c0)
225 t0, c0 = bits.Add64(lo, carry, 0)
226 carry, _ = bits.Add64(hi, c1, c0)
227
228 hi, lo = bits.Mul64(m, pLimbs[2])
229 lo, c0 = bits.Add64(lo, t2, 0)
230 hi, c1 = bits.Add64(hi, 0, c0)
231 t1, c0 = bits.Add64(lo, carry, 0)
232 carry, _ = bits.Add64(hi, c1, c0)
233
234 hi, lo = bits.Mul64(m, pLimbs[3])
235 lo, c0 = bits.Add64(lo, t3, 0)
236 hi, c1 = bits.Add64(hi, 0, c0)
237 t2, c0 = bits.Add64(lo, carry, 0)
238 carry, _ = bits.Add64(hi, c1, c0)
239 t3, c0 = bits.Add64(t4, carry, 0)
240 t4 = c0
241
242 // ---- i = 1 ----
243 hi, lo = bits.Mul64(a[1], b[0])
244 t0, c0 = bits.Add64(lo, t0, 0)
245 carry, _ = bits.Add64(hi, 0, c0)
246
247 hi, lo = bits.Mul64(a[1], b[1])
248 lo, c0 = bits.Add64(lo, t1, 0)
249 hi, c1 = bits.Add64(hi, 0, c0)
250 t1, c0 = bits.Add64(lo, carry, 0)
251 carry, _ = bits.Add64(hi, c1, c0)
252
253 hi, lo = bits.Mul64(a[1], b[2])
254 lo, c0 = bits.Add64(lo, t2, 0)
255 hi, c1 = bits.Add64(hi, 0, c0)
256 t2, c0 = bits.Add64(lo, carry, 0)
257 carry, _ = bits.Add64(hi, c1, c0)
258
259 hi, lo = bits.Mul64(a[1], b[3])
260 lo, c0 = bits.Add64(lo, t3, 0)
261 hi, c1 = bits.Add64(hi, 0, c0)
262 t3, c0 = bits.Add64(lo, carry, 0)
263 hi, _ = bits.Add64(hi, c1, c0)
264 t4 += hi
265
266 m = t0 * pPrime
267 hi, lo = bits.Mul64(m, pLimbs[0])
268 _, c0 = bits.Add64(t0, lo, 0)
269 carry, _ = bits.Add64(hi, 0, c0)
270
271 hi, lo = bits.Mul64(m, pLimbs[1])
272 lo, c0 = bits.Add64(lo, t1, 0)
273 hi, c1 = bits.Add64(hi, 0, c0)
274 t0, c0 = bits.Add64(lo, carry, 0)
275 carry, _ = bits.Add64(hi, c1, c0)
276
277 hi, lo = bits.Mul64(m, pLimbs[2])
278 lo, c0 = bits.Add64(lo, t2, 0)
279 hi, c1 = bits.Add64(hi, 0, c0)
280 t1, c0 = bits.Add64(lo, carry, 0)
281 carry, _ = bits.Add64(hi, c1, c0)
282
283 hi, lo = bits.Mul64(m, pLimbs[3])
284 lo, c0 = bits.Add64(lo, t3, 0)
285 hi, c1 = bits.Add64(hi, 0, c0)
286 t2, c0 = bits.Add64(lo, carry, 0)
287 carry, _ = bits.Add64(hi, c1, c0)
288 t3, c0 = bits.Add64(t4, carry, 0)
289 t4 = c0
290
291 // ---- i = 2 ----
292 hi, lo = bits.Mul64(a[2], b[0])
293 t0, c0 = bits.Add64(lo, t0, 0)
294 carry, _ = bits.Add64(hi, 0, c0)
295
296 hi, lo = bits.Mul64(a[2], b[1])
297 lo, c0 = bits.Add64(lo, t1, 0)
298 hi, c1 = bits.Add64(hi, 0, c0)
299 t1, c0 = bits.Add64(lo, carry, 0)
300 carry, _ = bits.Add64(hi, c1, c0)
301
302 hi, lo = bits.Mul64(a[2], b[2])
303 lo, c0 = bits.Add64(lo, t2, 0)
304 hi, c1 = bits.Add64(hi, 0, c0)
305 t2, c0 = bits.Add64(lo, carry, 0)
306 carry, _ = bits.Add64(hi, c1, c0)
307
308 hi, lo = bits.Mul64(a[2], b[3])
309 lo, c0 = bits.Add64(lo, t3, 0)
310 hi, c1 = bits.Add64(hi, 0, c0)
311 t3, c0 = bits.Add64(lo, carry, 0)
312 hi, _ = bits.Add64(hi, c1, c0)
313 t4 += hi
314
315 m = t0 * pPrime
316 hi, lo = bits.Mul64(m, pLimbs[0])
317 _, c0 = bits.Add64(t0, lo, 0)
318 carry, _ = bits.Add64(hi, 0, c0)
319
320 hi, lo = bits.Mul64(m, pLimbs[1])
321 lo, c0 = bits.Add64(lo, t1, 0)
322 hi, c1 = bits.Add64(hi, 0, c0)
323 t0, c0 = bits.Add64(lo, carry, 0)
324 carry, _ = bits.Add64(hi, c1, c0)
325
326 hi, lo = bits.Mul64(m, pLimbs[2])
327 lo, c0 = bits.Add64(lo, t2, 0)
328 hi, c1 = bits.Add64(hi, 0, c0)
329 t1, c0 = bits.Add64(lo, carry, 0)
330 carry, _ = bits.Add64(hi, c1, c0)
331
332 hi, lo = bits.Mul64(m, pLimbs[3])
333 lo, c0 = bits.Add64(lo, t3, 0)
334 hi, c1 = bits.Add64(hi, 0, c0)
335 t2, c0 = bits.Add64(lo, carry, 0)
336 carry, _ = bits.Add64(hi, c1, c0)
337 t3, c0 = bits.Add64(t4, carry, 0)
338 t4 = c0
339
340 // ---- i = 3 ----
341 hi, lo = bits.Mul64(a[3], b[0])
342 t0, c0 = bits.Add64(lo, t0, 0)
343 carry, _ = bits.Add64(hi, 0, c0)
344
345 hi, lo = bits.Mul64(a[3], b[1])
346 lo, c0 = bits.Add64(lo, t1, 0)
347 hi, c1 = bits.Add64(hi, 0, c0)
348 t1, c0 = bits.Add64(lo, carry, 0)
349 carry, _ = bits.Add64(hi, c1, c0)
350
351 hi, lo = bits.Mul64(a[3], b[2])
352 lo, c0 = bits.Add64(lo, t2, 0)
353 hi, c1 = bits.Add64(hi, 0, c0)
354 t2, c0 = bits.Add64(lo, carry, 0)
355 carry, _ = bits.Add64(hi, c1, c0)
356
357 hi, lo = bits.Mul64(a[3], b[3])
358 lo, c0 = bits.Add64(lo, t3, 0)
359 hi, c1 = bits.Add64(hi, 0, c0)
360 t3, c0 = bits.Add64(lo, carry, 0)
361 hi, _ = bits.Add64(hi, c1, c0)
362 t4 += hi
363
364 m = t0 * pPrime
365 hi, lo = bits.Mul64(m, pLimbs[0])
366 _, c0 = bits.Add64(t0, lo, 0)
367 carry, _ = bits.Add64(hi, 0, c0)
368
369 hi, lo = bits.Mul64(m, pLimbs[1])
370 lo, c0 = bits.Add64(lo, t1, 0)
371 hi, c1 = bits.Add64(hi, 0, c0)
372 t0, c0 = bits.Add64(lo, carry, 0)
373 carry, _ = bits.Add64(hi, c1, c0)
374
375 hi, lo = bits.Mul64(m, pLimbs[2])
376 lo, c0 = bits.Add64(lo, t2, 0)
377 hi, c1 = bits.Add64(hi, 0, c0)
378 t1, c0 = bits.Add64(lo, carry, 0)
379 carry, _ = bits.Add64(hi, c1, c0)
380
381 hi, lo = bits.Mul64(m, pLimbs[3])
382 lo, c0 = bits.Add64(lo, t3, 0)
383 hi, c1 = bits.Add64(hi, 0, c0)
384 t2, c0 = bits.Add64(lo, carry, 0)
385 carry, _ = bits.Add64(hi, c1, c0)
386 t3, c0 = bits.Add64(t4, carry, 0)
387 t4 = c0
388
389 r[0] = t0
390 r[1] = t1
391 r[2] = t2
392 r[3] = t3
393 feReduce(r)
394 }
395
396 // montSquareGeneric computes r = a^2 * R^{-1} mod P.
397 func montSquareGeneric(r, a *fe) {
398 montMulGeneric(r, a, a)
399 }
400
401 // feToMont converts a normal-form value to Montgomery form: aR mod P.
402 func feToMont(r, a *fe) {
403 montMul(r, a, &feR2)
404 }
405
406 // feFromMont converts a Montgomery-form value to normal form: a * R^{-1} mod P.
407 func feFromMont(r, a *fe) {
408 one := fe{1, 0, 0, 0}
409 montMul(r, a, &one)
410 }
411
412 // feFromSmall converts a small integer to Montgomery form.
413 func feFromSmall(r *fe, v int64) {
414 if v >= 0 {
415 *r = fe{uint64(v), 0, 0, 0}
416 } else {
417 *r = pLimbs
418 var borrow uint64
419 r[0], borrow = bits.Sub64(r[0], uint64(-v), 0)
420 r[1], borrow = bits.Sub64(r[1], 0, borrow)
421 r[2], borrow = bits.Sub64(r[2], 0, borrow)
422 r[3], _ = bits.Sub64(r[3], 0, borrow)
423 }
424 feToMont(r, r)
425 }
426
427 // fieldLen is the byte length of a serialized field element (27 bytes for 216-bit P).
428 const fieldLen = 27
429
430 // feFromBytes27 sets r from a 27-byte big-endian encoding and converts to Montgomery form.
431 // Returns false if the value >= P.
432 func feFromBytes27(r *fe, b []byte) bool {
433 if len(b) < 27 {
434 *r = feZero
435 return false
436 }
437 // 27 bytes = 216 bits. Big-endian layout:
438 // b[0..2] → low 24 bits of limb[3] (bits 192..215)
439 // b[3..10] → limb[2] (bits 128..191)
440 // b[11..18] → limb[1] (bits 64..127)
441 // b[19..26] → limb[0] (bits 0..63)
442 r[3] = uint64(b[0])<<16 | uint64(b[1])<<8 | uint64(b[2])
443 r[2] = uint64(b[3])<<56 | uint64(b[4])<<48 | uint64(b[5])<<40 | uint64(b[6])<<32 |
444 uint64(b[7])<<24 | uint64(b[8])<<16 | uint64(b[9])<<8 | uint64(b[10])
445 r[1] = uint64(b[11])<<56 | uint64(b[12])<<48 | uint64(b[13])<<40 | uint64(b[14])<<32 |
446 uint64(b[15])<<24 | uint64(b[16])<<16 | uint64(b[17])<<8 | uint64(b[18])
447 r[0] = uint64(b[19])<<56 | uint64(b[20])<<48 | uint64(b[21])<<40 | uint64(b[22])<<32 |
448 uint64(b[23])<<24 | uint64(b[24])<<16 | uint64(b[25])<<8 | uint64(b[26])
449
450 // Check >= P.
451 d0, borrow := bits.Sub64(r[0], pLimbs[0], 0)
452 d1, borrow := bits.Sub64(r[1], pLimbs[1], borrow)
453 _, borrow = bits.Sub64(r[2], pLimbs[2], borrow)
454 _, borrow = bits.Sub64(r[3], pLimbs[3], borrow)
455 _ = d0
456 _ = d1
457 if borrow == 0 {
458 *r = feZero
459 return false
460 }
461
462 feToMont(r, r)
463 return true
464 }
465
466 // feToBytes27 writes a Montgomery-form element to a 27-byte big-endian buffer.
467 func feToBytes27(b []byte, a *fe) {
468 var norm fe
469 feFromMont(&norm, a)
470 // limb[3] low 24 bits → b[0..2]
471 b[0] = byte(norm[3] >> 16)
472 b[1] = byte(norm[3] >> 8)
473 b[2] = byte(norm[3])
474 // limb[2] → b[3..10]
475 b[3] = byte(norm[2] >> 56)
476 b[4] = byte(norm[2] >> 48)
477 b[5] = byte(norm[2] >> 40)
478 b[6] = byte(norm[2] >> 32)
479 b[7] = byte(norm[2] >> 24)
480 b[8] = byte(norm[2] >> 16)
481 b[9] = byte(norm[2] >> 8)
482 b[10] = byte(norm[2])
483 // limb[1] → b[11..18]
484 b[11] = byte(norm[1] >> 56)
485 b[12] = byte(norm[1] >> 48)
486 b[13] = byte(norm[1] >> 40)
487 b[14] = byte(norm[1] >> 32)
488 b[15] = byte(norm[1] >> 24)
489 b[16] = byte(norm[1] >> 16)
490 b[17] = byte(norm[1] >> 8)
491 b[18] = byte(norm[1])
492 // limb[0] → b[19..26]
493 b[19] = byte(norm[0] >> 56)
494 b[20] = byte(norm[0] >> 48)
495 b[21] = byte(norm[0] >> 40)
496 b[22] = byte(norm[0] >> 32)
497 b[23] = byte(norm[0] >> 24)
498 b[24] = byte(norm[0] >> 16)
499 b[25] = byte(norm[0] >> 8)
500 b[26] = byte(norm[0])
501 }
502
503 // fePow computes r = base^exp where exp is a 256-bit exponent in little-endian limbs
504 // (plain integer, NOT Montgomery). base is in Montgomery form, result in Montgomery form.
505 func fePow(r, base *fe, exp *fe) {
506 var result fe
507 feSet(&result, &feOne)
508 var b fe
509 feSet(&b, base)
510
511 for i := 0; i < 4; i++ {
512 word := exp[i]
513 for bit := 0; bit < 64; bit++ {
514 if word&1 == 1 {
515 montMul(&result, &result, &b)
516 }
517 montSquare(&b, &b)
518 word >>= 1
519 }
520 }
521 feSet(r, &result)
522 }
523
524 // feInv computes r = a^{-1} in Montgomery form via Fermat: a^{P-2}.
525 func feInv(r, a *fe) {
526 fePow(r, a, &pMinus2)
527 }
528
529 // feSqrt computes r = sqrt(a) mod P in Montgomery form using Tonelli-Shanks.
530 // Returns false if a is not a quadratic residue.
531 //
532 // P - 1 = 2^2 * q where q is odd (tsS = 2).
533 // Non-residue: 5 (QNR by construction of P).
534 func feSqrt(r, a *fe) bool {
535 if feIsZero(a) == 1 {
536 *r = feZero
537 return true
538 }
539
540 var feFive fe
541 feFive[0] = 5
542 feToMont(&feFive, &feFive)
543
544 // t = a^q
545 var t fe
546 fePow(&t, a, &tsOddQ)
547
548 // rr = a^{(q+1)/2} = a * a^{(q-1)/2}
549 var rr, tmp fe
550 fePow(&tmp, a, &tsQm1h)
551 montMul(&rr, &tmp, a)
552
553 // c = 5^q
554 var c fe
555 fePow(&c, &feFive, &tsOddQ)
556
557 m := tsS
558
559 for {
560 if feEqual(&t, &feOne) == 1 {
561 feSet(r, &rr)
562 return true
563 }
564
565 var b fe
566 feSet(&b, &t)
567 i := 0
568 for j := 1; j < m; j++ {
569 montSquare(&b, &b)
570 if feEqual(&b, &feOne) == 1 {
571 i = j
572 break
573 }
574 }
575 if i == 0 {
576 return false
577 }
578
579 feSet(&b, &c)
580 for j := 0; j < m-i-1; j++ {
581 montSquare(&b, &b)
582 }
583
584 montMul(&rr, &rr, &b)
585 montSquare(&c, &b)
586 montMul(&t, &t, &c)
587 m = i
588 }
589 }
590
591 // feIsGtThirdP returns 1 if a (in NORMAL form, not Montgomery) is > (P-1)/3.
592 // Used for 3-way canonicalization of torus y-coordinate.
593 func feIsGtThirdP(a *fe) int {
594 for i := 3; i >= 0; i-- {
595 if a[i] > thirdPNorm[i] {
596 return 1
597 }
598 if a[i] < thirdPNorm[i] {
599 return 0
600 }
601 }
602 return 0 // equal → not strictly greater
603 }
604
605 // feIsGtHalfP returns 1 if a (in NORMAL form, not Montgomery) is > (P-1)/2.
606 func feIsGtHalfP(a *fe) int {
607 for i := 3; i >= 0; i-- {
608 if a[i] > halfPNorm[i] {
609 return 1
610 }
611 if a[i] < halfPNorm[i] {
612 return 0
613 }
614 }
615 return 0
616 }
617