add polynomial division

src/poly/tests.rs

@@ -1,4 +1,4 @@
-use super::flat::{Mono, Poly};
+use super::flat::{lex_cmp, Mono, Poly};
 use super::var::StaticVar;
 
 #[test]
@@ -99,3 +99,130 @@ fn test_poly_sub() {
     let b: Poly<StaticVar> = [(4, [("x", 2)]), (1, [("y", 1)])].into();
     assert_eq!(a - b, Poly::default());
 }
+
+#[test]
+fn test_lex_cmp() {
+    use std::cmp::Ordering;
+
+    let x2: Mono<StaticVar> = [("x", 2)].into();
+    let xy: Mono<StaticVar> = [("x", 1), ("y", 1)].into();
+    let y2: Mono<StaticVar> = [("y", 2)].into();
+    let x: Mono<StaticVar> = [("x", 1)].into();
+    let one: Mono<StaticVar> = Mono { term: vec![] };
+
+    // x² > xy (x exponent 2 vs 1)
+    assert_eq!(lex_cmp(&x2, &xy), Ordering::Greater);
+    // x > y² (x has higher priority)
+    assert_eq!(lex_cmp(&x, &y2), Ordering::Greater);
+    // xy > y² (x present in xy but not y²)
+    assert_eq!(lex_cmp(&xy, &y2), Ordering::Greater);
+    // 1 < x
+    assert_eq!(lex_cmp(&one, &x), Ordering::Less);
+    // reflexive
+    assert_eq!(lex_cmp(&x2, &x2), Ordering::Equal);
+}
+
+#[test]
+fn test_mono_div() {
+    // x² / x = x
+    let a: Mono<StaticVar> = [("x", 2)].into();
+    let b: Mono<StaticVar> = [("x", 1)].into();
+    assert_eq!(a.div(&b), Mono::from([("x", 1)]));
+
+    // x²y / xy = x
+    let a: Mono<StaticVar> = [("x", 2), ("y", 1)].into();
+    let b: Mono<StaticVar> = [("x", 1), ("y", 1)].into();
+    assert_eq!(a.div(&b), Mono::from([("x", 1)]));
+
+    // x²y / y = x²
+    let a: Mono<StaticVar> = [("x", 2), ("y", 1)].into();
+    let b: Mono<StaticVar> = [("y", 1)].into();
+    assert_eq!(a.div(&b), Mono::from([("x", 2)]));
+
+    // x / x = 1
+    let a: Mono<StaticVar> = [("x", 1)].into();
+    let b: Mono<StaticVar> = [("x", 1)].into();
+    assert_eq!(a.div(&b), Mono { term: vec![] });
+}
+
+fn make_const_poly(c: i32) -> Poly<StaticVar> {
+    Poly { mono: [(Mono { term: vec![] }, c)].into_iter().collect() }
+}
+
+fn verify_div_rem(f: Poly<StaticVar>, g: &Poly<StaticVar>, d: u32, q: Poly<StaticVar>, r: Poly<StaticVar>) {
+    // lc(g)^d * f == q * g + r
+    let (_, lc_g) = g.leading_term_lex().unwrap();
+    let lhs = lc_g.pow(d) * f;
+    let rhs = q * g.clone() + r;
+    assert_eq!(lhs, rhs);
+}
+
+#[test]
+fn test_div_rem() {
+    // x³ / x² = x, r=0, d=0
+    let f: Poly<StaticVar> = [(1, [("x", 3)])].into();
+    let g: Poly<StaticVar> = [(1, [("x", 2)])].into();
+    let expected_q: Poly<StaticVar> = [(1, [("x", 1)])].into();
+    let (d, q, r) = f.clone().div_rem(&g);
+    assert_eq!(q, expected_q);
+    assert!(r.is_zero());
+    verify_div_rem(f, &g, d, q, r);
+
+    // (x³ + x²y) / x² = x + y, r=0
+    // f = x²(x + y), g = x²
+    let f: Poly<StaticVar> = [
+        (1i32, Mono::from([("x", 3u32)])),
+        (1i32, Mono::from([("x", 2u32), ("y", 1u32)])),
+    ].into_iter().collect();
+    let g: Poly<StaticVar> = [(1, [("x", 2)])].into();
+    let expected_q: Poly<StaticVar> = [(1, [("x", 1)]), (1, [("y", 1)])].into();
+    let (d, q, r) = f.clone().div_rem(&g);
+    assert_eq!(q, expected_q);
+    assert!(r.is_zero());
+    verify_div_rem(f, &g, d, q, r);
+
+    // (x³ + y) / x² = x, r=y
+    // LT(x³) divisible by x², LT(y) is not
+    let f: Poly<StaticVar> = [(1, [("x", 3)]), (1, [("y", 1)])].into();
+    let g: Poly<StaticVar> = [(1, [("x", 2)])].into();
+    let expected_q: Poly<StaticVar> = [(1, [("x", 1)])].into();
+    let expected_r: Poly<StaticVar> = [(1, [("y", 1)])].into();
+    let (d, q, r) = f.clone().div_rem(&g);
+    assert_eq!(q, expected_q);
+    assert_eq!(r, expected_r);
+    verify_div_rem(f, &g, d, q, r);
+
+    // 3x² / (2x): 2 ∤ 3, needs pseudo-division
+    // 2¹ · 3x² = 3x · 2x  =>  d=1, q=3x, r=0
+    let f: Poly<StaticVar> = [(3, [("x", 2)])].into();
+    let g: Poly<StaticVar> = [(2, [("x", 1)])].into();
+    let expected_q: Poly<StaticVar> = [(3, [("x", 1)])].into();
+    let (d, q, r) = f.clone().div_rem(&g);
+    assert_eq!(d, 1);
+    assert_eq!(q, expected_q);
+    assert!(r.is_zero());
+    verify_div_rem(f, &g, d, q, r);
+
+    // 6xy / 2 = 3xy, r=0, d=0
+    let f: Poly<StaticVar> = [(6, [("x", 1), ("y", 1)])].into();
+    let g = make_const_poly(2);
+    let expected_q: Poly<StaticVar> = [(3, [("x", 1), ("y", 1)])].into();
+    let (d, q, r) = f.clone().div_rem(&g);
+    assert_eq!(d, 0);
+    assert_eq!(q, expected_q);
+    assert!(r.is_zero());
+    verify_div_rem(f, &g, d, q, r);
+
+    // (x² + xy) / (x + y) = x, r=0
+    // f = x(x + y), g = x + y
+    let f: Poly<StaticVar> = [
+        (1i32, Mono::from([("x", 2u32)])),
+        (1i32, Mono::from([("x", 1u32), ("y", 1u32)])),
+    ].into_iter().collect();
+    let g: Poly<StaticVar> = [(1, [("x", 1)]), (1, [("y", 1)])].into();
+    let expected_q: Poly<StaticVar> = [(1, [("x", 1)])].into();
+    let (d, q, r) = f.clone().div_rem(&g);
+    assert_eq!(q, expected_q);
+    assert!(r.is_zero());
+    verify_div_rem(f, &g, d, q, r);
+}