snarkvm_algorithms/fft/polynomial/
dense.rs

1// Copyright (c) 2019-2025 Provable Inc.
2// This file is part of the snarkVM library.
3
4// Licensed under the Apache License, Version 2.0 (the "License");
5// you may not use this file except in compliance with the License.
6// You may obtain a copy of the License at:
7
8// http://www.apache.org/licenses/LICENSE-2.0
9
10// Unless required by applicable law or agreed to in writing, software
11// distributed under the License is distributed on an "AS IS" BASIS,
12// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13// See the License for the specific language governing permissions and
14// limitations under the License.
15
16//! A polynomial represented in coefficient form.
17
18use super::PolyMultiplier;
19use crate::fft::{EvaluationDomain, Evaluations, Polynomial};
20use snarkvm_fields::{Field, PrimeField};
21use snarkvm_utilities::{cfg_iter_mut, serialize::*};
22
23use anyhow::Result;
24use num_traits::CheckedDiv;
25use rand::Rng;
26use std::{
27    fmt,
28    ops::{Add, AddAssign, Deref, DerefMut, Div, Mul, MulAssign, Neg, Sub, SubAssign},
29};
30
31use itertools::Itertools;
32
33#[cfg(not(feature = "serial"))]
34use rayon::prelude::*;
35
36/// Stores a polynomial in coefficient form.
37#[derive(Clone, PartialEq, Eq, Hash, Default, CanonicalSerialize, CanonicalDeserialize)]
38#[must_use]
39pub struct DensePolynomial<F: Field> {
40    /// The coefficient of `x^i` is stored at location `i` in `self.coeffs`.
41    pub coeffs: Vec<F>,
42}
43
44impl<F: Field> fmt::Debug for DensePolynomial<F> {
45    fn fmt(&self, f: &mut fmt::Formatter) -> Result<(), fmt::Error> {
46        for (i, coeff) in self.coeffs.iter().enumerate().filter(|(_, c)| !c.is_zero()) {
47            if i == 0 {
48                write!(f, "\n{coeff:?}",)?;
49            } else if i == 1 {
50                write!(f, " + \n{coeff:?} * x")?;
51            } else {
52                write!(f, " + \n{coeff:?} * x^{i}")?;
53            }
54        }
55        Ok(())
56    }
57}
58
59impl<F: Field> DensePolynomial<F> {
60    /// Returns the zero polynomial.
61    pub fn zero() -> Self {
62        Self { coeffs: Vec::new() }
63    }
64
65    /// Checks if the given polynomial is zero.
66    pub fn is_zero(&self) -> bool {
67        self.coeffs.is_empty() || self.coeffs.iter().all(|coeff| coeff.is_zero())
68    }
69
70    /// Constructs a new polynomial from a list of coefficients.
71    pub fn from_coefficients_slice(coeffs: &[F]) -> Self {
72        Self::from_coefficients_vec(coeffs.to_vec())
73    }
74
75    /// Constructs a new polynomial from a list of coefficients.
76    pub fn from_coefficients_vec(mut coeffs: Vec<F>) -> Self {
77        // While there are zeros at the end of the coefficient vector, pop them off.
78        while let Some(true) = coeffs.last().map(|c| c.is_zero()) {
79            coeffs.pop();
80        }
81        // Check that either the coefficients vec are empty or that the last coeff is
82        // non-zero.
83        assert!(coeffs.last().map_or(true, |coeff| !coeff.is_zero()));
84
85        Self { coeffs }
86    }
87
88    /// Returns the degree of the polynomial.
89    pub fn degree(&self) -> usize {
90        if self.is_zero() {
91            0
92        } else {
93            assert!(self.coeffs.last().map_or(false, |coeff| !coeff.is_zero()));
94            self.coeffs.len() - 1
95        }
96    }
97
98    /// Evaluates `self` at the given `point` in the field.
99    pub fn evaluate(&self, point: F) -> F {
100        if self.is_zero() {
101            return F::zero();
102        } else if point.is_zero() {
103            return self.coeffs[0];
104        }
105        let mut powers_of_point = Vec::with_capacity(1 + self.degree());
106        powers_of_point.push(F::one());
107        let mut cur = point;
108        for _ in 0..self.degree() {
109            powers_of_point.push(cur);
110            cur *= point;
111        }
112        let zero = F::zero();
113        let mapping = crate::cfg_into_iter!(powers_of_point).zip_eq(&self.coeffs).map(|(power, coeff)| power * coeff);
114        crate::cfg_reduce!(mapping, || zero, |a, b| a + b)
115    }
116
117    /// Outputs a univariate polynomial of degree `d` where each non-leading
118    /// coefficient is sampled uniformly at random from R and the leading
119    /// coefficient is sampled uniformly at random from among the non-zero
120    /// elements of R.
121    pub fn rand<R: Rng>(d: usize, rng: &mut R) -> Self {
122        let mut random_coeffs = (0..(d + 1)).map(|_| F::rand(rng)).collect_vec();
123        while random_coeffs[d].is_zero() {
124            // In the extremely unlikely event, sample again.
125            random_coeffs[d] = F::rand(rng);
126        }
127        Self::from_coefficients_vec(random_coeffs)
128    }
129
130    /// Returns the coefficients of `self`.
131    pub fn coeffs(&self) -> &[F] {
132        &self.coeffs
133    }
134
135    /// Perform a naive n^2 multiplication of `self` by `other`.
136    #[cfg(test)]
137    fn naive_mul(&self, other: &Self) -> Self {
138        if self.is_zero() || other.is_zero() {
139            DensePolynomial::zero()
140        } else {
141            let mut result = vec![F::zero(); self.degree() + other.degree() + 1];
142            for (i, self_coeff) in self.coeffs.iter().enumerate() {
143                for (j, other_coeff) in other.coeffs.iter().enumerate() {
144                    result[i + j] += *self_coeff * other_coeff;
145                }
146            }
147            DensePolynomial::from_coefficients_vec(result)
148        }
149    }
150}
151
152impl<F: PrimeField> DensePolynomial<F> {
153    /// Multiply `self` by the vanishing polynomial for the domain `domain`.
154    pub fn mul_by_vanishing_poly(&self, domain: EvaluationDomain<F>) -> DensePolynomial<F> {
155        let mut shifted = vec![F::zero(); domain.size()];
156        shifted.extend_from_slice(&self.coeffs);
157        crate::cfg_iter_mut!(shifted[..self.coeffs.len()]).zip_eq(&self.coeffs).for_each(|(s, c)| *s -= c);
158        DensePolynomial::from_coefficients_vec(shifted)
159    }
160
161    /// Divide `self` by the vanishing polynomial for the domain `domain`.
162    /// Returns the quotient and remainder of the division.
163    pub fn divide_by_vanishing_poly(
164        &self,
165        domain: EvaluationDomain<F>,
166    ) -> Result<(DensePolynomial<F>, DensePolynomial<F>)> {
167        let self_poly = Polynomial::from(self);
168        let vanishing_poly = Polynomial::from(domain.vanishing_polynomial());
169        self_poly.divide_with_q_and_r(&vanishing_poly)
170    }
171
172    /// Evaluate `self` over `domain`.
173    pub fn evaluate_over_domain_by_ref(&self, domain: EvaluationDomain<F>) -> Evaluations<F> {
174        let poly: Polynomial<'_, F> = self.into();
175        Polynomial::<F>::evaluate_over_domain(poly, domain)
176    }
177
178    /// Evaluate `self` over `domain`.
179    pub fn evaluate_over_domain(self, domain: EvaluationDomain<F>) -> Evaluations<F> {
180        let poly: Polynomial<'_, F> = self.into();
181        Polynomial::<F>::evaluate_over_domain(poly, domain)
182    }
183}
184
185impl<F: Field> From<super::SparsePolynomial<F>> for DensePolynomial<F> {
186    fn from(other: super::SparsePolynomial<F>) -> Self {
187        let mut result = vec![F::zero(); other.degree() + 1];
188        for (i, coeff) in other.coeffs() {
189            result[*i] = *coeff;
190        }
191        DensePolynomial::from_coefficients_vec(result)
192    }
193}
194
195impl<'a, F: Field> Add<&'a DensePolynomial<F>> for &'_ DensePolynomial<F> {
196    type Output = DensePolynomial<F>;
197
198    fn add(self, other: &'a DensePolynomial<F>) -> DensePolynomial<F> {
199        let mut result = if self.is_zero() {
200            other.clone()
201        } else if other.is_zero() {
202            self.clone()
203        } else if self.degree() >= other.degree() {
204            let mut result = self.clone();
205            // Zip safety: `result` and `other` could have different lengths.
206            cfg_iter_mut!(result.coeffs).zip(&other.coeffs).for_each(|(a, b)| *a += b);
207            result
208        } else {
209            let mut result = other.clone();
210            // Zip safety: `result` and `other` could have different lengths.
211            cfg_iter_mut!(result.coeffs).zip(&self.coeffs).for_each(|(a, b)| *a += b);
212            result
213        };
214        // If the leading coefficient ends up being zero, pop it off.
215        while let Some(true) = result.coeffs.last().map(|c| c.is_zero()) {
216            result.coeffs.pop();
217        }
218        result
219    }
220}
221
222impl<'a, F: Field> AddAssign<&'a DensePolynomial<F>> for DensePolynomial<F> {
223    fn add_assign(&mut self, other: &'a DensePolynomial<F>) {
224        if self.is_zero() {
225            self.coeffs.clear();
226            self.coeffs.extend_from_slice(&other.coeffs);
227        } else if other.is_zero() {
228            // return
229        } else if self.degree() >= other.degree() {
230            // Zip safety: `self` and `other` could have different lengths.
231            cfg_iter_mut!(self.coeffs).zip(&other.coeffs).for_each(|(a, b)| *a += b);
232        } else {
233            // Add the necessary number of zero coefficients.
234            self.coeffs.resize(other.coeffs.len(), F::zero());
235            // Zip safety: `self` and `other` have the same length.
236            cfg_iter_mut!(self.coeffs).zip(&other.coeffs).for_each(|(a, b)| *a += b);
237        }
238        // If the leading coefficient ends up being zero, pop it off.
239        while let Some(true) = self.coeffs.last().map(|c| c.is_zero()) {
240            self.coeffs.pop();
241        }
242    }
243}
244
245impl<'a, F: Field> AddAssign<&'a Polynomial<'a, F>> for DensePolynomial<F> {
246    fn add_assign(&mut self, other: &'a Polynomial<F>) {
247        match other {
248            Polynomial::Sparse(p) => *self += &Self::from(p.clone().into_owned()),
249            Polynomial::Dense(p) => *self += p.as_ref(),
250        }
251    }
252}
253
254impl<'a, F: Field> AddAssign<(F, &'a Polynomial<'a, F>)> for DensePolynomial<F> {
255    fn add_assign(&mut self, (f, other): (F, &'a Polynomial<F>)) {
256        match other {
257            Polynomial::Sparse(p) => *self += (f, &Self::from(p.clone().into_owned())),
258            Polynomial::Dense(p) => *self += (f, p.as_ref()),
259        }
260    }
261}
262
263impl<'a, F: Field> AddAssign<(F, &'a DensePolynomial<F>)> for DensePolynomial<F> {
264    #[allow(clippy::suspicious_op_assign_impl)]
265    fn add_assign(&mut self, (f, other): (F, &'a DensePolynomial<F>)) {
266        if self.is_zero() {
267            self.coeffs.clear();
268            self.coeffs.extend_from_slice(&other.coeffs);
269            self.coeffs.iter_mut().for_each(|c| *c *= &f);
270        } else if other.is_zero() {
271            // return
272        } else if self.degree() >= other.degree() {
273            // Zip safety: `self` and `other` could have different lengths.
274            cfg_iter_mut!(self.coeffs).zip(&other.coeffs).for_each(|(a, b)| {
275                *a += f * b;
276            });
277        } else {
278            // Add the necessary number of zero coefficients.
279            self.coeffs.resize(other.coeffs.len(), F::zero());
280            // Zip safety: `self` and `other` have the same length after the resize.
281            cfg_iter_mut!(self.coeffs).zip(&other.coeffs).for_each(|(a, b)| {
282                *a += f * b;
283            });
284        }
285        // If the leading coefficient ends up being zero, pop it off.
286        while let Some(true) = self.coeffs.last().map(|c| c.is_zero()) {
287            self.coeffs.pop();
288        }
289    }
290}
291
292impl<F: Field> Neg for DensePolynomial<F> {
293    type Output = DensePolynomial<F>;
294
295    #[inline]
296    fn neg(mut self) -> DensePolynomial<F> {
297        for coeff in &mut self.coeffs {
298            *coeff = -*coeff;
299        }
300        self
301    }
302}
303
304impl<'a, F: Field> Sub<&'a DensePolynomial<F>> for &'_ DensePolynomial<F> {
305    type Output = DensePolynomial<F>;
306
307    #[inline]
308    fn sub(self, other: &'a DensePolynomial<F>) -> DensePolynomial<F> {
309        let mut result = if self.is_zero() {
310            let mut result = other.clone();
311            for coeff in &mut result.coeffs {
312                *coeff = -(*coeff);
313            }
314            result
315        } else if other.is_zero() {
316            self.clone()
317        } else if self.degree() >= other.degree() {
318            let mut result = self.clone();
319            // Zip safety: `result` and `other` could have different degrees.
320            cfg_iter_mut!(result.coeffs).zip(&other.coeffs).for_each(|(a, b)| *a -= b);
321            result
322        } else {
323            let mut result = self.clone();
324            result.coeffs.resize(other.coeffs.len(), F::zero());
325            // Zip safety: `result` and `other` have the same length after the resize.
326            cfg_iter_mut!(result.coeffs).zip(&other.coeffs).for_each(|(a, b)| {
327                *a -= b;
328            });
329            result
330        };
331        // If the leading coefficient ends up being zero, pop it off.
332        while let Some(true) = result.coeffs.last().map(|c| c.is_zero()) {
333            result.coeffs.pop();
334        }
335        result
336    }
337}
338
339impl<'a, F: Field> SubAssign<&'a DensePolynomial<F>> for DensePolynomial<F> {
340    #[inline]
341    fn sub_assign(&mut self, other: &'a DensePolynomial<F>) {
342        if self.is_zero() {
343            self.coeffs.resize(other.coeffs.len(), F::zero());
344            for (i, coeff) in other.coeffs.iter().enumerate() {
345                self.coeffs[i] -= coeff;
346            }
347        } else if other.is_zero() {
348            // return
349        } else if self.degree() >= other.degree() {
350            // Zip safety: self and other could have different lengths.
351            cfg_iter_mut!(self.coeffs).zip(&other.coeffs).for_each(|(a, b)| *a -= b);
352        } else {
353            // Add the necessary number of zero coefficients.
354            self.coeffs.resize(other.coeffs.len(), F::zero());
355            // Zip safety: self and other have the same length after the resize.
356            cfg_iter_mut!(self.coeffs).zip(&other.coeffs).for_each(|(a, b)| *a -= b);
357        }
358        // If the leading coefficient ends up being zero, pop it off.
359        while let Some(true) = self.coeffs.last().map(|c| c.is_zero()) {
360            self.coeffs.pop();
361        }
362    }
363}
364
365impl<'a, F: Field> AddAssign<&'a super::SparsePolynomial<F>> for DensePolynomial<F> {
366    #[inline]
367    fn add_assign(&mut self, other: &'a super::SparsePolynomial<F>) {
368        if self.degree() < other.degree() {
369            self.coeffs.resize(other.degree() + 1, F::zero());
370        }
371        for (i, b) in other.coeffs() {
372            self.coeffs[*i] += b;
373        }
374        // If the leading coefficient ends up being zero, pop it off.
375        while let Some(true) = self.coeffs.last().map(|c| c.is_zero()) {
376            self.coeffs.pop();
377        }
378    }
379}
380
381impl<'a, F: Field> Sub<&'a super::SparsePolynomial<F>> for DensePolynomial<F> {
382    type Output = Self;
383
384    #[inline]
385    fn sub(mut self, other: &'a super::SparsePolynomial<F>) -> Self::Output {
386        if self.degree() < other.degree() {
387            self.coeffs.resize(other.degree() + 1, F::zero());
388        }
389        for (i, b) in other.coeffs() {
390            self.coeffs[*i] -= b;
391        }
392        // If the leading coefficient ends up being zero, pop it off.
393        while let Some(true) = self.coeffs.last().map(|c| c.is_zero()) {
394            self.coeffs.pop();
395        }
396        self
397    }
398}
399
400impl<'a, F: Field> Div<&'a DensePolynomial<F>> for &'_ DensePolynomial<F> {
401    type Output = DensePolynomial<F>;
402
403    /// This division can panic and ignores remainders
404    #[inline]
405    fn div(self, divisor: &'a DensePolynomial<F>) -> DensePolynomial<F> {
406        let a: Polynomial<_> = self.into();
407        let b: Polynomial<_> = divisor.into();
408        a.divide_with_q_and_r(&b).expect("division failed").0
409    }
410}
411
412impl<F: Field> Div<DensePolynomial<F>> for DensePolynomial<F> {
413    type Output = DensePolynomial<F>;
414
415    /// This division can panic and ignores remainders
416    #[inline]
417    fn div(self, divisor: DensePolynomial<F>) -> DensePolynomial<F> {
418        let a: Polynomial<_> = self.into();
419        let b: Polynomial<_> = divisor.into();
420        a.divide_with_q_and_r(&b).expect("division failed").0
421    }
422}
423
424impl<F: Field> CheckedDiv for DensePolynomial<F> {
425    #[inline]
426    fn checked_div(&self, divisor: &DensePolynomial<F>) -> Option<DensePolynomial<F>> {
427        let a: Polynomial<_> = self.into();
428        let b: Polynomial<_> = divisor.into();
429        match a.divide_with_q_and_r(&b) {
430            Ok((divisor, remainder)) => {
431                if remainder.is_zero() {
432                    Some(divisor)
433                } else {
434                    None
435                }
436            }
437            Err(_) => None,
438        }
439    }
440}
441
442/// Performs O(nlogn) multiplication of polynomials if F is smooth.
443impl<'a, F: PrimeField> Mul<&'a DensePolynomial<F>> for &'_ DensePolynomial<F> {
444    type Output = DensePolynomial<F>;
445
446    #[inline]
447    #[allow(clippy::suspicious_arithmetic_impl)]
448    fn mul(self, other: &'a DensePolynomial<F>) -> DensePolynomial<F> {
449        if self.is_zero() || other.is_zero() {
450            DensePolynomial::zero()
451        } else {
452            let mut m = PolyMultiplier::new();
453            m.add_polynomial_ref(self, "");
454            m.add_polynomial_ref(other, "");
455            m.multiply().unwrap()
456        }
457    }
458}
459
460/// Multiplies `self` by `other: F`.
461impl<F: Field> Mul<F> for DensePolynomial<F> {
462    type Output = Self;
463
464    #[inline]
465    fn mul(mut self, other: F) -> Self {
466        self.iter_mut().for_each(|c| *c *= other);
467        self
468    }
469}
470
471/// Multiplies `self` by `other: F`.
472impl<F: Field> Mul<F> for &'_ DensePolynomial<F> {
473    type Output = DensePolynomial<F>;
474
475    #[inline]
476    fn mul(self, other: F) -> Self::Output {
477        let result = self.clone();
478        result * other
479    }
480}
481
482/// Multiplies `self` by `other: F`.
483impl<F: Field> MulAssign<F> for DensePolynomial<F> {
484    #[allow(clippy::suspicious_arithmetic_impl)]
485    fn mul_assign(&mut self, other: F) {
486        cfg_iter_mut!(self).for_each(|c| *c *= other);
487    }
488}
489
490/// Multiplies `self` by `other: F`.
491impl<F: Field> std::iter::Sum for DensePolynomial<F> {
492    fn sum<I: Iterator<Item = Self>>(iter: I) -> Self {
493        iter.fold(DensePolynomial::zero(), |a, b| &a + &b)
494    }
495}
496
497impl<F: Field> Deref for DensePolynomial<F> {
498    type Target = [F];
499
500    fn deref(&self) -> &[F] {
501        &self.coeffs
502    }
503}
504
505impl<F: Field> DerefMut for DensePolynomial<F> {
506    fn deref_mut(&mut self) -> &mut [F] {
507        &mut self.coeffs
508    }
509}
510
511#[cfg(test)]
512mod tests {
513    use crate::fft::polynomial::*;
514    use num_traits::CheckedDiv;
515    use snarkvm_curves::bls12_377::Fr;
516    use snarkvm_fields::{Field, One, Zero};
517    use snarkvm_utilities::rand::{TestRng, Uniform};
518
519    use rand::RngCore;
520
521    #[test]
522    fn double_polynomials_random() {
523        let rng = &mut TestRng::default();
524        for degree in 0..70 {
525            let p = DensePolynomial::<Fr>::rand(degree, rng);
526            let p_double = &p + &p;
527            let p_quad = &p_double + &p_double;
528            assert_eq!(&(&(&p + &p) + &p) + &p, p_quad);
529        }
530    }
531
532    #[test]
533    fn add_polynomials() {
534        let rng = &mut TestRng::default();
535        for a_degree in 0..70 {
536            for b_degree in 0..70 {
537                let p1 = DensePolynomial::<Fr>::rand(a_degree, rng);
538                let p2 = DensePolynomial::<Fr>::rand(b_degree, rng);
539                let res1 = &p1 + &p2;
540                let res2 = &p2 + &p1;
541                assert_eq!(res1, res2);
542            }
543        }
544    }
545
546    #[test]
547    fn add_polynomials_with_mul() {
548        let rng = &mut TestRng::default();
549        for a_degree in 0..70 {
550            for b_degree in 0..70 {
551                let mut p1 = DensePolynomial::rand(a_degree, rng);
552                let p2 = DensePolynomial::rand(b_degree, rng);
553                let f = Fr::rand(rng);
554                let f_p2 = DensePolynomial::from_coefficients_vec(p2.coeffs.iter().map(|c| f * c).collect());
555                let res2 = &f_p2 + &p1;
556                p1 += (f, &p2);
557                let res1 = p1;
558                assert_eq!(res1, res2);
559            }
560        }
561    }
562
563    #[test]
564    fn sub_polynomials() {
565        let rng = &mut TestRng::default();
566        let p1 = DensePolynomial::<Fr>::rand(5, rng);
567        let p2 = DensePolynomial::<Fr>::rand(3, rng);
568        let res1 = &p1 - &p2;
569        let res2 = &p2 - &p1;
570        assert_eq!(&res1 + &p2, p1, "Subtraction should be inverse of addition!");
571        assert_eq!(res1, -res2, "p2 - p1 = -(p1 - p2)");
572    }
573
574    #[test]
575    fn divide_polynomials_fixed() {
576        let dividend = DensePolynomial::from_coefficients_slice(&[
577            "4".parse().unwrap(),
578            "8".parse().unwrap(),
579            "5".parse().unwrap(),
580            "1".parse().unwrap(),
581        ]);
582        let divisor = DensePolynomial::from_coefficients_slice(&[Fr::one(), Fr::one()]); // Construct a monic linear polynomial.
583        let result = dividend.checked_div(&divisor).unwrap();
584        let expected_result = DensePolynomial::from_coefficients_slice(&[
585            "4".parse().unwrap(),
586            "4".parse().unwrap(),
587            "1".parse().unwrap(),
588        ]);
589        assert_eq!(expected_result, result);
590    }
591
592    #[test]
593    #[allow(clippy::needless_borrow)]
594    fn divide_polynomials_random() {
595        let rng = &mut TestRng::default();
596
597        for a_degree in 0..70 {
598            for b_degree in 0..70 {
599                let dividend = DensePolynomial::<Fr>::rand(a_degree, rng);
600                let divisor = DensePolynomial::<Fr>::rand(b_degree, rng);
601                let (quotient, remainder) =
602                    Polynomial::divide_with_q_and_r(&(&dividend).into(), &(&divisor).into()).unwrap();
603                assert_eq!(dividend, &(&divisor * &quotient) + &remainder)
604            }
605        }
606    }
607
608    #[test]
609    fn evaluate_polynomials() {
610        let rng = &mut TestRng::default();
611        for a_degree in 0..70 {
612            let p = DensePolynomial::rand(a_degree, rng);
613            let point: Fr = Fr::from(10u64);
614            let mut total = Fr::zero();
615            for (i, coeff) in p.coeffs.iter().enumerate() {
616                total += point.pow([i as u64]) * coeff;
617            }
618            assert_eq!(p.evaluate(point), total);
619        }
620    }
621
622    #[test]
623    fn divide_poly_by_zero() {
624        let a = Polynomial::<Fr>::zero();
625        let b = Polynomial::<Fr>::zero();
626        assert!(a.divide_with_q_and_r(&b).is_err());
627    }
628
629    #[test]
630    fn mul_polynomials_random() {
631        let rng = &mut TestRng::default();
632        for a_degree in 0..70 {
633            for b_degree in 0..70 {
634                dbg!(a_degree);
635                dbg!(b_degree);
636                let a = DensePolynomial::<Fr>::rand(a_degree, rng);
637                let b = DensePolynomial::<Fr>::rand(b_degree, rng);
638                assert_eq!(&a * &b, a.naive_mul(&b))
639            }
640        }
641    }
642
643    #[test]
644    fn mul_polynomials_n_random() {
645        let rng = &mut TestRng::default();
646
647        let max_degree = 1 << 8;
648
649        for _ in 0..10 {
650            let mut multiplier = PolyMultiplier::new();
651            let a = DensePolynomial::<Fr>::rand(max_degree / 2, rng);
652            let mut mul_degree = a.degree() + 1;
653            multiplier.add_polynomial(a.clone(), "a");
654            let mut naive = a.clone();
655
656            // Include polynomials and evaluations
657            let num_polys = (rng.next_u32() as usize) % 8;
658            let num_evals = (rng.next_u32() as usize) % 4;
659            println!("\nnum_polys {num_polys}, num_evals {num_evals}");
660
661            for _ in 1..num_polys {
662                let degree = (rng.next_u32() as usize) % max_degree;
663                mul_degree += degree + 1;
664                println!("poly degree {degree}");
665                let a = DensePolynomial::<Fr>::rand(degree, rng);
666                naive = naive.naive_mul(&a);
667                multiplier.add_polynomial(a.clone(), "a");
668            }
669
670            // Add evaluations but don't overflow the domain
671            let mut eval_degree = mul_degree;
672            let domain = EvaluationDomain::new(mul_degree).unwrap();
673            println!("mul_degree {}, domain {}", mul_degree, domain.size());
674            for _ in 0..num_evals {
675                let a = DensePolynomial::<Fr>::rand(mul_degree / 8, rng);
676                eval_degree += a.len() + 1;
677                if eval_degree < domain.size() {
678                    println!("eval degree {eval_degree}");
679                    let mut a_evals = a.clone().coeffs;
680                    domain.fft_in_place(&mut a_evals);
681                    let a_evals = Evaluations::from_vec_and_domain(a_evals, domain);
682
683                    naive = naive.naive_mul(&a);
684                    multiplier.add_evaluation(a_evals, "a");
685                }
686            }
687
688            assert_eq!(multiplier.multiply().unwrap(), naive);
689        }
690    }
691
692    #[test]
693    fn mul_polynomials_corner_cases() {
694        let rng = &mut TestRng::default();
695
696        let a_degree = 70;
697
698        // Single polynomial
699        println!("Test single polynomial");
700        let a = DensePolynomial::<Fr>::rand(a_degree, rng);
701        let mut multiplier = PolyMultiplier::new();
702        multiplier.add_polynomial(a.clone(), "a");
703        assert_eq!(multiplier.multiply().unwrap(), a);
704
705        // Note PolyMultiplier doesn't support evaluations with no polynomials
706    }
707
708    #[test]
709    fn mul_by_vanishing_poly() {
710        let rng = &mut TestRng::default();
711        for size in 1..10 {
712            let domain = EvaluationDomain::new(1 << size).unwrap();
713            for degree in 0..70 {
714                let p = DensePolynomial::<Fr>::rand(degree, rng);
715                let ans1 = p.mul_by_vanishing_poly(domain);
716                let ans2 = &p * &domain.vanishing_polynomial().into();
717                assert_eq!(ans1, ans2);
718            }
719        }
720    }
721}