Skip to main content

oxiz_math/polynomial/
types.rs

1//! Type definitions for polynomial arithmetic.
2//!
3//! This module contains the core types used in polynomial representation:
4//! - VarPower: Variable with exponent
5//! - Monomial: Product of variables with exponents
6//! - Term: Coefficient multiplied by a monomial
7//! - MonomialOrder: Ordering for polynomial canonicalization
8
9#[allow(unused_imports)]
10use crate::prelude::*;
11use core::cmp::Ordering;
12use core::fmt;
13use core::hash::{Hash, Hasher};
14use num_rational::BigRational;
15use num_traits::{One, Zero};
16use smallvec::SmallVec;
17
18/// Variable identifier for polynomials.
19pub type Var = u32;
20
21/// Null variable constant (indicates no variable).
22pub const NULL_VAR: Var = u32::MAX;
23
24/// Power of a variable (variable, exponent).
25#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
26pub struct VarPower {
27    /// The variable identifier.
28    pub var: Var,
29    /// The exponent (power) of the variable.
30    pub power: u32,
31}
32
33impl VarPower {
34    /// Create a new variable power.
35    #[inline]
36    pub fn new(var: Var, power: u32) -> Self {
37        Self { var, power }
38    }
39}
40
41impl PartialOrd for VarPower {
42    fn partial_cmp(&self, other: &Self) -> Option<Ordering> {
43        Some(self.cmp(other))
44    }
45}
46
47impl Ord for VarPower {
48    fn cmp(&self, other: &Self) -> Ordering {
49        self.var.cmp(&other.var)
50    }
51}
52
53/// A monomial is a product of variables with exponents.
54/// Represented as a sorted list of (variable, power) pairs.
55/// The unit monomial (1) is represented as an empty list.
56#[derive(Clone, PartialEq, Eq)]
57pub struct Monomial {
58    /// Variables with their exponents, sorted by variable index.
59    pub(crate) vars: SmallVec<[VarPower; 4]>,
60    /// Cached total degree.
61    pub(crate) total_degree: u32,
62    /// Cached hash value.
63    pub(crate) hash: u64,
64}
65
66impl Monomial {
67    /// Create the unit monomial (1).
68    #[inline]
69    pub fn unit() -> Self {
70        Self {
71            vars: SmallVec::new(),
72            total_degree: 0,
73            hash: 0,
74        }
75    }
76
77    /// Create a monomial from a single variable with power 1.
78    #[inline]
79    pub fn from_var(var: Var) -> Self {
80        Self::from_var_power(var, 1)
81    }
82
83    /// Create a monomial from a single variable with a given power.
84    pub fn from_var_power(var: Var, power: u32) -> Self {
85        if power == 0 {
86            return Self::unit();
87        }
88        let mut vars = SmallVec::new();
89        vars.push(VarPower::new(var, power));
90        Self {
91            total_degree: power,
92            hash: compute_monomial_hash(&vars),
93            vars,
94        }
95    }
96
97    /// Create a monomial from a list of (variable, power) pairs.
98    /// The input doesn't need to be sorted or normalized.
99    pub fn from_powers(powers: impl IntoIterator<Item = (Var, u32)>) -> Self {
100        let mut var_powers: FxHashMap<Var, u32> = FxHashMap::default();
101        for (var, power) in powers {
102            if power > 0 {
103                *var_powers.entry(var).or_insert(0) += power;
104            }
105        }
106
107        let mut vars: SmallVec<[VarPower; 4]> = var_powers
108            .into_iter()
109            .filter(|(_, p)| *p > 0)
110            .map(|(v, p)| VarPower::new(v, p))
111            .collect();
112        vars.sort_by_key(|vp| vp.var);
113
114        let total_degree = vars.iter().map(|vp| vp.power).sum();
115        let hash = compute_monomial_hash(&vars);
116
117        Self {
118            vars,
119            total_degree,
120            hash,
121        }
122    }
123
124    /// Returns true if this is the unit monomial.
125    #[inline]
126    pub fn is_unit(&self) -> bool {
127        self.vars.is_empty()
128    }
129
130    /// Returns the total degree of the monomial.
131    #[inline]
132    pub fn total_degree(&self) -> u32 {
133        self.total_degree
134    }
135
136    /// Returns the number of variables in the monomial.
137    #[inline]
138    pub fn num_vars(&self) -> usize {
139        self.vars.len()
140    }
141
142    /// Returns the variable-power pairs.
143    #[inline]
144    pub fn vars(&self) -> &[VarPower] {
145        &self.vars
146    }
147
148    /// Returns the degree of a specific variable in this monomial.
149    pub fn degree(&self, var: Var) -> u32 {
150        self.vars
151            .iter()
152            .find(|vp| vp.var == var)
153            .map(|vp| vp.power)
154            .unwrap_or(0)
155    }
156
157    /// Returns the maximum variable in this monomial, or NULL_VAR if unit.
158    pub fn max_var(&self) -> Var {
159        self.vars.last().map(|vp| vp.var).unwrap_or(NULL_VAR)
160    }
161
162    /// Check if this monomial is univariate (contains at most one variable).
163    #[inline]
164    pub fn is_univariate(&self) -> bool {
165        self.vars.len() <= 1
166    }
167
168    /// Check if this monomial is linear (degree 0 or 1).
169    #[inline]
170    pub fn is_linear(&self) -> bool {
171        self.total_degree <= 1
172    }
173
174    /// Multiply two monomials.
175    pub fn mul(&self, other: &Monomial) -> Monomial {
176        if self.is_unit() {
177            return other.clone();
178        }
179        if other.is_unit() {
180            return self.clone();
181        }
182
183        let mut vars: SmallVec<[VarPower; 4]> = SmallVec::new();
184        let mut i = 0;
185        let mut j = 0;
186
187        while i < self.vars.len() && j < other.vars.len() {
188            match self.vars[i].var.cmp(&other.vars[j].var) {
189                Ordering::Less => {
190                    vars.push(self.vars[i]);
191                    i += 1;
192                }
193                Ordering::Greater => {
194                    vars.push(other.vars[j]);
195                    j += 1;
196                }
197                Ordering::Equal => {
198                    vars.push(VarPower::new(
199                        self.vars[i].var,
200                        self.vars[i].power + other.vars[j].power,
201                    ));
202                    i += 1;
203                    j += 1;
204                }
205            }
206        }
207        vars.extend_from_slice(&self.vars[i..]);
208        vars.extend_from_slice(&other.vars[j..]);
209
210        let total_degree = self.total_degree + other.total_degree;
211        let hash = compute_monomial_hash(&vars);
212
213        Monomial {
214            vars,
215            total_degree,
216            hash,
217        }
218    }
219
220    /// Check if other divides self. Returns the quotient if it does.
221    pub fn div(&self, other: &Monomial) -> Option<Monomial> {
222        if other.is_unit() {
223            return Some(self.clone());
224        }
225
226        let mut vars: SmallVec<[VarPower; 4]> = SmallVec::new();
227        let mut j = 0;
228
229        for vp in &self.vars {
230            if j < other.vars.len() && other.vars[j].var == vp.var {
231                if vp.power < other.vars[j].power {
232                    return None;
233                }
234                let new_power = vp.power - other.vars[j].power;
235                if new_power > 0 {
236                    vars.push(VarPower::new(vp.var, new_power));
237                }
238                j += 1;
239            } else if j < other.vars.len() && other.vars[j].var < vp.var {
240                return None;
241            } else {
242                vars.push(*vp);
243            }
244        }
245
246        if j < other.vars.len() {
247            return None;
248        }
249
250        let total_degree = vars.iter().map(|vp| vp.power).sum();
251        let hash = compute_monomial_hash(&vars);
252
253        Some(Monomial {
254            vars,
255            total_degree,
256            hash,
257        })
258    }
259
260    /// Compute the GCD of two monomials.
261    pub fn gcd(&self, other: &Monomial) -> Monomial {
262        if self.is_unit() || other.is_unit() {
263            return Monomial::unit();
264        }
265
266        let mut vars: SmallVec<[VarPower; 4]> = SmallVec::new();
267        let mut i = 0;
268        let mut j = 0;
269
270        while i < self.vars.len() && j < other.vars.len() {
271            match self.vars[i].var.cmp(&other.vars[j].var) {
272                Ordering::Less => {
273                    i += 1;
274                }
275                Ordering::Greater => {
276                    j += 1;
277                }
278                Ordering::Equal => {
279                    let min_power = self.vars[i].power.min(other.vars[j].power);
280                    if min_power > 0 {
281                        vars.push(VarPower::new(self.vars[i].var, min_power));
282                    }
283                    i += 1;
284                    j += 1;
285                }
286            }
287        }
288
289        let total_degree = vars.iter().map(|vp| vp.power).sum();
290        let hash = compute_monomial_hash(&vars);
291
292        Monomial {
293            vars,
294            total_degree,
295            hash,
296        }
297    }
298
299    /// Raise monomial to a power.
300    pub fn pow(&self, n: u32) -> Monomial {
301        if n == 0 {
302            return Monomial::unit();
303        }
304        if n == 1 {
305            return self.clone();
306        }
307
308        let vars: SmallVec<[VarPower; 4]> = self
309            .vars
310            .iter()
311            .map(|vp| VarPower::new(vp.var, vp.power * n))
312            .collect();
313        let total_degree = self.total_degree * n;
314        let hash = compute_monomial_hash(&vars);
315
316        Monomial {
317            vars,
318            total_degree,
319            hash,
320        }
321    }
322
323    /// Lexicographic comparison of monomials.
324    pub fn lex_cmp(&self, other: &Monomial) -> Ordering {
325        let mut i = 0;
326        let mut j = 0;
327
328        while i < self.vars.len() && j < other.vars.len() {
329            match self.vars[i].var.cmp(&other.vars[j].var) {
330                Ordering::Less => return Ordering::Greater,
331                Ordering::Greater => return Ordering::Less,
332                Ordering::Equal => match self.vars[i].power.cmp(&other.vars[j].power) {
333                    Ordering::Equal => {
334                        i += 1;
335                        j += 1;
336                    }
337                    ord => return ord,
338                },
339            }
340        }
341
342        if i < self.vars.len() {
343            Ordering::Greater
344        } else if j < other.vars.len() {
345            Ordering::Less
346        } else {
347            Ordering::Equal
348        }
349    }
350
351    /// Graded lexicographic comparison (total degree first, then lex).
352    pub fn grlex_cmp(&self, other: &Monomial) -> Ordering {
353        match self.total_degree.cmp(&other.total_degree) {
354            Ordering::Equal => self.lex_cmp(other),
355            ord => ord,
356        }
357    }
358
359    /// Graded reverse lexicographic comparison.
360    pub fn grevlex_cmp(&self, other: &Monomial) -> Ordering {
361        match self.total_degree.cmp(&other.total_degree) {
362            Ordering::Equal => {
363                // Reverse lex: compare from highest variable
364                let mut i = self.vars.len();
365                let mut j = other.vars.len();
366
367                while i > 0 && j > 0 {
368                    i -= 1;
369                    j -= 1;
370
371                    match self.vars[i].var.cmp(&other.vars[j].var) {
372                        Ordering::Less => return Ordering::Less,
373                        Ordering::Greater => return Ordering::Greater,
374                        Ordering::Equal => match self.vars[i].power.cmp(&other.vars[j].power) {
375                            Ordering::Equal => {}
376                            Ordering::Less => return Ordering::Greater,
377                            Ordering::Greater => return Ordering::Less,
378                        },
379                    }
380                }
381
382                if i > 0 {
383                    Ordering::Less
384                } else if j > 0 {
385                    Ordering::Greater
386                } else {
387                    Ordering::Equal
388                }
389            }
390            ord => ord,
391        }
392    }
393}
394
395fn compute_monomial_hash(vars: &[VarPower]) -> u64 {
396    use core::hash::{BuildHasher, BuildHasherDefault};
397    use rustc_hash::FxHasher;
398    let mut hasher = BuildHasherDefault::<FxHasher>::default().build_hasher();
399    for vp in vars {
400        vp.hash(&mut hasher);
401    }
402    hasher.finish()
403}
404
405impl Hash for Monomial {
406    fn hash<H: Hasher>(&self, state: &mut H) {
407        self.hash.hash(state);
408    }
409}
410
411impl fmt::Debug for Monomial {
412    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
413        if self.is_unit() {
414            write!(f, "1")
415        } else {
416            for (i, vp) in self.vars.iter().enumerate() {
417                if i > 0 {
418                    write!(f, "*")?;
419                }
420                if vp.power == 1 {
421                    write!(f, "x{}", vp.var)?;
422                } else {
423                    write!(f, "x{}^{}", vp.var, vp.power)?;
424                }
425            }
426            Ok(())
427        }
428    }
429}
430
431impl fmt::Display for Monomial {
432    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
433        fmt::Debug::fmt(self, f)
434    }
435}
436
437/// A term is a coefficient multiplied by a monomial.
438#[derive(Clone, PartialEq, Eq)]
439pub struct Term {
440    /// The coefficient of the term.
441    pub coeff: BigRational,
442    /// The monomial part of the term.
443    pub monomial: Monomial,
444}
445
446impl Term {
447    /// Create a new term.
448    #[inline]
449    pub fn new(coeff: BigRational, monomial: Monomial) -> Self {
450        Self { coeff, monomial }
451    }
452
453    /// Create a constant term.
454    #[inline]
455    pub fn constant(c: BigRational) -> Self {
456        Self::new(c, Monomial::unit())
457    }
458
459    /// Create a term from a single variable.
460    #[inline]
461    pub fn from_var(var: Var) -> Self {
462        Self::new(BigRational::one(), Monomial::from_var(var))
463    }
464
465    /// Check if this term is zero.
466    #[inline]
467    pub fn is_zero(&self) -> bool {
468        self.coeff.is_zero()
469    }
470
471    /// Check if this is a constant term.
472    #[inline]
473    pub fn is_constant(&self) -> bool {
474        self.monomial.is_unit()
475    }
476}
477
478impl fmt::Debug for Term {
479    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
480        if self.monomial.is_unit() {
481            write!(f, "{}", self.coeff)
482        } else if self.coeff.is_one() {
483            write!(f, "{:?}", self.monomial)
484        } else if self.coeff == -BigRational::one() {
485            write!(f, "-{:?}", self.monomial)
486        } else {
487            write!(f, "{}*{:?}", self.coeff, self.monomial)
488        }
489    }
490}
491
492/// Monomial ordering for polynomial canonicalization.
493#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
494pub enum MonomialOrder {
495    /// Lexicographic order.
496    Lex,
497    /// Graded lexicographic order.
498    #[default]
499    GrLex,
500    /// Graded reverse lexicographic order.
501    GRevLex,
502}
503
504impl MonomialOrder {
505    /// Compare two monomials using this ordering.
506    pub fn compare(&self, a: &Monomial, b: &Monomial) -> Ordering {
507        match self {
508            MonomialOrder::Lex => a.lex_cmp(b),
509            MonomialOrder::GrLex => a.grlex_cmp(b),
510            MonomialOrder::GRevLex => a.grevlex_cmp(b),
511        }
512    }
513}