1// Copyright 2020 ConsenSys AG
  2//
  3// Licensed under the Apache License, Version 2.0 (the "License");
  4// you may not use this file except in compliance with the License.
  5// You may obtain a copy of the License at
  6//
  7//     http://www.apache.org/licenses/LICENSE-2.0
  8//
  9// Unless required by applicable law or agreed to in writing, software
 10// distributed under the License is distributed on an "AS IS" BASIS,
 11// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 12// See the License for the specific language governing permissions and
 13// limitations under the License.
 14
 15// field modulus q =
 16//
 17// 21888242871839275222246405745257275088548364400416034343698204186575808495617
 18// Code generated by goff DO NOT EDIT
 19// goff version:  - build:
 20// Element are assumed to be in Montgomery form in all methods
 21
 22// Package ff (generated by goff) contains field arithmetics operations
 23package ff
 24
 25import (
 26	"crypto/rand"
 27	"encoding/binary"
 28	"io"
 29	"math/big"
 30	"math/bits"
 31	"sync"
 32
 33	"unsafe"
 34)
 35
 36// Element represents a field element stored on 4 words (uint64)
 37// Element are assumed to be in Montgomery form in all methods
 38type Element [4]uint64
 39
 40// ElementLimbs number of 64 bits words needed to represent Element
 41const ElementLimbs = 4
 42
 43// ElementBits number bits needed to represent Element
 44const ElementBits = 254
 45
 46// SetUint64 z = v, sets z LSB to v (non-Montgomery form) and convert z to Montgomery form
 47func (z *Element) SetUint64(v uint64) *Element {
 48	z[0] = v
 49	z[1] = 0
 50	z[2] = 0
 51	z[3] = 0
 52	return z.ToMont()
 53}
 54
 55// Set z = x
 56func (z *Element) Set(x *Element) *Element {
 57	z[0] = x[0]
 58	z[1] = x[1]
 59	z[2] = x[2]
 60	z[3] = x[3]
 61	return z
 62}
 63
 64// SetZero z = 0
 65func (z *Element) SetZero() *Element {
 66	z[0] = 0
 67	z[1] = 0
 68	z[2] = 0
 69	z[3] = 0
 70	return z
 71}
 72
 73// SetOne z = 1 (in Montgomery form)
 74func (z *Element) SetOne() *Element {
 75	z[0] = 12436184717236109307
 76	z[1] = 3962172157175319849
 77	z[2] = 7381016538464732718
 78	z[3] = 1011752739694698287
 79	return z
 80}
 81
 82// Neg z = q - x
 83func (z *Element) Neg(x *Element) *Element {
 84	if x.IsZero() {
 85		return z.SetZero()
 86	}
 87	var borrow uint64
 88	z[0], borrow = bits.Sub64(4891460686036598785, x[0], 0)
 89	z[1], borrow = bits.Sub64(2896914383306846353, x[1], borrow)
 90	z[2], borrow = bits.Sub64(13281191951274694749, x[2], borrow)
 91	z[3], _ = bits.Sub64(3486998266802970665, x[3], borrow)
 92	return z
 93}
 94
 95// Div z = x*y^-1 mod q
 96func (z *Element) Div(x, y *Element) *Element {
 97	var yInv Element
 98	yInv.Inverse(y)
 99	z.Mul(x, &yInv)
100	return z
101}
102
103// Equal returns z == x
104func (z *Element) Equal(x *Element) bool {
105	return (z[3] == x[3]) && (z[2] == x[2]) && (z[1] == x[1]) && (z[0] == x[0])
106}
107
108// IsZero returns z == 0
109func (z *Element) IsZero() bool {
110	return (z[3] | z[2] | z[1] | z[0]) == 0
111}
112
113// field modulus stored as big.Int
114var _elementModulusBigInt big.Int
115var onceelementModulus sync.Once
116
117func elementModulusBigInt() *big.Int {
118	onceelementModulus.Do(func() {
119		_elementModulusBigInt.SetString("21888242871839275222246405745257275088548364400416034343698204186575808495617", 10)
120	})
121	return &_elementModulusBigInt
122}
123
124// Inverse z = x^-1 mod q
125// Algorithm 16 in "Efficient Software-Implementation of Finite Fields with Applications to Cryptography"
126// if x == 0, sets and returns z = x
127func (z *Element) Inverse(x *Element) *Element {
128	if x.IsZero() {
129		return z.Set(x)
130	}
131
132	// initialize u = q
133	var u = Element{
134		4891460686036598785,
135		2896914383306846353,
136		13281191951274694749,
137		3486998266802970665,
138	}
139
140	// initialize s = r^2
141	var s = Element{
142		1997599621687373223,
143		6052339484930628067,
144		10108755138030829701,
145		150537098327114917,
146	}
147
148	// r = 0
149	r := Element{}
150
151	v := *x
152
153	var carry, borrow, t, t2 uint64
154	var bigger, uIsOne, vIsOne bool
155
156	for !uIsOne && !vIsOne {
157		for v[0]&1 == 0 {
158
159			// v = v >> 1
160			t2 = v[3] << 63
161			v[3] >>= 1
162			t = t2
163			t2 = v[2] << 63
164			v[2] = (v[2] >> 1) | t
165			t = t2
166			t2 = v[1] << 63
167			v[1] = (v[1] >> 1) | t
168			t = t2
169			v[0] = (v[0] >> 1) | t
170
171			if s[0]&1 == 1 {
172
173				// s = s + q
174				s[0], carry = bits.Add64(s[0], 4891460686036598785, 0)
175				s[1], carry = bits.Add64(s[1], 2896914383306846353, carry)
176				s[2], carry = bits.Add64(s[2], 13281191951274694749, carry)
177				s[3], _ = bits.Add64(s[3], 3486998266802970665, carry)
178
179			}
180
181			// s = s >> 1
182			t2 = s[3] << 63
183			s[3] >>= 1
184			t = t2
185			t2 = s[2] << 63
186			s[2] = (s[2] >> 1) | t
187			t = t2
188			t2 = s[1] << 63
189			s[1] = (s[1] >> 1) | t
190			t = t2
191			s[0] = (s[0] >> 1) | t
192
193		}
194		for u[0]&1 == 0 {
195
196			// u = u >> 1
197			t2 = u[3] << 63
198			u[3] >>= 1
199			t = t2
200			t2 = u[2] << 63
201			u[2] = (u[2] >> 1) | t
202			t = t2
203			t2 = u[1] << 63
204			u[1] = (u[1] >> 1) | t
205			t = t2
206			u[0] = (u[0] >> 1) | t
207
208			if r[0]&1 == 1 {
209
210				// r = r + q
211				r[0], carry = bits.Add64(r[0], 4891460686036598785, 0)
212				r[1], carry = bits.Add64(r[1], 2896914383306846353, carry)
213				r[2], carry = bits.Add64(r[2], 13281191951274694749, carry)
214				r[3], _ = bits.Add64(r[3], 3486998266802970665, carry)
215
216			}
217
218			// r = r >> 1
219			t2 = r[3] << 63
220			r[3] >>= 1
221			t = t2
222			t2 = r[2] << 63
223			r[2] = (r[2] >> 1) | t
224			t = t2
225			t2 = r[1] << 63
226			r[1] = (r[1] >> 1) | t
227			t = t2
228			r[0] = (r[0] >> 1) | t
229
230		}
231
232		// v >= u
233		bigger = !(v[3] < u[3] || (v[3] == u[3] && (v[2] < u[2] || (v[2] == u[2] && (v[1] < u[1] || (v[1] == u[1] && (v[0] < u[0])))))))
234
235		if bigger {
236
237			// v = v - u
238			v[0], borrow = bits.Sub64(v[0], u[0], 0)
239			v[1], borrow = bits.Sub64(v[1], u[1], borrow)
240			v[2], borrow = bits.Sub64(v[2], u[2], borrow)
241			v[3], _ = bits.Sub64(v[3], u[3], borrow)
242
243			// r >= s
244			bigger = !(r[3] < s[3] || (r[3] == s[3] && (r[2] < s[2] || (r[2] == s[2] && (r[1] < s[1] || (r[1] == s[1] && (r[0] < s[0])))))))
245
246			if bigger {
247
248				// s = s + q
249				s[0], carry = bits.Add64(s[0], 4891460686036598785, 0)
250				s[1], carry = bits.Add64(s[1], 2896914383306846353, carry)
251				s[2], carry = bits.Add64(s[2], 13281191951274694749, carry)
252				s[3], _ = bits.Add64(s[3], 3486998266802970665, carry)
253
254			}
255
256			// s = s - r
257			s[0], borrow = bits.Sub64(s[0], r[0], 0)
258			s[1], borrow = bits.Sub64(s[1], r[1], borrow)
259			s[2], borrow = bits.Sub64(s[2], r[2], borrow)
260			s[3], _ = bits.Sub64(s[3], r[3], borrow)
261
262		} else {
263
264			// u = u - v
265			u[0], borrow = bits.Sub64(u[0], v[0], 0)
266			u[1], borrow = bits.Sub64(u[1], v[1], borrow)
267			u[2], borrow = bits.Sub64(u[2], v[2], borrow)
268			u[3], _ = bits.Sub64(u[3], v[3], borrow)
269
270			// s >= r
271			bigger = !(s[3] < r[3] || (s[3] == r[3] && (s[2] < r[2] || (s[2] == r[2] && (s[1] < r[1] || (s[1] == r[1] && (s[0] < r[0])))))))
272
273			if bigger {
274
275				// r = r + q
276				r[0], carry = bits.Add64(r[0], 4891460686036598785, 0)
277				r[1], carry = bits.Add64(r[1], 2896914383306846353, carry)
278				r[2], carry = bits.Add64(r[2], 13281191951274694749, carry)
279				r[3], _ = bits.Add64(r[3], 3486998266802970665, carry)
280
281			}
282
283			// r = r - s
284			r[0], borrow = bits.Sub64(r[0], s[0], 0)
285			r[1], borrow = bits.Sub64(r[1], s[1], borrow)
286			r[2], borrow = bits.Sub64(r[2], s[2], borrow)
287			r[3], _ = bits.Sub64(r[3], s[3], borrow)
288
289		}
290		uIsOne = (u[0] == 1) && (u[3]|u[2]|u[1]) == 0
291		vIsOne = (v[0] == 1) && (v[3]|v[2]|v[1]) == 0
292	}
293
294	if uIsOne {
295		z.Set(&r)
296	} else {
297		z.Set(&s)
298	}
299
300	return z
301}
302
303// SetRandom sets z to a random element < q
304func (z *Element) SetRandom() *Element {
305	bytes := make([]byte, 32)
306	io.ReadFull(rand.Reader, bytes)
307	z[0] = binary.BigEndian.Uint64(bytes[0:8])
308	z[1] = binary.BigEndian.Uint64(bytes[8:16])
309	z[2] = binary.BigEndian.Uint64(bytes[16:24])
310	z[3] = binary.BigEndian.Uint64(bytes[24:32])
311	z[3] %= 3486998266802970665
312
313	// if z > q --> z -= q
314	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
315		var b uint64
316		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
317		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
318		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
319		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
320	}
321
322	return z
323}
324
325// Add z = x + y mod q
326func (z *Element) Add(x, y *Element) *Element {
327	var carry uint64
328
329	z[0], carry = bits.Add64(x[0], y[0], 0)
330	z[1], carry = bits.Add64(x[1], y[1], carry)
331	z[2], carry = bits.Add64(x[2], y[2], carry)
332	z[3], _ = bits.Add64(x[3], y[3], carry)
333
334	// if z > q --> z -= q
335	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
336		var b uint64
337		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
338		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
339		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
340		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
341	}
342	return z
343}
344
345// AddAssign z = z + x mod q
346func (z *Element) AddAssign(x *Element) *Element {
347	var carry uint64
348
349	z[0], carry = bits.Add64(z[0], x[0], 0)
350	z[1], carry = bits.Add64(z[1], x[1], carry)
351	z[2], carry = bits.Add64(z[2], x[2], carry)
352	z[3], _ = bits.Add64(z[3], x[3], carry)
353
354	// if z > q --> z -= q
355	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
356		var b uint64
357		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
358		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
359		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
360		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
361	}
362	return z
363}
364
365// Double z = x + x mod q, aka Lsh 1
366func (z *Element) Double(x *Element) *Element {
367	var carry uint64
368
369	z[0], carry = bits.Add64(x[0], x[0], 0)
370	z[1], carry = bits.Add64(x[1], x[1], carry)
371	z[2], carry = bits.Add64(x[2], x[2], carry)
372	z[3], _ = bits.Add64(x[3], x[3], carry)
373
374	// if z > q --> z -= q
375	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
376		var b uint64
377		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
378		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
379		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
380		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
381	}
382	return z
383}
384
385// Sub  z = x - y mod q
386func (z *Element) Sub(x, y *Element) *Element {
387	var b uint64
388	z[0], b = bits.Sub64(x[0], y[0], 0)
389	z[1], b = bits.Sub64(x[1], y[1], b)
390	z[2], b = bits.Sub64(x[2], y[2], b)
391	z[3], b = bits.Sub64(x[3], y[3], b)
392	if b != 0 {
393		var c uint64
394		z[0], c = bits.Add64(z[0], 4891460686036598785, 0)
395		z[1], c = bits.Add64(z[1], 2896914383306846353, c)
396		z[2], c = bits.Add64(z[2], 13281191951274694749, c)
397		z[3], _ = bits.Add64(z[3], 3486998266802970665, c)
398	}
399	return z
400}
401
402// SubAssign  z = z - x mod q
403func (z *Element) SubAssign(x *Element) *Element {
404	var b uint64
405	z[0], b = bits.Sub64(z[0], x[0], 0)
406	z[1], b = bits.Sub64(z[1], x[1], b)
407	z[2], b = bits.Sub64(z[2], x[2], b)
408	z[3], b = bits.Sub64(z[3], x[3], b)
409	if b != 0 {
410		var c uint64
411		z[0], c = bits.Add64(z[0], 4891460686036598785, 0)
412		z[1], c = bits.Add64(z[1], 2896914383306846353, c)
413		z[2], c = bits.Add64(z[2], 13281191951274694749, c)
414		z[3], _ = bits.Add64(z[3], 3486998266802970665, c)
415	}
416	return z
417}
418
419// Exp z = x^e mod q
420func (z *Element) Exp(x Element, e uint64) *Element {
421	if e == 0 {
422		return z.SetOne()
423	}
424
425	z.Set(&x)
426
427	l := bits.Len64(e) - 2
428	for i := l; i >= 0; i-- {
429		z.Square(z)
430		if e&(1<<uint(i)) != 0 {
431			z.MulAssign(&x)
432		}
433	}
434	return z
435}
436
437// FromMont converts z in place (i.e. mutates) from Montgomery to regular representation
438// sets and returns z = z * 1
439func (z *Element) FromMont() *Element {
440
441	// the following lines implement z = z * 1
442	// with a modified CIOS montgomery multiplication
443	{
444		// m = z[0]n'[0] mod W
445		m := z[0] * 14042775128853446655
446		C := madd0(m, 4891460686036598785, z[0])
447		C, z[0] = madd2(m, 2896914383306846353, z[1], C)
448		C, z[1] = madd2(m, 13281191951274694749, z[2], C)
449		C, z[2] = madd2(m, 3486998266802970665, z[3], C)
450		z[3] = C
451	}
452	{
453		// m = z[0]n'[0] mod W
454		m := z[0] * 14042775128853446655
455		C := madd0(m, 4891460686036598785, z[0])
456		C, z[0] = madd2(m, 2896914383306846353, z[1], C)
457		C, z[1] = madd2(m, 13281191951274694749, z[2], C)
458		C, z[2] = madd2(m, 3486998266802970665, z[3], C)
459		z[3] = C
460	}
461	{
462		// m = z[0]n'[0] mod W
463		m := z[0] * 14042775128853446655
464		C := madd0(m, 4891460686036598785, z[0])
465		C, z[0] = madd2(m, 2896914383306846353, z[1], C)
466		C, z[1] = madd2(m, 13281191951274694749, z[2], C)
467		C, z[2] = madd2(m, 3486998266802970665, z[3], C)
468		z[3] = C
469	}
470	{
471		// m = z[0]n'[0] mod W
472		m := z[0] * 14042775128853446655
473		C := madd0(m, 4891460686036598785, z[0])
474		C, z[0] = madd2(m, 2896914383306846353, z[1], C)
475		C, z[1] = madd2(m, 13281191951274694749, z[2], C)
476		C, z[2] = madd2(m, 3486998266802970665, z[3], C)
477		z[3] = C
478	}
479
480	// if z > q --> z -= q
481	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
482		var b uint64
483		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
484		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
485		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
486		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
487	}
488	return z
489}
490
491// ToMont converts z to Montgomery form
492// sets and returns z = z * r^2
493func (z *Element) ToMont() *Element {
494	var rSquare = Element{
495		1997599621687373223,
496		6052339484930628067,
497		10108755138030829701,
498		150537098327114917,
499	}
500	return z.MulAssign(&rSquare)
501}
502
503// ToRegular returns z in regular form (doesn't mutate z)
504func (z Element) ToRegular() Element {
505	return *z.FromMont()
506}
507
508// String returns the string form of an Element in Montgomery form
509func (z *Element) String() string {
510	var _z big.Int
511	return z.ToBigIntRegular(&_z).String()
512}
513
514// ToBigInt returns z as a big.Int in Montgomery form
515func (z *Element) ToBigInt(res *big.Int) *big.Int {
516	bits := (*[4]big.Word)(unsafe.Pointer(z))
517	return res.SetBits(bits[:])
518}
519
520// ToBigIntRegular returns z as a big.Int in regular form
521func (z Element) ToBigIntRegular(res *big.Int) *big.Int {
522	z.FromMont()
523	bits := (*[4]big.Word)(unsafe.Pointer(&z))
524	return res.SetBits(bits[:])
525}
526
527// SetBigInt sets z to v (regular form) and returns z in Montgomery form
528func (z *Element) SetBigInt(v *big.Int) *Element {
529	z.SetZero()
530
531	zero := big.NewInt(0)
532	q := elementModulusBigInt()
533
534	// copy input
535	vv := new(big.Int).Set(v)
536
537	// while v < 0, v+=q
538	for vv.Cmp(zero) == -1 {
539		vv.Add(vv, q)
540	}
541	// while v > q, v-=q
542	for vv.Cmp(q) == 1 {
543		vv.Sub(vv, q)
544	}
545	// if v == q, return 0
546	if vv.Cmp(q) == 0 {
547		return z
548	}
549	// v should
550	vBits := vv.Bits()
551	for i := 0; i < len(vBits); i++ {
552		z[i] = uint64(vBits[i])
553	}
554	return z.ToMont()
555}
556
557// SetString creates a big.Int with s (in base 10) and calls SetBigInt on z
558func (z *Element) SetString(s string) *Element {
559	x, ok := new(big.Int).SetString(s, 10)
560	if !ok {
561		panic("Element.SetString failed -> can't parse number in base10 into a big.Int")
562	}
563	return z.SetBigInt(x)
564}
565
566// Mul z = x * y mod q
567func (z *Element) Mul(x, y *Element) *Element {
568
569	var t [4]uint64
570	var c [3]uint64
571	{
572		// round 0
573		v := x[0]
574		c[1], c[0] = bits.Mul64(v, y[0])
575		m := c[0] * 14042775128853446655
576		c[2] = madd0(m, 4891460686036598785, c[0])
577		c[1], c[0] = madd1(v, y[1], c[1])
578		c[2], t[0] = madd2(m, 2896914383306846353, c[2], c[0])
579		c[1], c[0] = madd1(v, y[2], c[1])
580		c[2], t[1] = madd2(m, 13281191951274694749, c[2], c[0])
581		c[1], c[0] = madd1(v, y[3], c[1])
582		t[3], t[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
583	}
584	{
585		// round 1
586		v := x[1]
587		c[1], c[0] = madd1(v, y[0], t[0])
588		m := c[0] * 14042775128853446655
589		c[2] = madd0(m, 4891460686036598785, c[0])
590		c[1], c[0] = madd2(v, y[1], c[1], t[1])
591		c[2], t[0] = madd2(m, 2896914383306846353, c[2], c[0])
592		c[1], c[0] = madd2(v, y[2], c[1], t[2])
593		c[2], t[1] = madd2(m, 13281191951274694749, c[2], c[0])
594		c[1], c[0] = madd2(v, y[3], c[1], t[3])
595		t[3], t[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
596	}
597	{
598		// round 2
599		v := x[2]
600		c[1], c[0] = madd1(v, y[0], t[0])
601		m := c[0] * 14042775128853446655
602		c[2] = madd0(m, 4891460686036598785, c[0])
603		c[1], c[0] = madd2(v, y[1], c[1], t[1])
604		c[2], t[0] = madd2(m, 2896914383306846353, c[2], c[0])
605		c[1], c[0] = madd2(v, y[2], c[1], t[2])
606		c[2], t[1] = madd2(m, 13281191951274694749, c[2], c[0])
607		c[1], c[0] = madd2(v, y[3], c[1], t[3])
608		t[3], t[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
609	}
610	{
611		// round 3
612		v := x[3]
613		c[1], c[0] = madd1(v, y[0], t[0])
614		m := c[0] * 14042775128853446655
615		c[2] = madd0(m, 4891460686036598785, c[0])
616		c[1], c[0] = madd2(v, y[1], c[1], t[1])
617		c[2], z[0] = madd2(m, 2896914383306846353, c[2], c[0])
618		c[1], c[0] = madd2(v, y[2], c[1], t[2])
619		c[2], z[1] = madd2(m, 13281191951274694749, c[2], c[0])
620		c[1], c[0] = madd2(v, y[3], c[1], t[3])
621		z[3], z[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
622	}
623
624	// if z > q --> z -= q
625	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
626		var b uint64
627		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
628		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
629		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
630		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
631	}
632	return z
633}
634
635// MulAssign z = z * x mod q
636func (z *Element) MulAssign(x *Element) *Element {
637
638	var t [4]uint64
639	var c [3]uint64
640	{
641		// round 0
642		v := z[0]
643		c[1], c[0] = bits.Mul64(v, x[0])
644		m := c[0] * 14042775128853446655
645		c[2] = madd0(m, 4891460686036598785, c[0])
646		c[1], c[0] = madd1(v, x[1], c[1])
647		c[2], t[0] = madd2(m, 2896914383306846353, c[2], c[0])
648		c[1], c[0] = madd1(v, x[2], c[1])
649		c[2], t[1] = madd2(m, 13281191951274694749, c[2], c[0])
650		c[1], c[0] = madd1(v, x[3], c[1])
651		t[3], t[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
652	}
653	{
654		// round 1
655		v := z[1]
656		c[1], c[0] = madd1(v, x[0], t[0])
657		m := c[0] * 14042775128853446655
658		c[2] = madd0(m, 4891460686036598785, c[0])
659		c[1], c[0] = madd2(v, x[1], c[1], t[1])
660		c[2], t[0] = madd2(m, 2896914383306846353, c[2], c[0])
661		c[1], c[0] = madd2(v, x[2], c[1], t[2])
662		c[2], t[1] = madd2(m, 13281191951274694749, c[2], c[0])
663		c[1], c[0] = madd2(v, x[3], c[1], t[3])
664		t[3], t[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
665	}
666	{
667		// round 2
668		v := z[2]
669		c[1], c[0] = madd1(v, x[0], t[0])
670		m := c[0] * 14042775128853446655
671		c[2] = madd0(m, 4891460686036598785, c[0])
672		c[1], c[0] = madd2(v, x[1], c[1], t[1])
673		c[2], t[0] = madd2(m, 2896914383306846353, c[2], c[0])
674		c[1], c[0] = madd2(v, x[2], c[1], t[2])
675		c[2], t[1] = madd2(m, 13281191951274694749, c[2], c[0])
676		c[1], c[0] = madd2(v, x[3], c[1], t[3])
677		t[3], t[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
678	}
679	{
680		// round 3
681		v := z[3]
682		c[1], c[0] = madd1(v, x[0], t[0])
683		m := c[0] * 14042775128853446655
684		c[2] = madd0(m, 4891460686036598785, c[0])
685		c[1], c[0] = madd2(v, x[1], c[1], t[1])
686		c[2], z[0] = madd2(m, 2896914383306846353, c[2], c[0])
687		c[1], c[0] = madd2(v, x[2], c[1], t[2])
688		c[2], z[1] = madd2(m, 13281191951274694749, c[2], c[0])
689		c[1], c[0] = madd2(v, x[3], c[1], t[3])
690		z[3], z[2] = madd3(m, 3486998266802970665, c[0], c[2], c[1])
691	}
692
693	// if z > q --> z -= q
694	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
695		var b uint64
696		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
697		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
698		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
699		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
700	}
701	return z
702}
703
704// Square z = x * x mod q
705func (z *Element) Square(x *Element) *Element {
706
707	var p [4]uint64
708
709	var u, v uint64
710	{
711		// round 0
712		u, p[0] = bits.Mul64(x[0], x[0])
713		m := p[0] * 14042775128853446655
714		C := madd0(m, 4891460686036598785, p[0])
715		var t uint64
716		t, u, v = madd1sb(x[0], x[1], u)
717		C, p[0] = madd2(m, 2896914383306846353, v, C)
718		t, u, v = madd1s(x[0], x[2], t, u)
719		C, p[1] = madd2(m, 13281191951274694749, v, C)
720		_, u, v = madd1s(x[0], x[3], t, u)
721		p[3], p[2] = madd3(m, 3486998266802970665, v, C, u)
722	}
723	{
724		// round 1
725		m := p[0] * 14042775128853446655
726		C := madd0(m, 4891460686036598785, p[0])
727		u, v = madd1(x[1], x[1], p[1])
728		C, p[0] = madd2(m, 2896914383306846353, v, C)
729		var t uint64
730		t, u, v = madd2sb(x[1], x[2], p[2], u)
731		C, p[1] = madd2(m, 13281191951274694749, v, C)
732		_, u, v = madd2s(x[1], x[3], p[3], t, u)
733		p[3], p[2] = madd3(m, 3486998266802970665, v, C, u)
734	}
735	{
736		// round 2
737		m := p[0] * 14042775128853446655
738		C := madd0(m, 4891460686036598785, p[0])
739		C, p[0] = madd2(m, 2896914383306846353, p[1], C)
740		u, v = madd1(x[2], x[2], p[2])
741		C, p[1] = madd2(m, 13281191951274694749, v, C)
742		_, u, v = madd2sb(x[2], x[3], p[3], u)
743		p[3], p[2] = madd3(m, 3486998266802970665, v, C, u)
744	}
745	{
746		// round 3
747		m := p[0] * 14042775128853446655
748		C := madd0(m, 4891460686036598785, p[0])
749		C, z[0] = madd2(m, 2896914383306846353, p[1], C)
750		C, z[1] = madd2(m, 13281191951274694749, p[2], C)
751		u, v = madd1(x[3], x[3], p[3])
752		z[3], z[2] = madd3(m, 3486998266802970665, v, C, u)
753	}
754
755	// if z > q --> z -= q
756	if !(z[3] < 3486998266802970665 || (z[3] == 3486998266802970665 && (z[2] < 13281191951274694749 || (z[2] == 13281191951274694749 && (z[1] < 2896914383306846353 || (z[1] == 2896914383306846353 && (z[0] < 4891460686036598785))))))) {
757		var b uint64
758		z[0], b = bits.Sub64(z[0], 4891460686036598785, 0)
759		z[1], b = bits.Sub64(z[1], 2896914383306846353, b)
760		z[2], b = bits.Sub64(z[2], 13281191951274694749, b)
761		z[3], _ = bits.Sub64(z[3], 3486998266802970665, b)
762	}
763	return z
764}