Skip to main content

binius_circuits/bignum/
prime_field.rs

1// Copyright 2026 The Binius Developers
2// Copyright 2025 Irreducible Inc.
3use binius_core::word::Word;
4use binius_frontend::{CircuitBuilder, Wire};
5
6use super::{
7	BigUint, BigUintDivideHint, ModDivideHint, ModInverseHint, PseudoMersenneModReduce,
8	addsub::{add_with_carry_out, sub_with_borrow_out},
9	biguint_lt, num_biguint_from_u64_limbs, sub, textbook_mul, textbook_square,
10};
11
12/// A struct that implements prime field arithmetic over pseudo-Mersenne modulus.
13///
14/// Field elements are `BigUint`s consisting of `PseudoMersennePrimeField::limbs_len` limbs.
15/// It is assumed that all field elements are correctly represented (less than modulus).
16pub struct PseudoMersennePrimeField {
17	modulus: BigUint,
18	modulus_po2: usize,
19	modulus_subtrahend: BigUint,
20}
21
22impl PseudoMersennePrimeField {
23	/// Create a new pseudo-Mersenne prime field.
24	///
25	/// See [`PseudoMersenneModReduce`] for description of the parameters.
26	pub fn new(b: &CircuitBuilder, modulus_po2: usize, modulus_subtrahend: &[u64]) -> Self {
27		let modulus_subtrahend = num_biguint_from_u64_limbs(modulus_subtrahend);
28		let po2 = num_bigint::BigUint::from(2usize).pow(modulus_po2 as u32);
29		assert!(po2 > modulus_subtrahend, "2^modulus_po2 - modulus_subtrahend > 0");
30
31		let modulus = BigUint::new_constant(b, &(po2 - &modulus_subtrahend));
32		let modulus_subtrahend = BigUint::new_constant(b, &modulus_subtrahend);
33
34		Self {
35			modulus,
36			modulus_po2,
37			modulus_subtrahend,
38		}
39	}
40
41	/// Number of limbs in `BigUint`s representing field elements.
42	pub const fn limbs_len(&self) -> usize {
43		self.modulus.limbs.len()
44	}
45
46	/// Field modulus.
47	pub const fn modulus(&self) -> &BigUint {
48		&self.modulus
49	}
50
51	/// Field addition.
52	///
53	/// Equivalent formula: `(fe1 + fe2) % modulus`
54	pub fn add(&self, b: &CircuitBuilder, fe1: &BigUint, fe2: &BigUint) -> BigUint {
55		let l = self.limbs_len();
56		assert!(fe1.limbs.len() == l && fe2.limbs.len() == l);
57
58		// Both inputs are reduced, so the sum is below `2 * modulus` and the overflow past
59		// the top limb is a single bit. Carrying it in a wire avoids widening every limbwise
60		// operation below by an extra limb.
61		let (sum, carry) = add_with_carry_out(b, fe1, fe2);
62
63		// A carry means the sum reached `2^(64 * l)`, which is above the modulus, so the
64		// comparison only matters when no carry occurred: with a carry the wrapped sum is
65		// below the modulus. Carry and the comparison are therefore mutually exclusive, so
66		// XOR matches the OR on the MSB that `zero_unless` reads.
67		// TODO: consider nondeterminism
68		let need_reduction = b.bxor(carry, b.bnot(biguint_lt(b, &sum, &self.modulus)));
69
70		// Wrapping subtraction is exact here: without a carry the difference is non-negative,
71		// and with one it borrows out by exactly the `2^(64 * l)` the carry stands for.
72		let (reduced, _borrow) =
73			sub_with_borrow_out(b, &sum, &self.modulus.zero_unless(b, need_reduction));
74
75		reduced
76	}
77
78	/// Field subtraction.
79	///
80	/// Equivalent formula: `(fe1 - fe2) % modulus`
81	pub fn sub(&self, b: &CircuitBuilder, fe1: &BigUint, fe2: &BigUint) -> BigUint {
82		assert!(fe1.limbs.len() == self.limbs_len() && fe2.limbs.len() == self.limbs_len());
83		// NB: fe1 - fe2 = fe1 + modulus - fe2 <= 2*modulus - 1 (one subtraction still normalizes)
84		let fe2_add_inv = sub(b, &self.modulus, fe2);
85		self.add(b, fe1, &fe2_add_inv)
86	}
87
88	/// Field squaring.
89	///
90	/// Equivalent formula: `(fe ** 2) % modulus`
91	pub fn square(&self, b: &CircuitBuilder, fe: &BigUint) -> BigUint {
92		assert!(fe.limbs.len() == self.limbs_len());
93		self.reduce_product(b, &textbook_square(b, fe))
94	}
95
96	/// Field multiplication.
97	///
98	/// Equivalent formula: `(fe1 * fe2) % modulus`
99	/// Note: Both fe1 and fe2 may be greater or equal to modulus.
100	pub fn mul(&self, b: &CircuitBuilder, fe1: &BigUint, fe2: &BigUint) -> BigUint {
101		assert!(fe1.limbs.len() == self.limbs_len() && fe2.limbs.len() == self.limbs_len());
102		self.reduce_product(b, &textbook_mul(b, fe1, fe2))
103	}
104
105	fn reduce_product(&self, b: &CircuitBuilder, product: &BigUint) -> BigUint {
106		let (quotient, remainder) = BigUintDivideHint::call(b, &product.limbs, &self.modulus.limbs);
107
108		let zero = b.add_constant(Word::ZERO);
109
110		let quotient = BigUint { limbs: quotient };
111		let remainder = BigUint { limbs: remainder }.pad_limbs_to(self.limbs_len(), zero);
112
113		b.assert_true("remainder < modulus", biguint_lt(b, &remainder, &self.modulus));
114
115		// constraint: product == remainder + quotient * modulus
116		PseudoMersenneModReduce::new(
117			b,
118			product,
119			self.modulus_po2,
120			&self.modulus_subtrahend,
121			&quotient,
122			&remainder,
123		)
124		.constrain(b);
125
126		remainder
127	}
128
129	/// Field inverse.
130	///
131	/// Equivalent formula (for prime modulus): `(fe1 ** (modulus - 2)) % modulus`
132	/// The wire parameter `exists` is a boolean-wire signifying the existence of the inverse;
133	/// if `exists` is false, the modular reduction constraint is not applied. This is useful
134	/// for avoiding overconstraining in skipped parts of larger circuits.
135	pub fn inverse(&self, b: &CircuitBuilder, fe: &BigUint, exists: Wire) -> BigUint {
136		assert!(fe.limbs.len() == self.limbs_len());
137		let (quotient, inverse) = ModInverseHint::call(b, &fe.limbs, &self.modulus.limbs);
138
139		let zero = b.add_constant(Word::ZERO);
140
141		let quotient = BigUint { limbs: quotient };
142		let inverse = BigUint { limbs: inverse }.pad_limbs_to(self.limbs_len(), zero);
143		let one = BigUint::new_constant(b, &num_bigint::BigUint::from(1usize))
144			.pad_limbs_to(self.limbs_len(), zero);
145
146		let product = textbook_mul(b, &inverse, fe);
147
148		b.assert_true("inverse < modulus", biguint_lt(b, &inverse, &self.modulus));
149
150		// constraint: base * inverse = 1 + quotient * modulus
151		PseudoMersenneModReduce::new(
152			b,
153			&product,
154			self.modulus_po2,
155			&self.modulus_subtrahend,
156			&quotient,
157			&one,
158		)
159		.constrain_cond(b, exists);
160
161		inverse
162	}
163
164	/// Field division.
165	///
166	/// Equivalent formula (for prime modulus): `(dividend * divisor ** (modulus - 2)) % modulus`,
167	/// i.e. `dividend / divisor (mod modulus)`.
168	///
169	/// This collapses a modular inverse followed by a multiplication into a single modular
170	/// reduction, halving the multiplication cost relative to `mul(dividend, inverse(divisor))`.
171	/// The returned `slope` is constrained by `slope * divisor = dividend + quotient * modulus`,
172	/// which is `slope * divisor ≡ dividend (mod modulus)`.
173	///
174	/// The wire parameter `exists` is a boolean-wire signifying the existence of the quotient
175	/// (i.e. that `divisor` is invertible modulo the modulus); if `exists` is false, the modular
176	/// reduction constraint is not applied. This is useful for avoiding overconstraining in
177	/// skipped parts of larger circuits. The `slope < modulus` range check is unconditional, so
178	/// when `exists` is false the returned value is an unconstrained dummy `< modulus`.
179	///
180	/// # Precondition and incompleteness
181	///
182	/// `dividend` must be reduced (`dividend < modulus`): it plays the role of the remainder in
183	/// the reduction `slope * divisor = dividend + quotient * modulus`. If `dividend >= modulus`
184	/// there is no non-negative `quotient` satisfying that equation, so (when `exists` is true)
185	/// the constraint system has no satisfying witness and proof generation fails — the gadget is
186	/// *incomplete* for unreduced dividends. Callers must reduce the dividend beforehand; the
187	/// field helpers [`add`](Self::add), [`sub`](Self::sub), [`mul`](Self::mul) and
188	/// [`square`](Self::square) all return reduced values. (This matches the standard ECDSA
189	/// convention of reducing the message hash to a scalar in `[0, n)` before computing
190	/// `u1`/`u2`.)
191	pub fn div(
192		&self,
193		b: &CircuitBuilder,
194		dividend: &BigUint,
195		divisor: &BigUint,
196		exists: Wire,
197	) -> BigUint {
198		assert!(
199			dividend.limbs.len() == self.limbs_len() && divisor.limbs.len() == self.limbs_len()
200		);
201		let (quotient, slope) =
202			ModDivideHint::call(b, &dividend.limbs, &divisor.limbs, &self.modulus.limbs);
203
204		let zero = b.add_constant(Word::ZERO);
205
206		let quotient = BigUint { limbs: quotient };
207		let slope = BigUint { limbs: slope }.pad_limbs_to(self.limbs_len(), zero);
208
209		let product = textbook_mul(b, &slope, divisor);
210
211		b.assert_true("slope < modulus", biguint_lt(b, &slope, &self.modulus));
212
213		// constraint: slope * divisor = dividend + quotient * modulus
214		PseudoMersenneModReduce::new(
215			b,
216			&product,
217			self.modulus_po2,
218			&self.modulus_subtrahend,
219			&quotient,
220			dividend,
221		)
222		.constrain_cond(b, exists);
223
224		slope
225	}
226}