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