magnesia/algebra/
polynomial.rs

1use super::ops_with_ref::*;
2use super::ring::{One, Ring, Zero};
3use std::ops::{AddAssign, Mul, SubAssign};
4
5/**
6A polynomial of arbitrary order.
7
8# Example
9```
10# use magnesia::algebra::Polynomial;
11let p = Polynomial::from_coefficients(vec![1, 2, 3]); // p(x) = 1 + 2*x + 3*x^2
12assert_eq!(p.eval(2), 1 + 2*2 + 3*4);
13```
14*/
15#[derive(Debug, Clone, Eq, PartialEq)]
16pub struct Polynomial<T: Ring> {
17    a: Vec<T>,
18}
19
20impl<T: Ring> Polynomial<T> {
21    /// Creates a new `Polynomial` from a vector of coefficients.
22    ///
23    /// Note that the array must start with the coefficients that are
24    /// multiplied with the lower powers of `x` and ends with the coefficients
25    /// that are multiplied with the highest powers of `x`.
26    /// This may be somewhat unintuitive, since in maths polynomials are
27    /// usually written in reverse direction.
28    ///
29    /// # Example
30    /// ```
31    /// # use magnesia::algebra::Polynomial;
32    /// // p(x) = 1 + 2*x + 3*x^2
33    /// let p = Polynomial::from_coefficients(vec![1, 2, 3]);
34    /// ```
35    pub fn from_coefficients(coefficients: Vec<T>) -> Polynomial<T> {
36        Polynomial { a: coefficients }
37    }
38
39    /// Evaluates a polynomial given the reference to a value.
40    ///
41    /// # Example
42    /// ```
43    /// # use magnesia::algebra::Polynomial;
44    /// // p(x) = 1 + 2*x + 3*x^2
45    /// let p = Polynomial::from_coefficients(vec![1, 2, 3]);
46    /// let two = 2;
47    /// assert_eq!(p.eval_ref(&two), 1 + 2*two + 3*two*two);
48    /// ```
49    pub fn eval_ref(&self, x: &T) -> T {
50        let mut y = T::zero();
51        for a in self.a.iter().rev() {
52            y = y.mul_refs(x);
53            y += a;
54        }
55        y
56    }
57
58    /// Evaluates a polynomial given a value.
59    ///
60    /// # Example
61    /// ```
62    /// # use magnesia::algebra::Polynomial;
63    /// // p(x) = 1 + 2*x + 3*x^2
64    /// let p = Polynomial::from_coefficients(vec![1, 2, 3]);
65    /// assert_eq!(p.eval(2), 1 + 2*2 + 3*4);
66    /// ```
67    pub fn eval(&self, x: T) -> T {
68        self.eval_ref(&x)
69    }
70}
71
72impl<T: Ring> Zero for Polynomial<T> {
73    fn zero() -> Self {
74        Polynomial { a: Vec::new() }
75    }
76}
77
78#[test]
79fn test_zero_polynomial_for_ints() {
80    let p = Polynomial::zero();
81    assert_eq!(p.eval(5), 0);
82}
83
84impl<T: Ring> One for Polynomial<T> {
85    fn one() -> Self {
86        Polynomial { a: vec![T::one()] }
87    }
88}
89
90#[test]
91fn test_one_polynomial_for_ints() {
92    let p = Polynomial::one();
93    assert_eq!(p.eval(5), 1);
94}
95
96impl<T: Ring> std::ops::Add<Polynomial<T>> for Polynomial<T> {
97    type Output = Self;
98
99    fn add(mut self, rhs: Polynomial<T>) -> Self {
100        if self.a.len() < rhs.a.len() {
101            rhs.add(self)
102        } else {
103            for (i, x) in rhs.a.into_iter().enumerate() {
104                self.a[i] += &x;
105            }
106            self
107        }
108    }
109}
110
111#[test]
112fn test_add_polynomials_for_ints() {
113    let p = Polynomial::from_coefficients(vec![1, 2, 3]);
114    let q = Polynomial::from_coefficients(vec![2, 3, 4]);
115    let r = p + q;
116    assert_eq!(r, Polynomial::from_coefficients(vec![3, 5, 7]));
117}
118
119impl<T: Ring> std::ops::AddAssign for Polynomial<T> {
120    fn add_assign(&mut self, rhs: Self) {
121        let my_len = self.a.len();
122        if self.a.len() < rhs.a.len() {
123            self.a.reserve(rhs.a.len() - self.a.len());
124        }
125        for (i, x) in rhs.a.into_iter().enumerate() {
126            if i < my_len {
127                self.a[i] += &x;
128            } else {
129                self.a.push(x);
130            }
131        }
132    }
133}
134
135#[test]
136fn test_add_assign_for_polynomials_of_ints() {
137    let mut p = Polynomial::from_coefficients(vec![1, 2, 3]);
138    let q = Polynomial::from_coefficients(vec![2, 3, 4]);
139    p += q;
140    assert_eq!(p, Polynomial::from_coefficients(vec![3, 5, 7]));
141}
142
143impl<T: Ring> AddAssign<&Polynomial<T>> for Polynomial<T> {
144    fn add_assign(&mut self, rhs: &Self) {
145        if self.a.len() < rhs.a.len() {
146            self.a.resize_with(rhs.a.len(), || T::zero());
147        }
148        for i in 0..rhs.a.len() {
149            self.a[i] += &rhs.a[i];
150        }
151    }
152}
153
154#[test]
155fn test_add_assign_for_ref_polynomials_of_ints() {
156    let mut p = Polynomial::from_coefficients(vec![1, 2, 3]);
157    let q = Polynomial::from_coefficients(vec![2, 3, 4]);
158    p += &q;
159    assert_eq!(p, Polynomial::from_coefficients(vec![3, 5, 7]));
160}
161
162impl<T: Ring> std::ops::Sub for Polynomial<T> {
163    type Output = Self;
164
165    fn sub(self, rhs: Self) -> Self {
166        self + -rhs
167    }
168}
169
170#[test]
171fn test_sub_for_polynomials_of_int() {
172    let p = Polynomial::from_coefficients(vec![1, 2, 3]);
173    let q = Polynomial::from_coefficients(vec![2, 3, 5]);
174    let r = p - q;
175    assert_eq!(r, Polynomial::from_coefficients(vec![-1, -1, -2]));
176}
177
178impl<T: Ring> std::ops::SubAssign for Polynomial<T> {
179    fn sub_assign(&mut self, rhs: Self) {
180        *self += -rhs
181    }
182}
183
184#[test]
185fn test_sub_assign_for_polynomials_of_ints() {
186    let mut p = Polynomial::from_coefficients(vec![1, 2, 3]);
187    let q = Polynomial::from_coefficients(vec![2, 3, 5]);
188    p -= q;
189    assert_eq!(p, Polynomial::from_coefficients(vec![-1, -1, -2]));
190}
191
192impl<T: Ring> SubAssign<&Polynomial<T>> for Polynomial<T> {
193    fn sub_assign(&mut self, rhs: &Self) {
194        if self.a.len() < rhs.a.len() {
195            self.a.resize_with(rhs.a.len(), || T::zero());
196        }
197        for i in 0..rhs.a.len() {
198            self.a[i] -= &rhs.a[i];
199        }
200    }
201}
202
203#[test]
204fn test_sub_assign_for_ref_polynomials_of_ints() {
205    let mut p = Polynomial::from_coefficients(vec![1, 2, 3]);
206    let q = Polynomial::from_coefficients(vec![2, 3, 5]);
207    p -= &q;
208    assert_eq!(p, Polynomial::from_coefficients(vec![-1, -1, -2]));
209}
210
211impl<T: Ring> std::ops::Mul for Polynomial<T> {
212    type Output = Self;
213
214    fn mul(self, rhs: Self) -> Self {
215        self.mul_refs(&rhs)
216    }
217}
218
219#[test]
220fn test_mul_for_polynomials_of_ints() {
221    // p(x) = 1 + 2*x + 3*x^2
222    // q(x) = x
223    let p = Polynomial::from_coefficients(vec![1, 2, 3]);
224    let q = Polynomial::from_coefficients(vec![0, 1]);
225    let r = p * q;
226    assert_eq!(r, Polynomial::from_coefficients(vec![0, 1, 2, 3]));
227}
228
229impl<T: Ring> std::ops::MulAssign for Polynomial<T> {
230    fn mul_assign(&mut self, rhs: Self) {
231        *self = self.mul_refs(&rhs);
232    }
233}
234
235#[test]
236fn test_mul_assign_for_polynomials_of_ints() {
237    // p(x) = 1 + 2*x + 3*x^2
238    // q(x) = x
239    let mut p = Polynomial::from_coefficients(vec![1, 2, 3]);
240    let q = Polynomial::from_coefficients(vec![0, 1]);
241    p *= q;
242    assert_eq!(p, Polynomial::from_coefficients(vec![0, 1, 2, 3]));
243}
244
245impl<T: Ring> Mul for &Polynomial<T> {
246    type Output = Polynomial<T>;
247
248    fn mul(self, rhs: Self) -> Polynomial<T> {
249        if self.a.is_empty() || rhs.a.is_empty() {
250            return <Polynomial<T>>::zero();
251        }
252        let mut a = vec![T::zero(); self.a.len() + rhs.a.len() - 1];
253        for i in 0..self.a.len() {
254            for j in 0..rhs.a.len() {
255                let prod = self.a[i].mul_refs(&rhs.a[j]);
256                a[i + j] += &prod;
257            }
258        }
259        <Polynomial<T>>::from_coefficients(a)
260    }
261}
262
263#[test]
264fn test_mul_for_ref_polynomials_of_ints() {
265    // p(x) = 1 + 2*x + 3*x^2
266    // q(x) = x
267    let p = Polynomial::from_coefficients(vec![1, 2, 3]);
268    let q = Polynomial::from_coefficients(vec![0, 1]);
269    let r = &p * &q;
270    assert_eq!(r, Polynomial::from_coefficients(vec![0, 1, 2, 3]));
271}
272
273impl<T: Ring> std::ops::Neg for Polynomial<T> {
274    type Output = Self;
275
276    fn neg(self) -> Self {
277        Self {
278            a: self.a.into_iter().map(|x| -x).collect(),
279        }
280    }
281}
282
283#[test]
284fn test_neg_for_polynomials_of_ints() {
285    let p = Polynomial::from_coefficients(vec![1, 2, 3]);
286    let q = -p;
287    assert_eq!(q, Polynomial::from_coefficients(vec![-1, -2, -3]));
288}
289
290impl<T: Ring> Ring for Polynomial<T> {}