summaryrefslogtreecommitdiff
path: root/src/poly.rs
diff options
context:
space:
mode:
Diffstat (limited to 'src/poly.rs')
-rw-r--r--src/poly.rs148
1 files changed, 148 insertions, 0 deletions
diff --git a/src/poly.rs b/src/poly.rs
new file mode 100644
index 0000000..2128d38
--- /dev/null
+++ b/src/poly.rs
@@ -0,0 +1,148 @@
+use std::cmp::Ordering;
+
+pub type Number = f32;
+
+pub type Poly = Vec<Number>;
+
+pub fn eval_poly(p: &Poly, t: Number) -> Number {
+ let mut r = 0.0;
+ for i in 0..p.len() {
+ r += p[i] * t.powi(i as i32);
+ }
+ r
+}
+
+pub fn sub(a: &Poly, b: &Poly) -> Poly {
+ let mut r = a.clone();
+ for i in 0..r.len() {
+ r[i] -= b[i];
+ }
+ r
+}
+
+pub fn skewed_sum(a: Poly, b: Poly) -> Poly {
+ let mut r = a.clone();
+ for i in 0..r.len() - 1 {
+ r[i + 1] += b[i];
+ }
+ r.push(b[b.len() - 1]);
+ r
+}
+
+fn is_zero(p: &Poly) -> bool {
+ for i in 0..p.len() {
+ if p[i] != 0.0 {
+ return false
+ }
+ }
+ true
+}
+
+fn degree(p: &Poly) -> usize {
+ let mut i = p.len() - 1;
+ while p[i] == 0.0 && i > 0 {
+ i -= 1;
+ }
+ i
+}
+
+fn prod(p: &Poly, n: Number) -> Poly {
+ let mut r = p.clone();
+ for i in 0..r.len() {
+ r[i] *= n;
+ }
+ r
+}
+
+fn shift(p: &Poly, amount: usize) -> Poly {
+ let mut r = vec![0.0; p.len()];
+ for i in 0..p.len() {
+ if i + amount < r.len() {
+ r[i + amount] = p[i];
+ }
+ }
+ r
+}
+
+fn rem(a: &Poly, b: &Poly) -> Poly {
+ let an = a[degree(a)];
+ let bn = b[degree(b)];
+ match degree(a).cmp(&degree(b)) {
+ Ordering::Equal => sub(a, &prod(b, an / bn)),
+ Ordering::Greater => sub(a, &prod(&shift(b, (degree(a) - degree(b)) as usize), an / bn)),
+ Ordering::Less => vec![0.0; b.len()]
+ }
+}
+
+pub fn gcd(a: &Poly, b: &Poly) -> Poly {
+ let c = rem(a, b);
+ if is_zero(&c) {
+ b.clone()
+ } else {
+ gcd(b, &c)
+ }
+}
+
+#[cfg(test)]
+mod tests {
+ use super::*;
+
+ #[test]
+ fn poly_is_zero() {
+ let p = vec![0.0; 5];
+ assert_eq!(is_zero(&p), true);
+ }
+
+ #[test]
+ fn poly_is_not_zero() {
+ let mut p = vec![0.0; 5];
+ p[2] = 8.0;
+ assert_eq!(is_zero(&p), false);
+ }
+
+ #[test]
+ fn degree_is_five() {
+ let mut p = vec![0.0; 6];
+ p[5] = 2.0;
+ assert_eq!(degree(&p), 5);
+ }
+
+ #[test]
+ fn degree_is_zero() {
+ let p = vec![0.0; 6];
+ assert_eq!(degree(&p), 0);
+ }
+
+ #[test]
+ fn prod_test() {
+ let p = vec![1.0, 2.0, 3.0, 4.0, 5.0];
+ assert_eq!(prod(&p, 2.0), vec![2.0, 4.0, 6.0, 8.0, 10.0]);
+ }
+
+ #[test]
+ fn shift_test() {
+ let p = vec![1.0, 2.0, 3.0, 0.0, 0.0];
+ assert_eq!(shift(&p, 2), vec![0.0, 0.0, 1.0, 2.0, 3.0]);
+ }
+
+ #[test]
+ fn rem_equal() {
+ let a = vec![6.0, 7.0, 1.0];
+ let b = vec![-6.0, -5.0, 1.0];
+ assert_eq!(rem(&a, &b), vec![12.0, 12.0, 0.0]);
+ }
+
+ #[test]
+ fn rem_greater() {
+ let a = vec![-6.0, -5.0, 1.0];
+ let b = vec![12.0, 12.0, 0.0];
+ assert_eq!(rem(&a, &b), vec![-6.0, -6.0, 0.0]);
+ }
+
+ #[test]
+ fn gcd_test() {
+ let a = vec![6.0, 7.0, 1.0];
+ let b = vec![-6.0, -5.0, 1.0];
+ assert_eq!(gcd(&a, &b), vec![-6.0, -6.0, 0.0]);
+ }
+}