1package fri
  2
  3import (
  4	"fmt"
  5	"math"
  6	"math/big"
  7	"math/bits"
  8
  9	"github.com/consensys/gnark-crypto/field/goldilocks"
 10	"github.com/consensys/gnark/frontend"
 11	gl "github.com/succinctlabs/gnark-plonky2-verifier/goldilocks"
 12	"github.com/succinctlabs/gnark-plonky2-verifier/poseidon"
 13	"github.com/succinctlabs/gnark-plonky2-verifier/types"
 14	"github.com/succinctlabs/gnark-plonky2-verifier/variables"
 15)
 16
 17type Chip struct {
 18	api               frontend.API `gnark:"-"`
 19	gl                gl.Chip      `gnark:"-"`
 20	poseidonBN254Chip *poseidon.BN254Chip
 21	friParams         *types.FriParams `gnark:"-"`
 22}
 23
 24func NewChip(
 25	api frontend.API,
 26	friParams *types.FriParams,
 27) *Chip {
 28	poseidonBN254Chip := poseidon.NewBN254Chip(api)
 29	return &Chip{
 30		api:               api,
 31		poseidonBN254Chip: poseidonBN254Chip,
 32		friParams:         friParams,
 33		gl:                *gl.New(api),
 34	}
 35}
 36
 37func (f *Chip) assertLeadingZeros(powWitness gl.Variable, friConfig types.FriConfig) {
 38	// Asserts that powWitness'es big-endian bit representation has at least `leading_zeros` leading zeros.
 39	// Note that this is assuming that the Goldilocks field is being used.  Specfically that the
 40	// field is 64 bits long
 41	maxPowWitness := uint64(math.Pow(2, float64(64-friConfig.ProofOfWorkBits))) - 1
 42	reducedPowWitness := f.gl.Reduce(powWitness)
 43	f.api.AssertIsLessOrEqual(reducedPowWitness.Limb, frontend.Variable(maxPowWitness))
 44}
 45
 46func (f *Chip) fromOpeningsAndAlpha(
 47	openings *Openings,
 48	alpha gl.QuadraticExtensionVariable,
 49) []gl.QuadraticExtensionVariable {
 50	// One reduced opening for all openings evaluated at point Zeta.
 51	// Another one for all openings evaluated at point Zeta * Omega (which is only PlonkZsNext polynomial)
 52
 53	reducedOpenings := make([]gl.QuadraticExtensionVariable, 0, 2)
 54	for _, batch := range openings.Batches {
 55		reducedOpenings = append(reducedOpenings, f.gl.ReduceWithPowers(batch.Values, alpha))
 56	}
 57
 58	return reducedOpenings
 59}
 60
 61func (f *Chip) verifyMerkleProofToCapWithCapIndex(
 62	leafData []gl.Variable,
 63	leafIndexBits []frontend.Variable,
 64	capIndexBits []frontend.Variable,
 65	merkleCap variables.FriMerkleCap,
 66	proof *variables.FriMerkleProof,
 67) {
 68	currentDigest := f.poseidonBN254Chip.HashOrNoop(leafData)
 69	for i, sibling := range proof.Siblings {
 70		bit := leafIndexBits[i]
 71		// TODO: Don't need to do two hashes by using a trick that the plonky2 verifier circuit does
 72		// https://github.com/mir-protocol/plonky2/blob/973624f12d2d12d74422b3ea051358b9eaacb050/plonky2/src/gates/poseidon.rs#L298
 73		leftHash := f.poseidonBN254Chip.TwoToOne(sibling, currentDigest)
 74		rightHash := f.poseidonBN254Chip.TwoToOne(currentDigest, sibling)
 75		currentDigest = f.api.Select(bit, leftHash, rightHash)
 76	}
 77
 78	// We assume that the cap_height is 4.  Create two levels of the Lookup2 circuit
 79	if len(capIndexBits) != 4 || len(merkleCap) != 16 {
 80		errorMsg, _ := fmt.Printf(
 81			"capIndexBits length should be 4 and the merkleCap length should be 16.  Actual values (capIndexBits: %d, merkleCap: %d)\n",
 82			len(capIndexBits),
 83			len(merkleCap),
 84		)
 85		panic(errorMsg)
 86	}
 87
 88	const NUM_LEAF_LOOKUPS = 4
 89	var leafLookups [NUM_LEAF_LOOKUPS]poseidon.BN254HashOut
 90	// First create the "leaf" lookup2 circuits
 91	// The will use the least significant bits of the capIndexBits array
 92	for i := 0; i < NUM_LEAF_LOOKUPS; i++ {
 93		leafLookups[i] = f.api.Lookup2(
 94			capIndexBits[0], capIndexBits[1],
 95			merkleCap[i*NUM_LEAF_LOOKUPS], merkleCap[i*NUM_LEAF_LOOKUPS+1], merkleCap[i*NUM_LEAF_LOOKUPS+2], merkleCap[i*NUM_LEAF_LOOKUPS+3],
 96		)
 97	}
 98
 99	// Use the most 2 significant bits of the capIndexBits array for the "root" lookup
100	merkleCapEntry := f.api.Lookup2(capIndexBits[2], capIndexBits[3], leafLookups[0], leafLookups[1], leafLookups[2], leafLookups[3])
101	f.api.AssertIsEqual(currentDigest, merkleCapEntry)
102}
103
104func (f *Chip) verifyInitialProof(xIndexBits []frontend.Variable, proof *variables.FriInitialTreeProof, initialMerkleCaps []variables.FriMerkleCap, capIndexBits []frontend.Variable) {
105	if len(proof.EvalsProofs) != len(initialMerkleCaps) {
106		panic("length of eval proofs in fri proof should equal length of initial merkle caps")
107	}
108
109	for i := 0; i < len(initialMerkleCaps); i++ {
110		evals := proof.EvalsProofs[i].Elements
111		merkleProof := proof.EvalsProofs[i].MerkleProof
112		cap := initialMerkleCaps[i]
113		f.verifyMerkleProofToCapWithCapIndex(evals, xIndexBits, capIndexBits, cap, &merkleProof)
114	}
115}
116
117// / We decompose FRI query indices into bits without verifying that the decomposition given by
118// / the prover is the canonical one. In particular, if `x_index < 2^field_bits - p`, then the
119// / prover could supply the binary encoding of either `x_index` or `x_index + p`, since they are
120// / congruent mod `p`. However, this only occurs with probability
121// /     p_ambiguous = (2^field_bits - p) / p
122// / which is small for the field that we use in practice.
123// /
124// / In particular, the soundness error of one FRI query is roughly the codeword rate, which
125// / is much larger than this ambiguous-element probability given any reasonable parameters.
126// / Thus ambiguous elements contribute a negligible amount to soundness error.
127// /
128// / Here we compare the probabilities as a sanity check, to verify the claim above.
129func (f *Chip) assertNoncanonicalIndicesOK() {
130	numAmbiguousElems := uint64(math.MaxUint64) - goldilocks.Modulus().Uint64() + 1
131	queryError := f.friParams.Config.Rate()
132	pAmbiguous := float64(numAmbiguousElems) / float64(goldilocks.Modulus().Uint64())
133
134	// TODO:  Check that pAmbiguous value is the same as the one in plonky2 verifier
135	if pAmbiguous >= queryError*1e-5 {
136		panic("A non-negligible portion of field elements are in the range that permits non-canonical encodings. Need to do more analysis or enforce canonical encodings.")
137	}
138}
139
140func (f *Chip) expFromBitsConstBase(
141	base goldilocks.Element,
142	exponentBits []frontend.Variable,
143) gl.Variable {
144	product := gl.One()
145	for i, bit := range exponentBits {
146		// If the bit is on, we multiply product by base^pow.
147		// We can arithmetize this as:
148		//     product *= 1 + bit (base^pow - 1)
149		//     product = (base^pow - 1) product bit + product
150		pow := int64(1 << i)
151		basePow := goldilocks.NewElement(0)
152		basePow.Exp(base, big.NewInt(pow))
153		basePowVariable := gl.NewVariable(basePow.Uint64() - 1)
154		product = f.gl.Add(
155			f.gl.Mul(
156				f.gl.Mul(
157					basePowVariable,
158					product,
159				),
160				gl.NewVariable(bit),
161			),
162			product,
163		)
164	}
165	return product
166}
167
168func (f *Chip) calculateSubgroupX(
169	xIndexBits []frontend.Variable,
170	nLog uint64,
171) gl.Variable {
172	// Compute x from its index
173	// `subgroup_x` is `subgroup[x_index]`, i.e., the actual field element in the domain.
174	// TODO - Make these as global values
175	g := gl.NewVariable(gl.MULTIPLICATIVE_GROUP_GENERATOR.Uint64())
176	base := gl.PrimitiveRootOfUnity(nLog)
177
178	// Create a reverse list of xIndexBits
179	xIndexBitsRev := make([]frontend.Variable, 0)
180	for i := len(xIndexBits) - 1; i >= 0; i-- {
181		xIndexBitsRev = append(xIndexBitsRev, xIndexBits[i])
182	}
183
184	product := f.expFromBitsConstBase(base, xIndexBitsRev)
185
186	return f.gl.Mul(g, product)
187}
188
189func (f *Chip) friCombineInitial(
190	instance InstanceInfo,
191	proof variables.FriInitialTreeProof,
192	friAlpha gl.QuadraticExtensionVariable,
193	subgroupX_QE gl.QuadraticExtensionVariable,
194	precomputedReducedEval []gl.QuadraticExtensionVariable,
195) gl.QuadraticExtensionVariable {
196	sum := gl.ZeroExtension()
197
198	if len(instance.Batches) != len(precomputedReducedEval) {
199		panic("len(openings) != len(precomputedReducedEval)")
200	}
201
202	for i := 0; i < len(instance.Batches); i++ {
203		batch := instance.Batches[i]
204		reducedOpenings := precomputedReducedEval[i]
205
206		point := batch.Point
207		evals := make([]gl.QuadraticExtensionVariable, 0)
208		for _, polynomial := range batch.Polynomials {
209			evals = append(
210				evals,
211				gl.QuadraticExtensionVariable{
212					proof.EvalsProofs[polynomial.OracleIndex].Elements[polynomial.PolynomialInfo],
213					gl.Zero(),
214				},
215			)
216		}
217
218		reducedEvals := f.gl.ReduceWithPowers(evals, friAlpha)
219		numerator := f.gl.SubExtensionNoReduce(reducedEvals, reducedOpenings)
220		denominator := f.gl.SubExtension(subgroupX_QE, point)
221		sum = f.gl.MulExtension(f.gl.ExpExtension(friAlpha, uint64(len(evals))), sum)
222		sum = f.gl.MulAddExtension(
223			numerator,
224			f.gl.InverseExtension(denominator),
225			sum,
226		)
227	}
228
229	return sum
230}
231
232func (f *Chip) finalPolyEval(finalPoly variables.PolynomialCoeffs, point gl.QuadraticExtensionVariable) gl.QuadraticExtensionVariable {
233	ret := gl.ZeroExtension()
234	for i := len(finalPoly.Coeffs) - 1; i >= 0; i-- {
235		ret = f.gl.MulAddExtension(ret, point, finalPoly.Coeffs[i])
236	}
237	return ret
238}
239
240func (f *Chip) interpolate(
241	x gl.QuadraticExtensionVariable,
242	xPoints []gl.QuadraticExtensionVariable,
243	yPoints []gl.QuadraticExtensionVariable,
244	barycentricWeights []gl.QuadraticExtensionVariable,
245) gl.QuadraticExtensionVariable {
246	if len(xPoints) != len(yPoints) || len(xPoints) != len(barycentricWeights) {
247		panic("length of xPoints, yPoints, and barycentricWeights are inconsistent")
248	}
249
250	lX := gl.OneExtension()
251	for i := 0; i < len(xPoints); i++ {
252		lX = f.gl.SubMulExtension(x, xPoints[i], lX)
253	}
254
255	sum := gl.ZeroExtension()
256	for i := 0; i < len(xPoints); i++ {
257		sum = f.gl.AddExtension(
258			f.gl.MulExtension(
259				f.gl.DivExtension(
260					barycentricWeights[i],
261					f.gl.SubExtension(
262						x,
263						xPoints[i],
264					),
265				),
266				yPoints[i],
267			),
268			sum,
269		)
270	}
271
272	interpolation := f.gl.MulExtension(lX, sum)
273
274	returnField := interpolation
275	// Now check if x is already within the xPoints
276	for i := 0; i < len(xPoints); i++ {
277		returnField = f.gl.Lookup(
278			f.gl.IsZero(f.gl.SubExtension(x, xPoints[i])),
279			returnField,
280			yPoints[i],
281		)
282	}
283
284	return returnField
285}
286
287func (f *Chip) computeEvaluation(
288	x gl.Variable,
289	xIndexWithinCosetBits []frontend.Variable,
290	arityBits uint64,
291	evals []gl.QuadraticExtensionVariable,
292	beta gl.QuadraticExtensionVariable,
293) gl.QuadraticExtensionVariable {
294	arity := 1 << arityBits
295	if (len(evals)) != arity {
296		panic("len(evals) ! arity")
297	}
298	if arityBits > 8 {
299		panic("currently assuming that arityBits is <= 8")
300	}
301
302	g := gl.PrimitiveRootOfUnity(arityBits)
303	gInv := goldilocks.NewElement(0)
304	gInv.Exp(g, big.NewInt(int64(arity-1)))
305
306	// The evaluation vector needs to be reordered first.  Permute the evals array such that each
307	// element's new index is the bit reverse of it's original index.
308	// TODO:  Optimization - Since the size of the evals array should be constant (e.g. 2^arityBits),
309	//        we can just hard code the permutation.
310	permutedEvals := make([]gl.QuadraticExtensionVariable, len(evals))
311	for i := uint8(0); i < uint8(len(evals)); i++ {
312		newIndex := bits.Reverse8(i) >> arityBits
313		permutedEvals[newIndex] = evals[i]
314	}
315
316	// Want `g^(arity - rev_x_index_within_coset)` as in the out-of-circuit version. Compute it
317	// as `(g^-1)^rev_x_index_within_coset`.
318	revXIndexWithinCosetBits := make([]frontend.Variable, len(xIndexWithinCosetBits))
319	for i := 0; i < len(xIndexWithinCosetBits); i++ {
320		revXIndexWithinCosetBits[len(xIndexWithinCosetBits)-1-i] = xIndexWithinCosetBits[i]
321	}
322	start := f.expFromBitsConstBase(gInv, revXIndexWithinCosetBits)
323	cosetStart := f.gl.Mul(start, x)
324
325	xPoints := make([]gl.QuadraticExtensionVariable, len(evals))
326	yPoints := permutedEvals
327
328	// TODO: Make g_F a constant
329	g_F := gl.NewVariable(g.Uint64()).ToQuadraticExtension()
330	xPoints[0] = gl.QuadraticExtensionVariable{cosetStart, gl.Zero()}
331	for i := 1; i < len(evals); i++ {
332		xPoints[i] = f.gl.MulExtension(xPoints[i-1], g_F)
333	}
334
335	// TODO:  This is n^2.  Is there a way to do this better?
336	// Compute the barycentric weights
337	barycentricWeights := make([]gl.QuadraticExtensionVariable, len(xPoints))
338	for i := 0; i < len(xPoints); i++ {
339		barycentricWeights[i] = gl.OneExtension()
340		for j := 0; j < len(xPoints); j++ {
341			if i != j {
342				barycentricWeights[i] = f.gl.SubMulExtension(
343					xPoints[i],
344					xPoints[j],
345					barycentricWeights[i],
346				)
347			}
348		}
349		// Take the inverse of the barycentric weights
350		// TODO: Can provide a witness to this value
351		barycentricWeights[i] = f.gl.InverseExtension(barycentricWeights[i])
352	}
353
354	return f.interpolate(beta, xPoints, yPoints, barycentricWeights)
355}
356
357func (f *Chip) verifyQueryRound(
358	instance InstanceInfo,
359	challenges *variables.FriChallenges,
360	precomputedReducedEval []gl.QuadraticExtensionVariable,
361	initialMerkleCaps []variables.FriMerkleCap,
362	proof *variables.FriProof,
363	xIndex gl.Variable,
364	n uint64,
365	nLog uint64,
366	roundProof *variables.FriQueryRound,
367) {
368	f.assertNoncanonicalIndicesOK()
369	xIndex = f.gl.Reduce(xIndex)
370	xIndexBits := f.api.ToBinary(xIndex.Limb, 64)[0 : f.friParams.DegreeBits+f.friParams.Config.RateBits]
371	capIndexBits := xIndexBits[len(xIndexBits)-int(f.friParams.Config.CapHeight):]
372
373	f.verifyInitialProof(xIndexBits, &roundProof.InitialTreesProof, initialMerkleCaps, capIndexBits)
374
375	subgroupX := f.calculateSubgroupX(
376		xIndexBits,
377		nLog,
378	)
379
380	subgroupX_QE := subgroupX.ToQuadraticExtension()
381
382	oldEval := f.friCombineInitial(
383		instance,
384		roundProof.InitialTreesProof,
385		challenges.FriAlpha,
386		subgroupX_QE,
387		precomputedReducedEval,
388	)
389
390	for i, arityBits := range f.friParams.ReductionArityBits {
391		evals := roundProof.Steps[i].Evals
392
393		cosetIndexBits := xIndexBits[arityBits:]
394		xIndexWithinCosetBits := xIndexBits[:arityBits]
395
396		// Assumes that the arity bits will be 4.  That means that the range of
397		// xIndexWithCoset is [0,2^4-1].  This is based on plonky2's circuit recursive
398		// config:  https://github.com/mir-protocol/plonky2/blob/main/plonky2/src/plonk/circuit_data.rs#L63
399		// Will use a two levels tree of 4-selector gadgets.
400		if arityBits != 4 {
401			panic("assuming arity bits is 4")
402		}
403
404		const NUM_LEAF_LOOKUPS = 4
405		var leafLookups [NUM_LEAF_LOOKUPS]gl.QuadraticExtensionVariable
406		// First create the "leaf" lookup2 circuits
407		// The will use the least significant bits of the xIndexWithCosetBits array
408		for i := 0; i < NUM_LEAF_LOOKUPS; i++ {
409			leafLookups[i] = f.gl.Lookup2(
410				xIndexWithinCosetBits[0],
411				xIndexWithinCosetBits[1],
412				evals[i*NUM_LEAF_LOOKUPS],
413				evals[i*NUM_LEAF_LOOKUPS+1],
414				evals[i*NUM_LEAF_LOOKUPS+2],
415				evals[i*NUM_LEAF_LOOKUPS+3],
416			)
417		}
418
419		// Use the most 2 significant bits of the xIndexWithCosetBits array for the "root" lookup
420		newEval := f.gl.Lookup2(
421			xIndexWithinCosetBits[2],
422			xIndexWithinCosetBits[3],
423			leafLookups[0],
424			leafLookups[1],
425			leafLookups[2],
426			leafLookups[3],
427		)
428
429		f.gl.AssertIsEqual(newEval[0], oldEval[0])
430		f.gl.AssertIsEqual(newEval[1], oldEval[1])
431
432		oldEval = f.computeEvaluation(
433			subgroupX,
434			xIndexWithinCosetBits,
435			arityBits,
436			evals,
437			challenges.FriBetas[i],
438		)
439
440		// Convert evals (array of QE) to fields by taking their 0th degree coefficients
441		fieldEvals := make([]gl.Variable, 0, 2*len(evals))
442		for j := 0; j < len(evals); j++ {
443			fieldEvals = append(fieldEvals, evals[j][0])
444			fieldEvals = append(fieldEvals, evals[j][1])
445		}
446		f.verifyMerkleProofToCapWithCapIndex(
447			fieldEvals,
448			cosetIndexBits,
449			capIndexBits,
450			proof.CommitPhaseMerkleCaps[i],
451			&roundProof.Steps[i].MerkleProof,
452		)
453
454		// Update the point x to x^arity.
455		for j := uint64(0); j < arityBits; j++ {
456			subgroupX = f.gl.Mul(subgroupX, subgroupX)
457		}
458
459		xIndexBits = cosetIndexBits
460	}
461
462	subgroupX_QE = subgroupX.ToQuadraticExtension()
463	finalPolyEval := f.finalPolyEval(proof.FinalPoly, subgroupX_QE)
464
465	f.gl.AssertIsEqual(oldEval[0], finalPolyEval[0])
466	f.gl.AssertIsEqual(oldEval[1], finalPolyEval[1])
467}
468
469func (f *Chip) VerifyFriProof(
470	instance InstanceInfo,
471	openings Openings,
472	friChallenges *variables.FriChallenges,
473	initialMerkleCaps []variables.FriMerkleCap,
474	friProof *variables.FriProof,
475) {
476	// TODO:  Check fri config
477	/* if let Some(max_arity_bits) = params.max_arity_bits() {
478		self.check_recursion_config::<C>(max_arity_bits);
479	}
480
481	debug_assert_eq!(
482		params.final_poly_len(),
483		proof.final_poly.len(),
484		"Final polynomial has wrong degree."
485	); */
486
487	// Check POW
488
489	f.assertLeadingZeros(friChallenges.FriPowResponse, f.friParams.Config)
490
491	precomputedReducedEvals := f.fromOpeningsAndAlpha(&openings, friChallenges.FriAlpha)
492
493	// Size of the LDE domain.
494	nLog := f.friParams.DegreeBits + f.friParams.Config.RateBits
495	n := uint64(math.Pow(2, float64(nLog)))
496
497	if len(friChallenges.FriQueryIndices) != len(friProof.QueryRoundProofs) {
498		panic(fmt.Sprintf(
499			"Number of query indices (%d) should equal number of query round proofs (%d)",
500			len(friChallenges.FriQueryIndices),
501			len(friProof.QueryRoundProofs),
502		))
503	}
504
505	for idx, xIndex := range friChallenges.FriQueryIndices {
506		roundProof := friProof.QueryRoundProofs[idx]
507
508		f.verifyQueryRound(
509			instance,
510			friChallenges,
511			precomputedReducedEvals,
512			initialMerkleCaps,
513			friProof,
514			xIndex,
515			n,
516			nLog,
517			&roundProof,
518		)
519	}
520}