mathhook_core/algebra/groebner/
monomial_order.rs

1//! Monomial Orderings for Gröbner Basis Computation
2//!
3//! Implements standard monomial orderings: lexicographic, graded lexicographic,
4//! and graded reverse lexicographic. These orderings are fundamental to
5//! Gröbner basis algorithms.
6
7use crate::core::{Expression, Number, Symbol};
8use std::cmp::Ordering;
9
10/// Monomial ordering types
11#[derive(Debug, Clone, Copy, PartialEq, Eq)]
12pub enum MonomialOrder {
13    /// Lexicographic ordering (Lex)
14    ///
15    /// Compares monomials by exponents left to right
16    /// x^a1 y^a2 > x^b1 y^b2 if first non-equal exponent ai > bi
17    Lex,
18
19    /// Graded lexicographic ordering (Grlex)
20    ///
21    /// First compare total degree, then use lex for ties
22    Grlex,
23
24    /// Graded reverse lexicographic ordering (Grevlex)
25    ///
26    /// First compare total degree, then reverse lex for ties
27    Grevlex,
28}
29
30/// Trait for comparing monomials according to a specific ordering
31pub trait MonomialOrdering {
32    /// Compare two monomials using this ordering
33    ///
34    /// # Arguments
35    ///
36    /// * `mono1` - First monomial
37    /// * `mono2` - Second monomial
38    /// * `variables` - Ordered list of variables
39    ///
40    /// # Returns
41    ///
42    /// Ordering result (Less, Equal, Greater)
43    fn compare_monomials(
44        &self,
45        mono1: &Expression,
46        mono2: &Expression,
47        variables: &[Symbol],
48    ) -> Ordering;
49
50    /// Get the leading monomial of a polynomial
51    ///
52    /// # Arguments
53    ///
54    /// * `poly` - Polynomial expression
55    /// * `variables` - Ordered list of variables
56    ///
57    /// # Returns
58    ///
59    /// The leading monomial according to this ordering
60    fn leading_monomial(&self, poly: &Expression, variables: &[Symbol]) -> Expression;
61
62    /// Get the leading coefficient of a polynomial
63    ///
64    /// # Arguments
65    ///
66    /// * `poly` - Polynomial expression
67    /// * `variables` - Ordered list of variables
68    ///
69    /// # Returns
70    ///
71    /// The coefficient of the leading monomial
72    fn leading_coefficient(&self, poly: &Expression, variables: &[Symbol]) -> Expression;
73}
74
75impl MonomialOrdering for MonomialOrder {
76    fn compare_monomials(
77        &self,
78        mono1: &Expression,
79        mono2: &Expression,
80        variables: &[Symbol],
81    ) -> Ordering {
82        let exp1 = extract_exponents(mono1, variables);
83        let exp2 = extract_exponents(mono2, variables);
84
85        match self {
86            MonomialOrder::Lex => compare_lex(&exp1, &exp2),
87            MonomialOrder::Grlex => compare_grlex(&exp1, &exp2),
88            MonomialOrder::Grevlex => compare_grevlex(&exp1, &exp2),
89        }
90    }
91
92    fn leading_monomial(&self, poly: &Expression, variables: &[Symbol]) -> Expression {
93        match poly {
94            Expression::Add(terms) => {
95                let mut leading = terms[0].clone();
96                for term in &terms[1..] {
97                    if self.compare_monomials(term, &leading, variables) == Ordering::Greater {
98                        leading = term.clone();
99                    }
100                }
101                extract_monomial_part(&leading)
102            }
103            _ => extract_monomial_part(poly),
104        }
105    }
106
107    fn leading_coefficient(&self, poly: &Expression, variables: &[Symbol]) -> Expression {
108        match poly {
109            Expression::Add(terms) => {
110                let mut leading = terms[0].clone();
111                for term in &terms[1..] {
112                    if self.compare_monomials(term, &leading, variables) == Ordering::Greater {
113                        leading = term.clone();
114                    }
115                }
116                extract_coefficient(&leading)
117            }
118            _ => extract_coefficient(poly),
119        }
120    }
121}
122
123/// Extract exponents of a monomial as a vector
124///
125/// # Arguments
126///
127/// * `mono` - Monomial expression (possibly with coefficient)
128/// * `variables` - Ordered list of variables
129///
130/// # Returns
131///
132/// Vector of exponents in the same order as variables
133fn extract_exponents(mono: &Expression, variables: &[Symbol]) -> Vec<i64> {
134    let mut exponents = vec![0i64; variables.len()];
135    let mono_part = extract_monomial_part(mono);
136
137    match mono_part {
138        Expression::Symbol(s) => {
139            if let Some(idx) = variables.iter().position(|v| v == &s) {
140                exponents[idx] = 1;
141            }
142        }
143        Expression::Pow(base, exp) => {
144            if let Expression::Symbol(s) = base.as_ref() {
145                if let Some(idx) = variables.iter().position(|v| v == s) {
146                    // Extract integer from Number variant
147                    if let Expression::Number(Number::Integer(i)) = exp.as_ref() {
148                        exponents[idx] = *i;
149                    }
150                }
151            }
152        }
153        Expression::Mul(factors) => {
154            // Dereference Box<Vec<Expression>> to iterate
155            for factor in factors.iter() {
156                let factor_exp = extract_exponents(factor, variables);
157                for (i, e) in factor_exp.iter().enumerate() {
158                    exponents[i] += e;
159                }
160            }
161        }
162        _ => {}
163    }
164
165    exponents
166}
167
168/// Extract the monomial part (without coefficient) of a term
169fn extract_monomial_part(term: &Expression) -> Expression {
170    match term {
171        Expression::Mul(factors) => {
172            let mut mono_factors = Vec::new();
173            // Dereference Box<Vec<Expression>> to iterate
174            for factor in factors.iter() {
175                if !is_numeric(factor) {
176                    mono_factors.push(factor.clone());
177                }
178            }
179            if mono_factors.is_empty() {
180                Expression::integer(1)
181            } else if mono_factors.len() == 1 {
182                mono_factors[0].clone()
183            } else {
184                Expression::mul(mono_factors)
185            }
186        }
187        _ => {
188            if is_numeric(term) {
189                Expression::integer(1)
190            } else {
191                term.clone()
192            }
193        }
194    }
195}
196
197/// Extract the coefficient of a term
198fn extract_coefficient(term: &Expression) -> Expression {
199    match term {
200        Expression::Mul(factors) => {
201            let mut coeffs = Vec::new();
202            // Dereference Box<Vec<Expression>> to iterate
203            for factor in factors.iter() {
204                if is_numeric(factor) {
205                    coeffs.push(factor.clone());
206                }
207            }
208            if coeffs.is_empty() {
209                Expression::integer(1)
210            } else if coeffs.len() == 1 {
211                coeffs[0].clone()
212            } else {
213                Expression::mul(coeffs)
214            }
215        }
216        Expression::Number(_) => term.clone(),
217        _ => Expression::integer(1),
218    }
219}
220
221/// Check if expression is a numeric constant
222fn is_numeric(expr: &Expression) -> bool {
223    matches!(expr, Expression::Number(_))
224}
225
226/// Lexicographic comparison of exponent vectors
227fn compare_lex(exp1: &[i64], exp2: &[i64]) -> Ordering {
228    for (e1, e2) in exp1.iter().zip(exp2.iter()) {
229        match e1.cmp(e2) {
230            Ordering::Equal => continue,
231            other => return other,
232        }
233    }
234    Ordering::Equal
235}
236
237/// Graded lexicographic comparison
238fn compare_grlex(exp1: &[i64], exp2: &[i64]) -> Ordering {
239    let deg1: i64 = exp1.iter().sum();
240    let deg2: i64 = exp2.iter().sum();
241
242    match deg1.cmp(&deg2) {
243        Ordering::Equal => compare_lex(exp1, exp2),
244        other => other,
245    }
246}
247
248/// Graded reverse lexicographic comparison
249fn compare_grevlex(exp1: &[i64], exp2: &[i64]) -> Ordering {
250    let deg1: i64 = exp1.iter().sum();
251    let deg2: i64 = exp2.iter().sum();
252
253    match deg1.cmp(&deg2) {
254        Ordering::Equal => {
255            for (e1, e2) in exp1.iter().zip(exp2.iter()).rev() {
256                match e2.cmp(e1) {
257                    Ordering::Equal => continue,
258                    other => return other,
259                }
260            }
261            Ordering::Equal
262        }
263        other => other,
264    }
265}
266
267#[cfg(test)]
268mod tests {
269    use super::*;
270    use crate::symbol;
271
272    #[test]
273    fn test_extract_exponents_simple() {
274        let x = symbol!(x);
275        let y = symbol!(y);
276        let vars = vec![x.clone(), y.clone()];
277
278        let mono = Expression::mul(vec![
279            Expression::pow(Expression::symbol(x.clone()), Expression::integer(2)),
280            Expression::symbol(y.clone()),
281        ]);
282
283        let exp = extract_exponents(&mono, &vars);
284        assert_eq!(exp, vec![2, 1]);
285    }
286
287    #[test]
288    fn test_lex_ordering() {
289        let x = symbol!(x);
290        let y = symbol!(y);
291        let vars = vec![x.clone(), y.clone()];
292
293        let mono1 = Expression::pow(Expression::symbol(x.clone()), Expression::integer(2));
294        let mono2 = Expression::pow(Expression::symbol(y.clone()), Expression::integer(3));
295
296        let order = MonomialOrder::Lex;
297        let cmp = order.compare_monomials(&mono1, &mono2, &vars);
298
299        assert_eq!(cmp, Ordering::Greater);
300    }
301
302    #[test]
303    fn test_grlex_ordering() {
304        let x = symbol!(x);
305        let y = symbol!(y);
306        let vars = vec![x.clone(), y.clone()];
307
308        let mono1 = Expression::pow(Expression::symbol(x.clone()), Expression::integer(1));
309        let mono2 = Expression::pow(Expression::symbol(y.clone()), Expression::integer(2));
310
311        let order = MonomialOrder::Grlex;
312        let cmp = order.compare_monomials(&mono1, &mono2, &vars);
313
314        assert_eq!(cmp, Ordering::Less);
315    }
316
317    #[test]
318    fn test_grevlex_ordering() {
319        let x = symbol!(x);
320        let y = symbol!(y);
321        let vars = vec![x.clone(), y.clone()];
322
323        let mono1 = Expression::mul(vec![
324            Expression::symbol(x.clone()),
325            Expression::symbol(y.clone()),
326        ]);
327        let mono2 = Expression::pow(Expression::symbol(x.clone()), Expression::integer(2));
328
329        let order = MonomialOrder::Grevlex;
330        let cmp = order.compare_monomials(&mono1, &mono2, &vars);
331
332        assert_eq!(cmp, Ordering::Less);
333    }
334
335    #[test]
336    fn test_leading_monomial() {
337        let x = symbol!(x);
338        let y = symbol!(y);
339        let vars = vec![x.clone(), y.clone()];
340
341        let poly = Expression::add(vec![
342            Expression::pow(Expression::symbol(x.clone()), Expression::integer(2)),
343            Expression::mul(vec![Expression::integer(3), Expression::symbol(y.clone())]),
344        ]);
345
346        let order = MonomialOrder::Lex;
347        let lm = order.leading_monomial(&poly, &vars);
348
349        assert_eq!(
350            lm,
351            Expression::pow(Expression::symbol(x), Expression::integer(2))
352        );
353    }
354}