1// This package implements efficient Golidlocks arithmetic operations within Gnark. We do not use
  2// the emulated field arithmetic API, because it is too slow for our purposes. Instead, we use
  3// an efficient reduction method that leverages the fact that the modulus is a simple
  4// linear combination of powers of two.
  5package goldilocks
  6
  7// In general, methods whose name do not contain `NoReduce` can be used without any extra mental
  8// overhead. These methods act exactly as you would expect a normal field would operate.
  9//
 10// However, if you want to aggressively optimize the number of constraints in your circuit, it can
 11// be very beneficial to use the no reduction methods and keep track of the maximum number of bits
 12// your computation uses.
 13
 14// This implementation is based on the following plonky2 implementation of Goldilocks
 15// Available here: https://github.com/0xPolygonZero/plonky2/blob/main/field/src/goldilocks_field.rs#L70
 16
 17import (
 18	"fmt"
 19	"math"
 20	"math/big"
 21	"os"
 22
 23	"github.com/consensys/gnark-crypto/field/goldilocks"
 24	"github.com/consensys/gnark/constraint/solver"
 25	"github.com/consensys/gnark/frontend"
 26	"github.com/consensys/gnark/std/math/emulated"
 27	"github.com/consensys/gnark/std/rangecheck"
 28)
 29
 30// The multiplicative group generator of the field.
 31var MULTIPLICATIVE_GROUP_GENERATOR goldilocks.Element = goldilocks.NewElement(7)
 32
 33// The two adicity of the field.
 34var TWO_ADICITY uint64 = 32
 35
 36// The power of two generator of the field.
 37var POWER_OF_TWO_GENERATOR goldilocks.Element = goldilocks.NewElement(1753635133440165772)
 38
 39// The modulus of the field.
 40var MODULUS *big.Int = emulated.Goldilocks{}.Modulus()
 41
 42// The threshold maximum number of bits at which we must reduce the element.
 43var REDUCE_NB_BITS_THRESHOLD uint8 = 254 - 64
 44
 45// The number of bits to use for range checks on inner products of field elements.
 46var RANGE_CHECK_NB_BITS int = 140
 47
 48// Registers the hint functions with the solver.
 49func init() {
 50	solver.RegisterHint(MulAddHint)
 51	solver.RegisterHint(ReduceHint)
 52	solver.RegisterHint(InverseHint)
 53	solver.RegisterHint(SplitLimbsHint)
 54}
 55
 56// A type alias used to represent Goldilocks field elements.
 57type Variable struct {
 58	Limb frontend.Variable
 59}
 60
 61// Creates a new Goldilocks field element from an existing variable. Assumes that the element is
 62// already reduced.
 63func NewVariable(x frontend.Variable) Variable {
 64	return Variable{Limb: x}
 65}
 66
 67// The zero element in the Golidlocks field.
 68func Zero() Variable {
 69	return NewVariable(0)
 70}
 71
 72// The one element in the Goldilocks field.
 73func One() Variable {
 74	return NewVariable(1)
 75}
 76
 77// The negative one element in the Goldilocks field.
 78func NegOne() Variable {
 79	return NewVariable(MODULUS.Uint64() - 1)
 80}
 81
 82// The chip used for Goldilocks field operations.
 83type Chip struct {
 84	api          frontend.API
 85	rangeChecker frontend.Rangechecker
 86}
 87
 88// Creates a new Goldilocks Chip.
 89func New(api frontend.API) *Chip {
 90	use_bit_decomp := os.Getenv("USE_BIT_DECOMPOSITION_RANGE_CHECK")
 91
 92	var rangeChecker frontend.Rangechecker
 93
 94	// If USE_BIT_DECOMPOSITION_RANGE_CHECK is not set, then use the std.rangecheck New function
 95	if use_bit_decomp == "" {
 96		rangeChecker = rangecheck.New(api)
 97	} else {
 98		rangeChecker = bitDecompChecker{api: api}
 99	}
100	return &Chip{api: api, rangeChecker: rangeChecker}
101}
102
103// Adds two field elements such that x + y = z within the Golidlocks field.
104func (p *Chip) Add(a Variable, b Variable) Variable {
105	return p.MulAdd(a, NewVariable(1), b)
106}
107
108// Adds two field elements such that x + y = z within the Golidlocks field without reducing.
109func (p *Chip) AddNoReduce(a Variable, b Variable) Variable {
110	return NewVariable(p.api.Add(a.Limb, b.Limb))
111}
112
113// Subtracts two field elements such that x + y = z within the Golidlocks field.
114func (p *Chip) Sub(a Variable, b Variable) Variable {
115	return p.MulAdd(b, NewVariable(MODULUS.Uint64()-1), a)
116}
117
118// Subtracts two field elements such that x + y = z within the Golidlocks field without reducing.
119func (p *Chip) SubNoReduce(a Variable, b Variable) Variable {
120	return NewVariable(p.api.Add(a.Limb, p.api.Mul(b.Limb, MODULUS.Uint64()-1)))
121}
122
123// Multiplies two field elements such that x * y = z within the Golidlocks field.
124func (p *Chip) Mul(a Variable, b Variable) Variable {
125	return p.MulAdd(a, b, Zero())
126}
127
128// Multiplies two field elements such that x * y = z within the Golidlocks field without reducing.
129func (p *Chip) MulNoReduce(a Variable, b Variable) Variable {
130	return NewVariable(p.api.Mul(a.Limb, b.Limb))
131}
132
133// Multiplies two field elements and adds a field element such that x * y + z = c within the
134// Golidlocks field.
135func (p *Chip) MulAdd(a Variable, b Variable, c Variable) Variable {
136	result, err := p.api.Compiler().NewHint(MulAddHint, 2, a.Limb, b.Limb, c.Limb)
137	if err != nil {
138		panic(err)
139	}
140
141	quotient := NewVariable(result[0])
142	remainder := NewVariable(result[1])
143
144	lhs := p.api.MulAcc(c.Limb, a.Limb, b.Limb)
145	rhs := p.api.MulAcc(remainder.Limb, MODULUS, quotient.Limb)
146	p.api.AssertIsEqual(lhs, rhs)
147
148	p.RangeCheck(quotient)
149	p.RangeCheck(remainder)
150	return remainder
151}
152
153// Multiplies two field elements and adds a field element such that x * y + z = c within the
154// Golidlocks field without reducing.
155func (p *Chip) MulAddNoReduce(a Variable, b Variable, c Variable) Variable {
156	return NewVariable(p.api.MulAcc(c.Limb, a.Limb, b.Limb))
157}
158
159// The hint used to compute MulAdd.
160func MulAddHint(_ *big.Int, inputs []*big.Int, results []*big.Int) error {
161	if len(inputs) != 3 {
162		panic("MulAddHint expects 3 input operands")
163	}
164
165	for _, operand := range inputs {
166		if operand.Cmp(MODULUS) >= 0 {
167			panic(fmt.Sprintf("%s is not in the field", operand.String()))
168		}
169	}
170
171	product := new(big.Int).Mul(inputs[0], inputs[1])
172	sum := new(big.Int).Add(product, inputs[2])
173	quotient := new(big.Int).Div(sum, MODULUS)
174	remainder := new(big.Int).Rem(sum, MODULUS)
175
176	results[0] = quotient
177	results[1] = remainder
178
179	return nil
180}
181
182// Reduces a field element x such that x % MODULUS = y.
183func (p *Chip) Reduce(x Variable) Variable {
184	// Witness a `quotient` and `remainder` such that:
185	//
186	// 		MODULUS * quotient + remainder = x
187	//
188	// Must check that offset \in [0, MODULUS) and carry \in [0, 2^RANGE_CHECK_NB_BITS) to ensure
189	// that this computation does not overflow. We use 2^RANGE_CHECK_NB_BITS to reduce the cost of the range check
190	//
191	// In other words, we assume that we at most compute a a dot product with dimension at most RANGE_CHECK_NB_BITS - 128.
192
193	result, err := p.api.Compiler().NewHint(ReduceHint, 2, x.Limb)
194	if err != nil {
195		panic(err)
196	}
197
198	quotient := result[0]
199	p.rangeChecker.Check(quotient, RANGE_CHECK_NB_BITS)
200
201	remainder := NewVariable(result[1])
202	p.RangeCheck(remainder)
203	return remainder
204}
205
206// Reduces a field element x such that x % MODULUS = y.
207func (p *Chip) ReduceWithMaxBits(x Variable, maxNbBits uint64) Variable {
208	// Witness a `quotient` and `remainder` such that:
209	//
210	// 		MODULUS * quotient + remainder = x
211	//
212	// Must check that remainder \in [0, MODULUS) and quotient \in [0, 2^maxNbBits) to ensure that this
213	// computation does not overflow.
214
215	result, err := p.api.Compiler().NewHint(ReduceHint, 2, x.Limb)
216	if err != nil {
217		panic(err)
218	}
219
220	quotient := result[0]
221	p.rangeChecker.Check(quotient, int(maxNbBits))
222
223	remainder := NewVariable(result[1])
224	p.RangeCheck(remainder)
225	return remainder
226}
227
228// The hint used to compute Reduce.
229func ReduceHint(_ *big.Int, inputs []*big.Int, results []*big.Int) error {
230	if len(inputs) != 1 {
231		panic("ReduceHint expects 1 input operand")
232	}
233	input := inputs[0]
234	quotient := new(big.Int).Div(input, MODULUS)
235	remainder := new(big.Int).Rem(input, MODULUS)
236	results[0] = quotient
237	results[1] = remainder
238	return nil
239}
240
241// Computes the inverse of a field element x such that x * x^-1 = 1.
242func (p *Chip) Inverse(x Variable) Variable {
243	result, err := p.api.Compiler().NewHint(InverseHint, 1, x.Limb)
244	if err != nil {
245		panic(err)
246	}
247
248	inverse := NewVariable(result[0])
249	product := p.Mul(inverse, x)
250	p.api.AssertIsEqual(product.Limb, frontend.Variable(1))
251	return inverse
252}
253
254// The hint used to compute Inverse.
255func InverseHint(_ *big.Int, inputs []*big.Int, results []*big.Int) error {
256	if len(inputs) != 1 {
257		panic("InverseHint expects 1 input operand")
258	}
259
260	input := inputs[0]
261	if input.Cmp(MODULUS) == 0 || input.Cmp(MODULUS) == 1 {
262		panic("Input is not in the field")
263	}
264
265	inputGl := goldilocks.NewElement(input.Uint64())
266	resultGl := goldilocks.NewElement(0)
267	resultGl.Inverse(&inputGl)
268
269	result := big.NewInt(0)
270	results[0] = resultGl.BigInt(result)
271
272	return nil
273}
274
275// Computes a field element raised to some power.
276func (p *Chip) Exp(x Variable, k *big.Int) Variable {
277	if k.IsUint64() && k.Uint64() == 0 {
278		return One()
279	}
280
281	e := k
282	if k.Sign() == -1 {
283		panic("Unsupported negative exponent. Need to implement inversion.")
284	}
285
286	z := x
287	for i := e.BitLen() - 2; i >= 0; i-- {
288		z = p.Mul(z, z)
289		if e.Bit(i) == 1 {
290			z = p.Mul(z, x)
291		}
292	}
293
294	return z
295}
296
297// The hint used to split a GoldilocksVariable into 2 32 bit limbs.
298func SplitLimbsHint(_ *big.Int, inputs []*big.Int, results []*big.Int) error {
299	if len(inputs) != 1 {
300		panic("SplitLimbsHint expects 1 input operand")
301	}
302
303	// The Goldilocks field element
304	input := inputs[0]
305
306	if input.Cmp(MODULUS) == 0 || input.Cmp(MODULUS) == 1 {
307		return fmt.Errorf("input is not in the field")
308	}
309
310	two_32 := big.NewInt(int64(math.Pow(2, 32)))
311
312	// The most significant bits
313	results[0] = new(big.Int).Quo(input, two_32)
314	// The least significant bits
315	results[1] = new(big.Int).Rem(input, two_32)
316
317	return nil
318}
319
320// Range checks a field element x to be less than the Golidlocks modulus 2 ^ 64 - 2 ^ 32 + 1.
321func (p *Chip) RangeCheck(x Variable) {
322	// The Goldilocks' modulus is 2^64 - 2^32 + 1, which is:
323	//
324	// 		1111111111111111111111111111111100000000000000000000000000000001
325	//
326	// in big endian binary. This function will first verify that x is at most 64 bits wide. Then it
327	// checks that if the bits[0:31] (in big-endian) are all 1, then bits[32:64] are all zero.
328
329	result, err := p.api.Compiler().NewHint(SplitLimbsHint, 2, x.Limb)
330	if err != nil {
331		panic(err)
332	}
333
334	// We check that this is a valid decomposition of the Goldilock's element and range-check each limb.
335	mostSigLimb := result[0]
336	leastSigLimb := result[1]
337	p.api.AssertIsEqual(
338		p.api.Add(
339			p.api.Mul(mostSigLimb, uint64(math.Pow(2, 32))),
340			leastSigLimb,
341		),
342		x.Limb,
343	)
344	p.rangeChecker.Check(mostSigLimb, 32)
345	p.rangeChecker.Check(leastSigLimb, 32)
346
347	// If the most significant bits are all 1, then we need to check that the least significant bits are all zero
348	// in order for element to be less than the Goldilock's modulus.
349	// Otherwise, we don't need to do any checks, since we already know that the element is less than the Goldilocks modulus.
350	shouldCheck := p.api.IsZero(p.api.Sub(mostSigLimb, uint64(math.Pow(2, 32))-1))
351	p.api.AssertIsEqual(
352		p.api.Select(
353			shouldCheck,
354			leastSigLimb,
355			frontend.Variable(0),
356		),
357		frontend.Variable(0),
358	)
359}
360
361func (p *Chip) AssertIsEqual(x, y Variable) {
362	p.api.AssertIsEqual(x.Limb, y.Limb)
363}
364
365// Computes the n'th primitive root of unity for the Goldilocks field.
366func PrimitiveRootOfUnity(nLog uint64) goldilocks.Element {
367	if nLog > TWO_ADICITY {
368		panic("nLog is greater than TWO_ADICITY")
369	}
370	res := goldilocks.NewElement(POWER_OF_TWO_GENERATOR.Uint64())
371	for i := 0; i < int(TWO_ADICITY-nLog); i++ {
372		res.Square(&res)
373	}
374	return res
375}
376
377func TwoAdicSubgroup(nLog uint64) []goldilocks.Element {
378	if nLog > TWO_ADICITY {
379		panic("nLog is greater than GOLDILOCKS_TWO_ADICITY")
380	}
381
382	var res []goldilocks.Element
383	rootOfUnity := PrimitiveRootOfUnity(nLog)
384	res = append(res, goldilocks.NewElement(1))
385
386	for i := 0; i < (1 << nLog); i++ {
387		lastElement := res[len(res)-1]
388		res = append(res, *lastElement.Mul(&lastElement, &rootOfUnity))
389	}
390
391	return res
392}