Skip to main content

feanor_math/rings/multivariate/
mod.rs

1use crate::homomorphism::Homomorphism;
2use crate::ring::*;
3use crate::seq::VectorFn;
4use crate::wrapper::RingElementWrapper;
5
6/// Contains an implementation [`multivariate_impl::MultivariatePolyRingImpl`] of a multivariate
7/// polynomial ring.
8pub mod multivariate_impl;
9
10use std::any::Any;
11use std::cmp::{Ordering, max};
12
13/// Type of coefficients of multivariate polynomials of the given ring.
14pub type PolyCoeff<P> = El<<<P as RingStore>::Type as RingExtension>::BaseRing>;
15
16/// Type of monomials of multivariate polynomials of the given ring.
17pub type PolyMonomial<P> = <<P as RingStore>::Type as MultivariatePolyRing>::Monomial;
18
19/// Trait for multivariate polynomial rings.
20pub trait MultivariatePolyRing: RingExtension {
21    type Monomial;
22    type TermIter<'a>: Iterator<Item = (&'a El<Self::BaseRing>, &'a Self::Monomial)>
23    where
24        Self: 'a;
25
26    /// Returns the number of variables of this polynomial ring, i.e. the transcendence degree
27    /// of the base ring.
28    fn indeterminate_count(&self) -> usize;
29
30    /// Returns the monomial `Xi`, where `Xi` is the `i`-th generator of this ring.
31    fn indeterminate(&self, i: usize) -> Self::Monomial {
32        assert!(i < self.indeterminate_count());
33        self.create_monomial((0..self.indeterminate_count()).map(|j| if i == j { 1 } else { 0 }))
34    }
35
36    /// Creates a monomial with the given exponents.
37    ///
38    /// Note that when building a polynomial, the most convenient method is usually
39    /// to use [`MultivariatePolyRingStore::with_wrapped_indeterminates()`].
40    ///
41    /// # Example
42    /// ```rust
43    /// # use feanor_math::ring::*;
44    /// # use feanor_math::primitive_int::*;
45    /// # use feanor_math::rings::multivariate::*;
46    /// # use feanor_math::rings::multivariate::multivariate_impl::*;
47    /// let poly_ring = MultivariatePolyRingImpl::new(StaticRing::<i64>::RING, 3);
48    /// let x_as_monomial = poly_ring.create_monomial([1, 0, 0]);
49    /// let x_as_poly = poly_ring.create_term(1, x_as_monomial);
50    /// assert_eq!("X0", format!("{}", poly_ring.format(&x_as_poly)));
51    /// ```
52    fn create_monomial<I>(&self, exponents: I) -> Self::Monomial
53    where
54        I: IntoIterator<Item = usize>,
55        I::IntoIter: ExactSizeIterator;
56
57    /// Multiplies the given polynomial with the given monomial.
58    fn mul_assign_monomial(&self, f: &mut Self::Element, monomial: Self::Monomial);
59
60    /// Returns the coefficient corresponding to the given monomial in the given polynomial.
61    /// If the polynomial does not contain a term with that monomial, zero is returned.
62    fn coefficient_at<'a>(&'a self, f: &'a Self::Element, m: &Self::Monomial) -> &'a El<Self::BaseRing>;
63
64    /// Returns the power of the `var_index`-th variable in the given monomial.
65    /// In other words, this maps `X1^i1 ... Xm^im` to `i(var_index)`.
66    fn exponent_at(&self, m: &Self::Monomial, var_index: usize) -> usize;
67
68    /// Writes the powers of each variable in the given monomial into the given
69    /// output slice.
70    ///
71    /// This is equivalent to performing `out[i] = self.exponent_at(m, i)` for
72    /// every `i` in `0..self.indeterminate_count()`.
73    fn expand_monomial_to(&self, m: &Self::Monomial, out: &mut [usize]) {
74        assert_eq!(out.len(), self.indeterminate_count());
75        for (i, out) in out.iter_mut().enumerate() {
76            *out = self.exponent_at(m, i);
77        }
78    }
79
80    /// Returns an iterator over all nonzero terms of the given polynomial.
81    fn terms<'a>(&'a self, f: &'a Self::Element) -> Self::TermIter<'a>;
82
83    /// Creates a new single-term polynomial.
84    fn create_term(&self, coeff: El<Self::BaseRing>, monomial: Self::Monomial) -> Self::Element {
85        let mut result = self.from(coeff);
86        self.mul_assign_monomial(&mut result, monomial);
87        return result;
88    }
89
90    /// Returns the **L**eading **T**erm of `f`, i.e. the term whose monomial is largest w.r.t. the
91    /// given order.
92    fn LT<'a, O: MonomialOrder>(
93        &'a self,
94        f: &'a Self::Element,
95        order: O,
96    ) -> Option<(&'a El<Self::BaseRing>, &'a Self::Monomial)> {
97        self.terms(f).max_by(|l, r| order.compare(RingRef::new(self), l.1, r.1))
98    }
99
100    /// Returns the term of `f` whose monomial is largest (w.r.t. the given order) among all
101    /// monomials smaller than `lt_than`.
102    fn largest_term_lt<'a, O: MonomialOrder>(
103        &'a self,
104        f: &'a Self::Element,
105        order: O,
106        lt_than: &Self::Monomial,
107    ) -> Option<(&'a El<Self::BaseRing>, &'a Self::Monomial)> {
108        self.terms(f)
109            .filter(|(_, m)| order.compare(RingRef::new(self), m, lt_than) == Ordering::Less)
110            .max_by(|l, r| order.compare(RingRef::new(self), l.1, r.1))
111    }
112
113    fn add_assign_from_terms<I>(&self, lhs: &mut Self::Element, rhs: I)
114    where
115        I: IntoIterator<Item = (El<Self::BaseRing>, Self::Monomial)>,
116    {
117        let self_ring = RingRef::new(self);
118        self.add_assign(lhs, self_ring.sum(rhs.into_iter().map(|(c, m)| self.create_term(c, m))));
119    }
120
121    /// Applies the given homomorphism `R -> S` to each coefficient of the given polynomial
122    /// in `R[X1, ..., Xm]` to produce a monomial in `S[X1, ..., Xm]`.
123    fn map_terms<P, H>(&self, from: &P, el: &P::Element, hom: H) -> Self::Element
124    where
125        P: ?Sized + MultivariatePolyRing,
126        H: Homomorphism<<P::BaseRing as RingStore>::Type, <Self::BaseRing as RingStore>::Type>,
127    {
128        assert!(self.base_ring().get_ring() == hom.codomain().get_ring());
129        assert!(from.base_ring().get_ring() == hom.domain().get_ring());
130        assert_eq!(self.indeterminate_count(), from.indeterminate_count());
131        let mut exponents_storage = (0..self.indeterminate_count()).map(|_| 0).collect::<Vec<_>>();
132        return RingRef::new(self).from_terms(from.terms(el).map(|(c, m)| {
133            from.expand_monomial_to(m, &mut exponents_storage);
134            (hom.map_ref(c), self.create_monomial(exponents_storage.iter().copied()))
135        }));
136    }
137
138    fn clone_monomial(&self, mon: &Self::Monomial) -> Self::Monomial {
139        self.create_monomial((0..self.indeterminate_count()).map(|i| self.exponent_at(mon, i)))
140    }
141
142    /// Returns a list of all variables appearing in the given polynomial.
143    /// Associated with each variable is the highest degree in which it
144    /// appears in some term.
145    ///
146    /// # Example
147    /// ```rust
148    /// # use feanor_math::primitive_int::*;
149    /// # use feanor_math::ring::*;
150    /// # use feanor_math::rings::multivariate::*;
151    /// # use feanor_math::rings::multivariate::multivariate_impl::*;
152    /// let poly_ring = MultivariatePolyRingImpl::new(StaticRing::<i64>::RING, 2);
153    /// let [f, g] =
154    ///     poly_ring.with_wrapped_indeterminates(|[X, Y]| [1 + X + X.pow_ref(2) * Y, X.pow_ref(3)]);
155    /// assert_eq!(vec![(0, 2), (1, 1)], poly_ring.appearing_indeterminates(&f));
156    /// assert_eq!(vec![(0, 3)], poly_ring.appearing_indeterminates(&g));
157    /// ```
158    fn appearing_indeterminates(&self, f: &Self::Element) -> Vec<(usize, usize)> {
159        let mut result = (0..self.indeterminate_count()).map(|_| 0).collect::<Vec<_>>();
160        for (_, m) in self.terms(f) {
161            for (i, result) in result.iter_mut().enumerate() {
162                *result = max(*result, self.exponent_at(m, i));
163            }
164        }
165        return result.into_iter().enumerate().filter(|(_, e)| *e > 0).collect();
166    }
167
168    /// Multiplies two monomials.
169    fn monomial_mul(&self, lhs: Self::Monomial, rhs: &Self::Monomial) -> Self::Monomial {
170        self.create_monomial(
171            (0..self.indeterminate_count()).map(|i| self.exponent_at(&lhs, i) + self.exponent_at(rhs, i)),
172        )
173    }
174
175    /// Returns the degree of a monomial, i.e. the sum of the exponents of all variables.
176    fn monomial_deg(&self, mon: &Self::Monomial) -> usize {
177        (0..self.indeterminate_count()).map(|i| self.exponent_at(mon, i)).sum()
178    }
179
180    /// Returns the least common multiple of two monomials.
181    fn monomial_lcm(&self, lhs: Self::Monomial, rhs: &Self::Monomial) -> Self::Monomial {
182        self.create_monomial(
183            (0..self.indeterminate_count()).map(|i| max(self.exponent_at(&lhs, i), self.exponent_at(rhs, i))),
184        )
185    }
186
187    /// Computes the quotient of two monomials.
188    ///
189    /// If `lhs` does not divide `rhs`, this returns `Result::Err` with the monomial
190    /// `lhs / gcd(rhs, lhs)`.
191    fn monomial_div(&self, lhs: Self::Monomial, rhs: &Self::Monomial) -> Result<Self::Monomial, Self::Monomial> {
192        let mut failed = false;
193        let result = self.create_monomial((0..self.indeterminate_count()).map(|i| {
194            if let Some(res) = self.exponent_at(&lhs, i).checked_sub(self.exponent_at(rhs, i)) {
195                res
196            } else {
197                failed = true;
198                0
199            }
200        }));
201        if failed { Err(result) } else { Ok(result) }
202    }
203
204    /// Evaluates the given polynomial at the given values.
205    ///
206    /// # Example
207    /// ```rust
208    /// # use feanor_math::ring::*;
209    /// # use feanor_math::rings::multivariate::*;
210    /// # use feanor_math::rings::multivariate::multivariate_impl::*;
211    /// # use feanor_math::primitive_int::*;
212    /// # use feanor_math::seq::*;
213    /// let poly_ring = MultivariatePolyRingImpl::new(StaticRing::<i64>::RING, 2);
214    /// let [f] = poly_ring.with_wrapped_indeterminates(|[X, Y]| [1 + X + X.pow_ref(2) * Y]);
215    /// assert_eq!(
216    ///     1 + 5 + 5 * 5 * 8,
217    ///     poly_ring.evaluate(
218    ///         &f,
219    ///         [5, 8].clone_ring_els(StaticRing::<i64>::RING),
220    ///         &poly_ring.base_ring().identity()
221    ///     )
222    /// );
223    /// ```
224    fn evaluate<R, V, H>(&self, f: &Self::Element, value: V, hom: H) -> R::Element
225    where
226        R: ?Sized + RingBase,
227        H: Homomorphism<<Self::BaseRing as RingStore>::Type, R>,
228        V: VectorFn<R::Element>,
229    {
230        assert_eq!(self.indeterminate_count(), value.len());
231        assert!(hom.domain().get_ring() == self.base_ring().get_ring());
232        hom.codomain().sum(self.terms(f).map(|(c, m)| {
233            hom.mul_map(
234                hom.codomain().prod(
235                    (0..self.indeterminate_count()).map(|i| hom.codomain().pow(value.at(i), self.exponent_at(m, i))),
236                ),
237                hom.domain().clone_el(c),
238            )
239        }))
240    }
241
242    /// Replaces the given indeterminate in the given polynomial by the value `val`.
243    ///
244    /// Conceptually, this is similar to [`MultivariatePolyRing::evaluate()`], but less general,
245    /// which can allow a faster implementation sometimes. In particular, this only replaces a
246    /// single indeterminate, and does not change the ring.
247    fn specialize(&self, f: &Self::Element, var: usize, val: &Self::Element) -> Self::Element {
248        assert!(var < self.indeterminate_count());
249        let mut parts = Vec::new();
250        for (c, m) in self.terms(f) {
251            while self.exponent_at(m, var) >= parts.len() {
252                parts.push(Vec::new());
253            }
254            let new_m = self.create_monomial(
255                (0..self.indeterminate_count()).map(|i| if i == var { 0 } else { self.exponent_at(m, i) }),
256            );
257            parts[self.exponent_at(m, var)].push((self.base_ring().clone_el(c), new_m));
258        }
259        if let Some(first) = parts.pop() {
260            let mut current = self.zero();
261            self.add_assign_from_terms(&mut current, first);
262            while let Some(new) = parts.pop() {
263                let mut next = self.zero();
264                self.add_assign_from_terms(&mut next, new);
265                self.mul_assign_ref(&mut current, val);
266                self.add_assign(&mut current, next);
267            }
268            return current;
269        } else {
270            return self.zero();
271        }
272    }
273}
274
275/// [`RingStore`] for [`MultivariatePolyRing`]
276pub trait MultivariatePolyRingStore: RingStore
277where
278    Self::Type: MultivariatePolyRing,
279{
280    delegate! { MultivariatePolyRing, fn indeterminate_count(&self) -> usize }
281    delegate! { MultivariatePolyRing, fn indeterminate(&self, i: usize) -> <Self::Type as MultivariatePolyRing>::Monomial }
282    delegate! { MultivariatePolyRing, fn create_term(&self, coeff: PolyCoeff<Self>, monomial: PolyMonomial<Self>) -> El<Self> }
283    delegate! { MultivariatePolyRing, fn exponent_at(&self, m: &PolyMonomial<Self>, var_index: usize) -> usize }
284    delegate! { MultivariatePolyRing, fn expand_monomial_to(&self, m: &PolyMonomial<Self>, out: &mut [usize]) -> () }
285    delegate! { MultivariatePolyRing, fn monomial_mul(&self, lhs: PolyMonomial<Self>, rhs: &PolyMonomial<Self>) -> PolyMonomial<Self> }
286    delegate! { MultivariatePolyRing, fn monomial_lcm(&self, lhs: PolyMonomial<Self>, rhs: &PolyMonomial<Self>) -> PolyMonomial<Self> }
287    delegate! { MultivariatePolyRing, fn monomial_div(&self, lhs: PolyMonomial<Self>, rhs: &PolyMonomial<Self>) -> Result<PolyMonomial<Self>, PolyMonomial<Self>> }
288    delegate! { MultivariatePolyRing, fn monomial_deg(&self, val: &PolyMonomial<Self>) -> usize }
289    delegate! { MultivariatePolyRing, fn mul_assign_monomial(&self, f: &mut El<Self>, monomial: PolyMonomial<Self>) -> () }
290    delegate! { MultivariatePolyRing, fn appearing_indeterminates(&self, f: &El<Self>) -> Vec<(usize, usize)> }
291    delegate! { MultivariatePolyRing, fn specialize(&self, f: &El<Self>, var: usize, val: &El<Self>) -> El<Self> }
292
293    fn expand_monomial(&self, m: &PolyMonomial<Self>) -> Vec<usize> {
294        let mut result = (0..self.indeterminate_count()).map(|_| 0).collect::<Vec<_>>();
295        self.expand_monomial_to(m, &mut result);
296        return result;
297    }
298
299    /// Returns the term of `f` whose monomial is largest (w.r.t. the given order) among all
300    /// monomials smaller than `lt_than`.
301    fn largest_term_lt<'a, O: MonomialOrder>(
302        &'a self,
303        f: &'a El<Self>,
304        order: O,
305        lt_than: &PolyMonomial<Self>,
306    ) -> Option<(&'a PolyCoeff<Self>, &'a PolyMonomial<Self>)> {
307        self.get_ring().largest_term_lt(f, order, lt_than)
308    }
309
310    /// Returns the **L**eading **T**erm of `f`, i.e. the term whose monomial is largest w.r.t. the
311    /// given order.
312    fn LT<'a, O: MonomialOrder>(
313        &'a self,
314        f: &'a El<Self>,
315        order: O,
316    ) -> Option<(&'a PolyCoeff<Self>, &'a PolyMonomial<Self>)> {
317        self.get_ring().LT(f, order)
318    }
319
320    /// Creates a new monomial with given exponents.
321    ///
322    /// For details, see [`MultivariatePolyRing::create_monomial()`].
323    fn create_monomial<I>(&self, exponents: I) -> PolyMonomial<Self>
324    where
325        I: IntoIterator<Item = usize>,
326        I::IntoIter: ExactSizeIterator,
327    {
328        self.get_ring().create_monomial(exponents)
329    }
330
331    fn clone_monomial(&self, mon: &PolyMonomial<Self>) -> PolyMonomial<Self> { self.get_ring().clone_monomial(mon) }
332
333    /// Returns the coefficient of a polynomial corresponding to a monomial.
334    ///
335    /// For details, see [`MultivariatePolyRing::coefficient_at()`].
336    fn coefficient_at<'a>(&'a self, f: &'a El<Self>, m: &PolyMonomial<Self>) -> &'a PolyCoeff<Self> {
337        self.get_ring().coefficient_at(f, m)
338    }
339
340    /// Returns an iterator over all nonzero terms of the polynomial.
341    ///
342    /// For details, see [`MultivariatePolyRing::terms()`].
343    fn terms<'a>(&'a self, f: &'a El<Self>) -> <Self::Type as MultivariatePolyRing>::TermIter<'a> {
344        self.get_ring().terms(f)
345    }
346
347    /// Creates a new polynomial by summing up all the given terms.
348    fn from_terms<I>(&self, terms: I) -> El<Self>
349    where
350        I: IntoIterator<Item = (PolyCoeff<Self>, PolyMonomial<Self>)>,
351    {
352        let mut result = self.zero();
353        self.get_ring().add_assign_from_terms(&mut result, terms);
354        return result;
355    }
356
357    /// Evaluates the polynomial at the given values.
358    ///
359    /// For details, see [`MultivariatePolyRing::evaluate()`].
360    fn evaluate<R, V, H>(&self, f: &El<Self>, value: V, hom: H) -> R::Element
361    where
362        R: ?Sized + RingBase,
363        H: Homomorphism<<<Self::Type as RingExtension>::BaseRing as RingStore>::Type, R>,
364        V: VectorFn<R::Element>,
365    {
366        self.get_ring().evaluate(f, value, hom)
367    }
368
369    /// Returns the homomorphism `R[X1, ..., Xm] -> S[X1, ..., Xm]` that is induced by
370    /// applying the given homomorphism `R -> S` coefficient-wise.
371    fn into_lifted_hom<P, H>(self, from: P, hom: H) -> CoefficientHom<P, Self, H>
372    where
373        P: RingStore,
374        P::Type: MultivariatePolyRing,
375        H: Homomorphism<
376                <<P::Type as RingExtension>::BaseRing as RingStore>::Type,
377                <<Self::Type as RingExtension>::BaseRing as RingStore>::Type,
378            >,
379    {
380        CoefficientHom { from, to: self, hom }
381    }
382
383    /// Returns the homomorphism `R[X1, ..., Xm] -> S[X1, ..., Xm]` that is induced by
384    /// applying the given homomorphism `R -> S` coefficient-wise.
385    ///
386    /// If the ownership of this ring should be transferred to the homomorphism, consider
387    /// using [`MultivariatePolyRingStore::into_lifted_hom()`].
388    fn lifted_hom<P, H>(&self, from: P, hom: H) -> CoefficientHom<P, &Self, H>
389    where
390        P: RingStore,
391        P::Type: MultivariatePolyRing,
392        H: Homomorphism<
393                <<P::Type as RingExtension>::BaseRing as RingStore>::Type,
394                <<Self::Type as RingExtension>::BaseRing as RingStore>::Type,
395            >,
396    {
397        self.into_lifted_hom(from, hom)
398    }
399
400    /// Invokes the function with a wrapped version of the indeterminates of this poly ring.
401    /// Use for convenient creation of polynomials.
402    ///
403    /// Note however that [`MultivariatePolyRingStore::from_terms()`] might be more performant.
404    ///
405    /// # Example
406    /// ```rust
407    /// use feanor_math::assert_el_eq;
408    /// use feanor_math::homomorphism::*;
409    /// use feanor_math::ring::*;
410    /// use feanor_math::rings::multivariate::multivariate_impl::*;
411    /// use feanor_math::rings::multivariate::*;
412    /// use feanor_math::rings::zn::zn_64::*;
413    /// let base_ring = Zn::new(7);
414    /// let poly_ring = MultivariatePolyRingImpl::new(base_ring, 3);
415    /// let f_version1 = poly_ring.from_terms(
416    ///     [
417    ///         (
418    ///             base_ring.int_hom().map(3),
419    ///             poly_ring.create_monomial([0, 0, 0]),
420    ///         ),
421    ///         (
422    ///             base_ring.int_hom().map(2),
423    ///             poly_ring.create_monomial([0, 1, 1]),
424    ///         ),
425    ///         (base_ring.one(), poly_ring.create_monomial([2, 0, 0])),
426    ///     ]
427    ///     .into_iter(),
428    /// );
429    /// let f_version2 = poly_ring
430    ///     .with_wrapped_indeterminates_dyn(|[x, y, z]| [3 + 2 * y * z + x.pow_ref(2)])
431    ///     .into_iter()
432    ///     .next()
433    ///     .unwrap();
434    /// let [f_version3] =
435    ///     poly_ring.with_wrapped_indeterminates(|[x, y, z]| [3 + 2 * y * z + x.pow_ref(2)]);
436    /// assert_el_eq!(poly_ring, f_version1, f_version2);
437    /// ```
438    #[stability::unstable(feature = "enable")]
439    fn with_wrapped_indeterminates_dyn<'a, F, T, const N: usize>(&'a self, f: F) -> Vec<El<Self>>
440    where
441        F: FnOnce([&RingElementWrapper<&'a Self>; N]) -> T,
442        T: IntoIterator<Item = RingElementWrapper<&'a Self>>,
443    {
444        assert_eq!(self.indeterminate_count(), N);
445        let wrapped_indets: [_; N] = std::array::from_fn(|i| {
446            RingElementWrapper::new(self, self.create_term(self.base_ring().one(), self.indeterminate(i)))
447        });
448        f(std::array::from_fn(|i| &wrapped_indets[i]))
449            .into_iter()
450            .map(|f| f.unwrap())
451            .collect()
452    }
453
454    /// Same as [`MultivariatePolyRingStore::with_wrapped_indeterminates_dyn()`], but returns result
455    /// as an array to allow pattern matching.
456    ///
457    /// # Example
458    /// ```rust
459    /// use feanor_math::assert_el_eq;
460    /// use feanor_math::homomorphism::*;
461    /// use feanor_math::ring::*;
462    /// use feanor_math::rings::multivariate::multivariate_impl::*;
463    /// use feanor_math::rings::multivariate::*;
464    /// use feanor_math::rings::zn::zn_64::*;
465    /// let base_ring = Zn::new(7);
466    /// let poly_ring = MultivariatePolyRingImpl::new(base_ring, 3);
467    /// let f_version1 = poly_ring.from_terms(
468    ///     [
469    ///         (
470    ///             base_ring.int_hom().map(3),
471    ///             poly_ring.create_monomial([0, 0, 0]),
472    ///         ),
473    ///         (
474    ///             base_ring.int_hom().map(2),
475    ///             poly_ring.create_monomial([0, 1, 1]),
476    ///         ),
477    ///         (base_ring.one(), poly_ring.create_monomial([2, 0, 0])),
478    ///     ]
479    ///     .into_iter(),
480    /// );
481    /// let f_version2 = poly_ring
482    ///     .with_wrapped_indeterminates_dyn(|[x, y, z]| [3 + 2 * y * z + x.pow_ref(2)])
483    ///     .into_iter()
484    ///     .next()
485    ///     .unwrap();
486    /// let [f_version3] =
487    ///     poly_ring.with_wrapped_indeterminates(|[x, y, z]| [3 + 2 * y * z + x.pow_ref(2)]);
488    /// assert_el_eq!(poly_ring, f_version1, f_version2);
489    /// ```
490    #[stability::unstable(feature = "enable")]
491    fn with_wrapped_indeterminates<'a, F, const N: usize, const M: usize>(&'a self, f: F) -> [El<Self>; M]
492    where
493        F: FnOnce([&RingElementWrapper<&'a Self>; N]) -> [RingElementWrapper<&'a Self>; M],
494    {
495        assert_eq!(self.indeterminate_count(), N);
496        let wrapped_indets: [_; N] = std::array::from_fn(|i| {
497            RingElementWrapper::new(self, self.create_term(self.base_ring().one(), self.indeterminate(i)))
498        });
499        let mut result_it = f(std::array::from_fn(|i| &wrapped_indets[i]))
500            .into_iter()
501            .map(|f| f.unwrap());
502        let result = std::array::from_fn(|_| result_it.next().unwrap());
503        debug_assert!(result_it.next().is_none());
504        return result;
505    }
506}
507
508impl<P> MultivariatePolyRingStore for P
509where
510    P: RingStore,
511    P::Type: MultivariatePolyRing,
512{
513}
514
515/// Trait for monomial orders, i.e. orderings on the monomials `X1^i1 ... Xm^im` of a multivariate
516/// polynomial ring that are compatible with multiplication. These are also sometimes called term
517/// orders.
518///
519/// To be precise, a total order on monomials is a monomial order, if for all `u, v, w` with `u <=
520/// v` also `uw <= vw`.
521pub trait MonomialOrder: Clone {
522    fn compare<P>(&self, ring: P, lhs: &PolyMonomial<P>, rhs: &PolyMonomial<P>) -> Ordering
523    where
524        P: RingStore,
525        P::Type: MultivariatePolyRing;
526
527    /// Checks whether two monomials are equal.
528    ///
529    /// This may be faster than [`MonomialOrder::compare()`], but clearly must have
530    /// a compatible behavior.
531    fn eq_mon<P>(&self, ring: P, lhs: &PolyMonomial<P>, rhs: &PolyMonomial<P>) -> bool
532    where
533        P: RingStore,
534        P::Type: MultivariatePolyRing,
535    {
536        self.compare(ring, lhs, rhs) == Ordering::Equal
537    }
538
539    /// Whether this order is the same as the given other order, i.e. [`MonomialOrder::compare()`]
540    /// gives the same output on all inputs.
541    ///
542    /// Many monomial orders are likely to be implemented as zero-sized types with only a single
543    /// instance. In this case, the default implementation is sufficient.
544    fn is_same<O>(&self, rhs: &O) -> bool
545    where
546        O: MonomialOrder,
547    {
548        assert!(self.as_any().is_some());
549        assert!(std::mem::size_of::<Self>() == 0);
550        if let Some(rhs_as_any) = rhs.as_any() {
551            self.as_any().unwrap().type_id() == rhs_as_any.type_id()
552        } else {
553            false
554        }
555    }
556
557    /// Upcasts this reference to `&dyn Any`, which is sometimes required to compare monomial order
558    /// objects of different types.
559    fn as_any(&self) -> Option<&dyn Any>;
560}
561
562/// Trait for [`MonomialOrder`]s that are graded, i.e. for `v, w` with `deg(v) < deg(w)`
563/// they always satisfy `v < w`.
564pub trait GradedMonomialOrder: MonomialOrder {}
565
566/// The graded reverse lexicographic order. This is the most important monomial order, since
567/// computing a Groebner basis w.r.t. this order is usually more efficient than for other orders.
568/// Also sometimes referred to as "grevlex".
569///
570/// Concretely, this is the ordering of monomials we get by first comparing monomial degrees, and
571/// in case of a tie reverse the outcome of a lexicographic comparison, using a reversed order
572/// of variables.
573///
574/// # Example
575/// ```rust
576/// # use feanor_math::ring::*;
577/// # use feanor_math::rings::multivariate::*;
578/// # use feanor_math::rings::multivariate::multivariate_impl::*;
579/// # use feanor_math::primitive_int::*;
580/// # use std::cmp::Ordering;
581/// let poly_ring = MultivariatePolyRingImpl::new(StaticRing::<i64>::RING, 3);
582/// let monomials_descending = [
583///     [2, 0, 0], // x1^2
584///     [1, 1, 0], // x1 x2
585///     [0, 2, 0], // x2^2
586///     [1, 0, 1], // x1 x3
587///     [0, 1, 1], // x2 x3
588///     [0, 0, 2], // x3^2
589/// ]
590/// .into_iter()
591/// .map(|m| poly_ring.create_monomial(m))
592/// .collect::<Vec<_>>();
593/// for i in 1..6 {
594///     assert!(
595///         DegRevLex.compare(
596///             &poly_ring,
597///             &monomials_descending[i - 1],
598///             &monomials_descending[i]
599///         ) == Ordering::Greater
600///     );
601/// }
602/// ```
603#[derive(Clone, Copy, PartialEq, Eq, Debug)]
604pub struct DegRevLex;
605
606impl MonomialOrder for DegRevLex {
607    fn as_any(&self) -> Option<&dyn Any> { Some(self as &dyn Any) }
608
609    fn compare<P>(&self, ring: P, lhs: &PolyMonomial<P>, rhs: &PolyMonomial<P>) -> Ordering
610    where
611        P: RingStore,
612        P::Type: MultivariatePolyRing,
613    {
614        let lhs_deg = ring.monomial_deg(lhs);
615        let rhs_deg = ring.monomial_deg(rhs);
616        if lhs_deg < rhs_deg {
617            return Ordering::Less;
618        } else if lhs_deg > rhs_deg {
619            return Ordering::Greater;
620        } else {
621            for i in (0..ring.indeterminate_count()).rev() {
622                if ring.exponent_at(lhs, i) > ring.exponent_at(rhs, i) {
623                    return Ordering::Less;
624                } else if ring.exponent_at(lhs, i) < ring.exponent_at(rhs, i) {
625                    return Ordering::Greater;
626                }
627            }
628            return Ordering::Equal;
629        }
630    }
631}
632
633impl GradedMonomialOrder for DegRevLex {}
634
635/// Lexicographic ordering of monomials.
636#[derive(Clone, Copy, PartialEq, Eq, Debug)]
637pub struct Lex;
638
639impl MonomialOrder for Lex {
640    fn as_any(&self) -> Option<&dyn Any> { Some(self as &dyn Any) }
641
642    fn compare<P>(&self, ring: P, lhs: &PolyMonomial<P>, rhs: &PolyMonomial<P>) -> Ordering
643    where
644        P: RingStore,
645        P::Type: MultivariatePolyRing,
646    {
647        for i in 0..ring.indeterminate_count() {
648            match ring.exponent_at(lhs, i).cmp(&ring.exponent_at(rhs, i)) {
649                Ordering::Less => {
650                    return Ordering::Less;
651                }
652                Ordering::Greater => {
653                    return Ordering::Greater;
654                }
655                Ordering::Equal => {}
656            }
657        }
658        return Ordering::Equal;
659    }
660}
661
662pub struct CoefficientHom<PFrom, PTo, H>
663where
664    PFrom: RingStore,
665    PTo: RingStore,
666    PFrom::Type: MultivariatePolyRing,
667    PTo::Type: MultivariatePolyRing,
668    H: Homomorphism<
669            <<PFrom::Type as RingExtension>::BaseRing as RingStore>::Type,
670            <<PTo::Type as RingExtension>::BaseRing as RingStore>::Type,
671        >,
672{
673    from: PFrom,
674    to: PTo,
675    hom: H,
676}
677
678impl<PFrom, PTo, H> Homomorphism<PFrom::Type, PTo::Type> for CoefficientHom<PFrom, PTo, H>
679where
680    PFrom: RingStore,
681    PTo: RingStore,
682    PFrom::Type: MultivariatePolyRing,
683    PTo::Type: MultivariatePolyRing,
684    H: Homomorphism<
685            <<PFrom::Type as RingExtension>::BaseRing as RingStore>::Type,
686            <<PTo::Type as RingExtension>::BaseRing as RingStore>::Type,
687        >,
688{
689    type DomainStore = PFrom;
690    type CodomainStore = PTo;
691
692    fn codomain(&self) -> &Self::CodomainStore { &self.to }
693
694    fn domain(&self) -> &Self::DomainStore { &self.from }
695
696    fn map(&self, x: <PFrom::Type as RingBase>::Element) -> <PTo::Type as RingBase>::Element { self.map_ref(&x) }
697
698    fn map_ref(&self, x: &<PFrom::Type as RingBase>::Element) -> <PTo::Type as RingBase>::Element {
699        self.to.get_ring().map_terms(self.from.get_ring(), x, &self.hom)
700    }
701}
702
703pub mod generic_impls {
704    use std::fmt::{Formatter, Result};
705
706    use super::*;
707
708    #[stability::unstable(feature = "enable")]
709    pub fn print<P>(ring: P, poly: &El<P>, out: &mut Formatter, env: EnvBindingStrength) -> Result
710    where
711        P: RingStore,
712        P::Type: MultivariatePolyRing,
713    {
714        if ring.is_zero(poly) {
715            ring.base_ring()
716                .get_ring()
717                .dbg_within(&ring.base_ring().zero(), out, env)?;
718        } else {
719            if env >= EnvBindingStrength::Product {
720                write!(out, "(")?;
721            }
722
723            let mut print_term = |c: &PolyCoeff<P>, m: &PolyMonomial<P>, with_plus: bool| {
724                if with_plus {
725                    write!(out, " + ")?;
726                }
727                let is_constant_term = ring.monomial_deg(m) == 0;
728                if !ring.base_ring().is_one(c) || is_constant_term {
729                    ring.base_ring().get_ring().dbg_within(
730                        c,
731                        out,
732                        if is_constant_term {
733                            EnvBindingStrength::Sum
734                        } else {
735                            EnvBindingStrength::Product
736                        },
737                    )?;
738                    if !is_constant_term {
739                        write!(out, " * ")?;
740                    }
741                }
742                let mut needs_space = false;
743                for i in 0..ring.indeterminate_count() {
744                    if ring.exponent_at(m, i) > 0 {
745                        if needs_space {
746                            write!(out, " * ")?;
747                        }
748                        write!(out, "X{}", i)?;
749                        needs_space = true;
750                    }
751                    if ring.exponent_at(m, i) > 1 {
752                        write!(out, "^{}", ring.exponent_at(m, i))?;
753                    }
754                }
755                return Ok::<(), std::fmt::Error>(());
756            };
757
758            for (i, (c, m)) in ring.terms(poly).enumerate() {
759                print_term(c, m, i != 0)?;
760            }
761            if env >= EnvBindingStrength::Product {
762                write!(out, ")")?;
763            }
764        }
765
766        return Ok(());
767    }
768}
769
770#[cfg(any(test, feature = "generic_tests"))]
771pub mod generic_tests {
772    use super::*;
773    use crate::seq::*;
774
775    #[stability::unstable(feature = "enable")]
776    pub fn test_poly_ring_axioms<P: RingStore, I: Iterator<Item = PolyCoeff<P>>>(
777        ring: P,
778        interesting_base_ring_elements: I,
779    ) where
780        P::Type: MultivariatePolyRing,
781    {
782        let elements = interesting_base_ring_elements.collect::<Vec<_>>();
783        let n = ring.indeterminate_count();
784        let base_ring = ring.base_ring();
785
786        // test multiplication of variables
787        for i in 0..n {
788            for j in 0..n {
789                let xi = ring.create_term(base_ring.one(), ring.indeterminate(i));
790                let xj = ring.create_term(base_ring.one(), ring.indeterminate(j));
791                let xixj = ring.create_term(
792                    base_ring.one(),
793                    ring.create_monomial((0..n).map(|k| {
794                        if k == i && k == j {
795                            2
796                        } else if k == j || k == i {
797                            1
798                        } else {
799                            0
800                        }
801                    })),
802                );
803                assert_el_eq!(ring, xixj, ring.mul(xi, xj));
804            }
805        }
806
807        // test monomial operations
808        for i in 0..n {
809            for j in 0..n {
810                let xi = ring.indeterminate(i);
811                let xj = ring.indeterminate(j);
812                let xixj_lcm = ring.create_monomial((0..n).map(|k| if k == j || k == i { 1 } else { 0 }));
813                assert_el_eq!(
814                    ring,
815                    ring.create_term(base_ring.one(), xixj_lcm),
816                    ring.create_term(base_ring.one(), ring.monomial_lcm(xi, &xj))
817                );
818            }
819        }
820
821        // all monomials should be different
822        for i in 0..n {
823            for j in 0..n {
824                let xi = ring.create_term(
825                    base_ring.one(),
826                    ring.create_monomial((0..n).map(|k| if k == i { 1 } else { 0 })),
827                );
828                let xj = ring.create_term(
829                    base_ring.one(),
830                    ring.create_monomial((0..n).map(|k| if k == j { 1 } else { 0 })),
831                );
832                assert!((i == j) == ring.eq_el(&xi, &xj));
833            }
834        }
835
836        // monomials shouldn't be zero divisors
837        for i in 0..n {
838            for a in &elements {
839                let xi = ring.create_term(base_ring.one(), ring.indeterminate(i));
840                assert!(base_ring.is_zero(a) == ring.is_zero(&ring.inclusion().mul_ref_map(&xi, a)));
841            }
842        }
843
844        // test add_assign_from_terms
845        for i in 0..n {
846            let xi = ring.indeterminate(i);
847            let mut a = ring.create_term(base_ring.int_hom().map(3), ring.create_monomial((0..n).map(|_| 0)));
848            let terms_with_multiples = [
849                (base_ring.one(), ring.clone_monomial(&xi)),
850                (base_ring.one(), ring.clone_monomial(&xi)),
851                (base_ring.one(), ring.create_monomial((0..n).map(|_| 0))),
852                (base_ring.one(), ring.create_monomial((0..n).map(|_| 0))),
853                (base_ring.one(), ring.clone_monomial(&xi)),
854                (base_ring.one(), ring.create_monomial((0..n).map(|_| 0))),
855                (base_ring.one(), ring.clone_monomial(&xi)),
856                (base_ring.one(), ring.create_monomial((0..n).map(|_| 0))),
857            ];
858            ring.get_ring().add_assign_from_terms(&mut a, terms_with_multiples);
859            assert_el_eq!(
860                &ring,
861                ring.from_terms([
862                    (base_ring.int_hom().map(7), ring.create_monomial((0..n).map(|_| 0))),
863                    (base_ring.int_hom().map(4), xi),
864                ]),
865                a
866            );
867        }
868
869        if n >= 2 {
870            let one = ring.create_monomial((0..n).map(|_| 0));
871            let x0 = ring.indeterminate(0);
872            let x1 = ring.indeterminate(1);
873            let x0_v = ring.create_term(base_ring.one(), ring.clone_monomial(&x0));
874            let x1_v = ring.create_term(base_ring.one(), ring.clone_monomial(&x1));
875            let x0_2 = ring.create_monomial((0..n).map(|k| if k == 0 { 2 } else { 0 }));
876            let x0_3 = ring.create_monomial((0..n).map(|k| if k == 0 { 3 } else { 0 }));
877            let x0_4 = ring.create_monomial((0..n).map(|k| if k == 0 { 4 } else { 0 }));
878            let x0x1 = ring.create_monomial((0..n).map(|k| if k == 0 || k == 1 { 1 } else { 0 }));
879            let x0_3x1 = ring.create_monomial((0..n).map(|k| {
880                if k == 0 {
881                    3
882                } else if k == 1 {
883                    1
884                } else {
885                    0
886                }
887            }));
888            let x1_2 = ring.create_monomial((0..n).map(|k| if k == 1 { 2 } else { 0 }));
889
890            // test product
891            for a in &elements {
892                for b in &elements {
893                    for c in &elements {
894                        let f = ring.add(
895                            ring.inclusion().mul_ref_map(&x0_v, a),
896                            ring.inclusion().mul_ref_map(&x1_v, b),
897                        );
898                        let g = ring.add(ring.inclusion().mul_ref_map(&x0_v, c), ring.clone_el(&x1_v));
899                        let h = ring.from_terms(
900                            [
901                                (base_ring.mul_ref(a, c), ring.clone_monomial(&x0_2)),
902                                (
903                                    base_ring.add_ref_snd(base_ring.mul_ref(b, c), a),
904                                    ring.clone_monomial(&x0x1),
905                                ),
906                                (base_ring.clone_el(b), ring.clone_monomial(&x1_2)),
907                            ]
908                            .into_iter(),
909                        );
910                        assert_el_eq!(ring, h, ring.mul(f, g));
911                    }
912                }
913            }
914
915            // test sum
916            for a in &elements {
917                for b in &elements {
918                    for c in &elements {
919                        let f = ring.from_terms([
920                            (base_ring.clone_el(a), ring.clone_monomial(&one)),
921                            (base_ring.clone_el(c), ring.clone_monomial(&x0_2)),
922                            (base_ring.one(), ring.clone_monomial(&x0x1)),
923                        ]);
924                        let g = ring.from_terms([
925                            (base_ring.clone_el(b), ring.clone_monomial(&x0)),
926                            (base_ring.one(), ring.clone_monomial(&x0_2)),
927                        ]);
928                        let h = ring.from_terms([
929                            (base_ring.clone_el(a), ring.clone_monomial(&one)),
930                            (base_ring.clone_el(b), ring.clone_monomial(&x0)),
931                            (base_ring.add_ref_fst(c, base_ring.one()), ring.clone_monomial(&x0_2)),
932                            (base_ring.one(), ring.clone_monomial(&x0x1)),
933                        ]);
934                        assert_el_eq!(ring, h, ring.add(f, g));
935                    }
936                }
937            }
938
939            // test mul_assign_monomial
940            for a in &elements {
941                for b in &elements {
942                    for c in &elements {
943                        let mut f = ring.from_terms([
944                            (base_ring.clone_el(a), ring.clone_monomial(&one)),
945                            (base_ring.clone_el(b), ring.clone_monomial(&x0)),
946                            (base_ring.clone_el(c), ring.clone_monomial(&x0_2)),
947                            (base_ring.one(), ring.clone_monomial(&x0x1)),
948                        ]);
949                        let h = ring.from_terms([
950                            (base_ring.clone_el(a), ring.clone_monomial(&x0_2)),
951                            (base_ring.clone_el(b), ring.clone_monomial(&x0_3)),
952                            (base_ring.clone_el(c), ring.clone_monomial(&x0_4)),
953                            (base_ring.one(), ring.clone_monomial(&x0_3x1)),
954                        ]);
955                        let m = ring.clone_monomial(&x0_2);
956                        ring.mul_assign_monomial(&mut f, m);
957                        assert_el_eq!(ring, h, f);
958                    }
959                }
960            }
961
962            // test evaluate
963            for a in &elements {
964                for b in &elements {
965                    let f = ring.from_terms([
966                        (base_ring.int_hom().map(3), ring.clone_monomial(&one)),
967                        (base_ring.int_hom().map(10), ring.clone_monomial(&x0)),
968                        (base_ring.neg_one(), ring.clone_monomial(&x0_2)),
969                        (base_ring.one(), ring.clone_monomial(&x0x1)),
970                    ]);
971                    let expected = <_ as RingStore>::sum(
972                        base_ring,
973                        [
974                            base_ring.int_hom().map(3),
975                            base_ring.int_hom().mul_ref_map(a, &10),
976                            base_ring.negate(base_ring.pow(base_ring.clone_el(a), 2)),
977                            base_ring.mul_ref(a, b),
978                        ],
979                    );
980                    let values = [base_ring.clone_el(a), base_ring.clone_el(b)]
981                        .into_iter()
982                        .chain((0..(ring.indeterminate_count() - 2)).map(|_| base_ring.zero()))
983                        .collect::<Vec<_>>();
984                    assert_el_eq!(
985                        &base_ring,
986                        &expected,
987                        &ring.evaluate(&f, values.into_clone_ring_els(base_ring), &base_ring.identity())
988                    );
989                }
990            }
991        }
992    }
993}