winter_air/air/
divisor.rs

1// Copyright (c) Facebook, Inc. and its affiliates.
2//
3// This source code is licensed under the MIT license found in the
4// LICENSE file in the root directory of this source tree.
5
6use alloc::vec::Vec;
7use core::fmt::{Display, Formatter};
8
9use math::{FieldElement, StarkField};
10
11use crate::air::Assertion;
12
13// CONSTRAINT DIVISOR
14// ================================================================================================
15/// The denominator portion of boundary and transition constraints.
16///
17/// A divisor is described by a combination of a sparse polynomial, which describes the numerator
18/// of the divisor and a set of exemption points, which describe the denominator of the divisor.
19/// The numerator polynomial is described as multiplication of tuples where each tuple encodes
20/// an expression $(x^a - b)$. The exemption points encode expressions $(x - a)$.
21///
22/// For example divisor $(x^a - 1) \cdot (x^b - 2) / (x - 3)$ can be represented as:
23/// numerator: `[(a, 1), (b, 2)]`, exemptions: `[3]`.
24///
25/// A divisor cannot be instantiated directly, and instead must be created either for an
26/// [Assertion] or for a transition constraint.
27#[derive(Clone, Debug, PartialEq, Eq)]
28pub struct ConstraintDivisor<B: StarkField> {
29    pub(super) numerator: Vec<(usize, B)>,
30    pub(super) exemptions: Vec<B>,
31}
32
33impl<B: StarkField> ConstraintDivisor<B> {
34    // CONSTRUCTORS
35    // --------------------------------------------------------------------------------------------
36
37    /// Returns a new divisor instantiated from the provided parameters.
38    fn new(numerator: Vec<(usize, B)>, exemptions: Vec<B>) -> Self {
39        ConstraintDivisor { numerator, exemptions }
40    }
41
42    /// Builds a divisor for transition constraints.
43    ///
44    /// For transition constraints, the divisor polynomial $z(x)$ is always the same:
45    ///
46    /// $$ z(x) = \frac{x^n - 1}{ \prod_{i=1}^k (x - g^{n-i})} $$
47    ///
48    /// where, $n$ is the length of the execution trace, $g$ is the generator of the trace domain,
49    /// and $k$ is the number of exemption points. The default value for $k$ is $1$.
50    ///
51    /// The above divisor specifies that transition constraints must hold on all steps of the
52    /// constraint enforcement domain except for the last $k$ steps.
53    pub fn from_transition(
54        constraint_enforcement_domain_size: usize,
55        num_exemptions: usize,
56    ) -> Self {
57        let exemptions = (constraint_enforcement_domain_size - num_exemptions
58            ..constraint_enforcement_domain_size)
59            .map(|step| get_trace_domain_value_at::<B>(constraint_enforcement_domain_size, step))
60            .collect();
61        Self::new(vec![(constraint_enforcement_domain_size, B::ONE)], exemptions)
62    }
63
64    /// Builds a divisor for a boundary constraint described by the assertion.
65    ///
66    /// For boundary constraints, the divisor polynomial is defined as:
67    ///
68    /// $$
69    /// z(x) = x^k - g^{a \cdot k}
70    /// $$
71    ///
72    /// where $g$ is the generator of the trace domain, $k$ is the number of asserted steps, and
73    /// $a$ is the step offset in the trace domain. Specifically:
74    /// * For an assertion against a single step, the polynomial is $(x - g^a)$, where $a$ is the
75    ///   step on which the assertion should hold.
76    /// * For an assertion against a sequence of steps which fall on powers of two, it is
77    ///   $(x^k - 1)$ where $k$ is the number of asserted steps.
78    /// * For assertions against a sequence of steps which repeat with a period that is a power
79    ///   of two but don't fall exactly on steps which are powers of two (e.g. 1, 9, 17, ... )
80    ///   it is $(x^k - g^{a \cdot k})$, where $a$ is the number of steps by which the assertion steps
81    ///   deviate from a power of two, and $k$ is the number of asserted steps. This is equivalent to
82    ///   $(x - g^a) \cdot (x - g^{a + j}) \cdot (x - g^{a + 2 \cdot j}) ... (x - g^{a + (k  - 1) \cdot j})$,
83    ///   where $j$ is the length of interval between asserted steps (e.g. 8).
84    ///
85    /// # Panics
86    /// Panics of the specified `trace_length` is inconsistent with the specified `assertion`.
87    pub fn from_assertion<E>(assertion: &Assertion<E>, trace_length: usize) -> Self
88    where
89        E: FieldElement<BaseField = B>,
90    {
91        let num_steps = assertion.get_num_steps(trace_length);
92        if assertion.first_step == 0 {
93            Self::new(vec![(num_steps, B::ONE)], vec![])
94        } else {
95            let trace_offset = num_steps * assertion.first_step;
96            let offset = get_trace_domain_value_at::<B>(trace_length, trace_offset);
97            Self::new(vec![(num_steps, offset)], vec![])
98        }
99    }
100
101    // PUBLIC ACCESSORS
102    // --------------------------------------------------------------------------------------------
103
104    /// Returns the numerator portion of this constraint divisor.
105    pub fn numerator(&self) -> &[(usize, B)] {
106        &self.numerator
107    }
108
109    /// Returns exemption points (the denominator portion) of this constraints divisor.
110    pub fn exemptions(&self) -> &[B] {
111        &self.exemptions
112    }
113
114    /// Returns the degree of the divisor polynomial
115    pub fn degree(&self) -> usize {
116        let numerator_degree = self.numerator.iter().fold(0, |degree, term| degree + term.0);
117        let denominator_degree = self.exemptions.len();
118        numerator_degree - denominator_degree
119    }
120
121    // EVALUATORS
122    // --------------------------------------------------------------------------------------------
123    /// Evaluates the divisor polynomial at the provided `x` coordinate.
124    pub fn evaluate_at<E: FieldElement<BaseField = B>>(&self, x: E) -> E {
125        // compute the numerator value
126        let mut numerator = E::ONE;
127        for (degree, constant) in self.numerator.iter() {
128            let v = x.exp((*degree as u32).into());
129            let v = v - E::from(*constant);
130            numerator *= v;
131        }
132
133        // compute the denominator value
134        let denominator = self.evaluate_exemptions_at(x);
135
136        numerator / denominator
137    }
138
139    /// Evaluates the denominator of this divisor (the exemption points) at the provided `x`
140    /// coordinate.
141    #[inline(always)]
142    pub fn evaluate_exemptions_at<E: FieldElement<BaseField = B>>(&self, x: E) -> E {
143        self.exemptions.iter().fold(E::ONE, |r, &e| r * (x - E::from(e)))
144    }
145}
146
147impl<B: StarkField> Display for ConstraintDivisor<B> {
148    fn fmt(&self, f: &mut Formatter) -> core::fmt::Result {
149        for (degree, offset) in self.numerator.iter() {
150            write!(f, "(x^{degree} - {offset})")?;
151        }
152        if !self.exemptions.is_empty() {
153            write!(f, " / ")?;
154            for x in self.exemptions.iter() {
155                write!(f, "(x - {x})")?;
156            }
157        }
158        Ok(())
159    }
160}
161
162// HELPER FUNCTIONS
163// ================================================================================================
164
165/// Returns g^step, where g is the generator of trace domain.
166fn get_trace_domain_value_at<B: StarkField>(trace_length: usize, step: usize) -> B {
167    debug_assert!(step < trace_length, "step must be in the trace domain [0, {trace_length})");
168    let g = B::get_root_of_unity(trace_length.ilog2());
169    g.exp((step as u64).into())
170}
171
172// TESTS
173// ================================================================================================
174
175#[cfg(test)]
176mod tests {
177    use math::{fields::f128::BaseElement, polynom};
178
179    use super::*;
180
181    #[test]
182    fn constraint_divisor_degree() {
183        // single term numerator
184        let div = ConstraintDivisor::new(vec![(4, BaseElement::ONE)], vec![]);
185        assert_eq!(4, div.degree());
186
187        // multi-term numerator
188        let div = ConstraintDivisor::new(
189            vec![(4, BaseElement::ONE), (2, BaseElement::new(2)), (3, BaseElement::new(3))],
190            vec![],
191        );
192        assert_eq!(9, div.degree());
193
194        // multi-term numerator with exemption points
195        let div = ConstraintDivisor::new(
196            vec![(4, BaseElement::ONE), (2, BaseElement::new(2)), (3, BaseElement::new(3))],
197            vec![BaseElement::ONE, BaseElement::new(2)],
198        );
199        assert_eq!(7, div.degree());
200    }
201
202    #[test]
203    fn constraint_divisor_evaluation() {
204        // single term numerator: (x^4 - 1)
205        let div = ConstraintDivisor::new(vec![(4, BaseElement::ONE)], vec![]);
206        assert_eq!(BaseElement::new(15), div.evaluate_at(BaseElement::new(2)));
207
208        // multi-term numerator: (x^4 - 1) * (x^2 - 2) * (x^3 - 3)
209        let div = ConstraintDivisor::new(
210            vec![(4, BaseElement::ONE), (2, BaseElement::new(2)), (3, BaseElement::new(3))],
211            vec![],
212        );
213        let expected = BaseElement::new(15) * BaseElement::new(2) * BaseElement::new(5);
214        assert_eq!(expected, div.evaluate_at(BaseElement::new(2)));
215
216        // multi-term numerator with exemption points:
217        // (x^4 - 1) * (x^2 - 2) * (x^3 - 3) / ((x - 1) * (x - 2))
218        let div = ConstraintDivisor::new(
219            vec![(4, BaseElement::ONE), (2, BaseElement::new(2)), (3, BaseElement::new(3))],
220            vec![BaseElement::ONE, BaseElement::new(2)],
221        );
222        let expected = BaseElement::new(255) * BaseElement::new(14) * BaseElement::new(61)
223            / BaseElement::new(6);
224        assert_eq!(expected, div.evaluate_at(BaseElement::new(4)));
225    }
226
227    #[test]
228    fn constraint_divisor_equivalence() {
229        let n = 8_usize;
230        let g = BaseElement::get_root_of_unity(n.trailing_zeros());
231        let k = 4_u32;
232        let j = n as u32 / k;
233
234        // ----- periodic assertion divisor, no offset --------------------------------------------
235
236        // create a divisor for assertion which repeats every 2 steps starting at step 0
237        let assertion = Assertion::periodic(0, 0, j as usize, BaseElement::ONE);
238        let divisor = ConstraintDivisor::from_assertion(&assertion, n);
239
240        // z(x) = x^4 - 1 = (x - 1) * (x - g^2) * (x - g^4) * (x - g^6)
241        let poly = polynom::mul(
242            &polynom::mul(
243                &[-BaseElement::ONE, BaseElement::ONE],
244                &[-g.exp(j.into()), BaseElement::ONE],
245            ),
246            &polynom::mul(
247                &[-g.exp((2 * j).into()), BaseElement::ONE],
248                &[-g.exp((3 * j).into()), BaseElement::ONE],
249            ),
250        );
251
252        for i in 0..n {
253            let expected = polynom::eval(&poly, g.exp((i as u32).into()));
254            let actual = divisor.evaluate_at(g.exp((i as u32).into()));
255            assert_eq!(expected, actual);
256            if i.is_multiple_of(j as usize) {
257                assert_eq!(BaseElement::ZERO, actual);
258            }
259        }
260
261        // ----- periodic assertion divisor, with offset ------------------------------------------
262
263        // create a divisor for assertion which repeats every 2 steps starting at step 1
264        let offset = 1_u32;
265        let assertion = Assertion::periodic(0, offset as usize, j as usize, BaseElement::ONE);
266        let divisor = ConstraintDivisor::from_assertion(&assertion, n);
267        assert_eq!(ConstraintDivisor::new(vec![(k as usize, g.exp(k.into()))], vec![]), divisor);
268
269        // z(x) = x^4 - g^4 = (x - g) * (x - g^3) * (x - g^5) * (x - g^7)
270        let poly = polynom::mul(
271            &polynom::mul(
272                &[-g.exp(offset.into()), BaseElement::ONE],
273                &[-g.exp((offset + j).into()), BaseElement::ONE],
274            ),
275            &polynom::mul(
276                &[-g.exp((offset + 2 * j).into()), BaseElement::ONE],
277                &[-g.exp((offset + 3 * j).into()), BaseElement::ONE],
278            ),
279        );
280
281        for i in 0..n {
282            let expected = polynom::eval(&poly, g.exp((i as u32).into()));
283            let actual = divisor.evaluate_at(g.exp((i as u32).into()));
284            assert_eq!(expected, actual);
285            if i % (j as usize) == offset as usize {
286                assert_eq!(BaseElement::ZERO, actual);
287            }
288        }
289
290        // create a divisor for assertion which repeats every 4 steps starting at step 3
291        let offset = 3_u32;
292        let k = 2_u32;
293        let j = n as u32 / k;
294        let assertion = Assertion::periodic(0, offset as usize, j as usize, BaseElement::ONE);
295        let divisor = ConstraintDivisor::from_assertion(&assertion, n);
296        assert_eq!(
297            ConstraintDivisor::new(vec![(k as usize, g.exp((offset * k).into()))], vec![]),
298            divisor
299        );
300
301        // z(x) = x^2 - g^6 = (x - g^3) * (x - g^7)
302        let poly = polynom::mul(
303            &[-g.exp(offset.into()), BaseElement::ONE],
304            &[-g.exp((offset + j).into()), BaseElement::ONE],
305        );
306
307        for i in 0..n {
308            let expected = polynom::eval(&poly, g.exp((i as u32).into()));
309            let actual = divisor.evaluate_at(g.exp((i as u32).into()));
310            assert_eq!(expected, actual);
311            if i % (j as usize) == offset as usize {
312                assert_eq!(BaseElement::ZERO, actual);
313            }
314        }
315    }
316}
winter_air/air/divisor.rs

winter_air/air/
divisor.rs
