Skip to main content

ps_ecc/polynomial/methods/
div_rem_inplace.rs

1use crate::{
2    finite_field::{inv, mul, sub},
3    Polynomial, PolynomialDivError,
4};
5
6impl Polynomial {
7    /// In-place polynomial division: `self` becomes the remainder, quotient is output.
8    ///
9    /// After this method, `self` contains the remainder and `quotient` contains
10    /// the quotient of dividing the original `self` by `divisor`.
11    ///
12    /// # Errors
13    ///
14    /// Returns `PolynomialDivError::ZeroDivisor` if the divisor is the zero polynomial.
15    /// Returns `PolynomialDivError::GFError` if field inversion fails.
16    #[inline]
17    pub fn div_rem_inplace(
18        &mut self,
19        divisor: &Self,
20        quotient: &mut Self,
21    ) -> Result<(), PolynomialDivError> {
22        let divisor_deg = divisor.degree as usize;
23
24        // Check for zero divisor
25        if divisor_deg == 0 && divisor.coefficients[0] == 0 {
26            return Err(PolynomialDivError::ZeroDivisor);
27        }
28
29        // Precompute inverse of leading coefficient
30        let divisor_lc_inv = inv(divisor.coefficients[divisor_deg])?.get();
31
32        // Clear quotient
33        quotient.coefficients.fill(0);
34        quotient.degree = 0;
35
36        let mut deg = self.degree as usize;
37
38        // If dividend degree < divisor degree, remainder is dividend, quotient is 0
39        if deg < divisor_deg {
40            return Ok(());
41        }
42
43        // Set quotient degree
44        #[allow(clippy::cast_possible_truncation)]
45        {
46            quotient.degree = (deg - divisor_deg) as u8;
47        }
48
49        while deg >= divisor_deg {
50            let lead_coef = self.coefficients[deg];
51
52            if lead_coef == 0 {
53                if deg == 0 {
54                    break;
55                }
56
57                deg -= 1;
58                continue;
59            }
60
61            let ratio = mul(lead_coef, divisor_lc_inv);
62            let shift = deg - divisor_deg;
63
64            quotient.coefficients[shift] = ratio;
65            self.coefficients[deg] = 0;
66
67            // Subtract ratio * divisor from self (in GF(2^8), sub == xor)
68            for i in 0..divisor_deg {
69                self.coefficients[shift + i] = sub(
70                    self.coefficients[shift + i],
71                    mul(ratio, divisor.coefficients[i]),
72                );
73            }
74
75            if deg == 0 {
76                break;
77            }
78
79            deg -= 1;
80        }
81
82        // Update degrees
83        self.trim_degree();
84        quotient.trim_degree();
85
86        Ok(())
87    }
88}
89
90#[cfg(test)]
91#[allow(clippy::expect_used)]
92mod tests {
93    use crate::{finite_field::mul, Polynomial};
94
95    #[test]
96    fn div_rem_inplace_basic() {
97        // (x^2 + 1) / (x + 1) in GF(256)
98        // = (x + 1) with remainder 0
99        let mut dividend = Polynomial::try_from(&[1u8, 0, 1][..]).expect("valid polynomial");
100
101        let divisor = Polynomial::try_from(&[1u8, 1][..]).expect("valid polynomial");
102        let mut quotient = Polynomial::default();
103
104        dividend
105            .div_rem_inplace(&divisor, &mut quotient)
106            .expect("division succeeds");
107
108        assert_eq!(quotient.coefficients(), &[1, 1]);
109        assert_eq!(dividend.coefficients(), &[0]);
110    }
111
112    #[test]
113    fn div_rem_inplace_with_remainder() {
114        // (x^2 + x + 1) / (x + 1) should have a remainder
115        let mut dividend = Polynomial::try_from(&[1u8, 1, 1][..]).expect("valid polynomial");
116
117        let divisor = Polynomial::try_from(&[1u8, 1][..]).expect("valid polynomial");
118
119        let original = Polynomial::try_from(&[1u8, 1, 1][..]).expect("valid polynomial");
120
121        let mut quotient = Polynomial::default();
122
123        dividend
124            .div_rem_inplace(&divisor, &mut quotient)
125            .expect("division succeeds");
126
127        // Verify: quotient * divisor + remainder = original
128        let mut reconstructed = Polynomial::default();
129
130        // Multiply quotient * divisor
131        for i in 0..=quotient.degree as usize {
132            for j in 0..=divisor.degree as usize {
133                reconstructed.coefficients[i + j] ^=
134                    mul(quotient.coefficients[i], divisor.coefficients[j]);
135            }
136        }
137
138        reconstructed.degree = quotient.degree + divisor.degree;
139        reconstructed.trim_degree();
140
141        // Add remainder
142        for i in 0..=dividend.degree as usize {
143            reconstructed.coefficients[i] ^= dividend.coefficients[i];
144        }
145
146        reconstructed.trim_degree();
147
148        assert_eq!(reconstructed, original);
149    }
150
151    #[test]
152    fn div_rem_inplace_smaller_dividend() {
153        // dividend degree < divisor degree => remainder = dividend, quotient = 0
154        let mut dividend = Polynomial::try_from(&[5u8, 3][..]).expect("valid polynomial");
155
156        let divisor = Polynomial::try_from(&[1u8, 2, 3][..]).expect("valid polynomial");
157
158        let original = Polynomial::try_from(&[5u8, 3][..]).expect("valid polynomial");
159        let mut quotient = Polynomial::default();
160
161        dividend
162            .div_rem_inplace(&divisor, &mut quotient)
163            .expect("division succeeds");
164
165        assert_eq!(quotient.coefficients(), &[0]);
166        assert_eq!(dividend, original);
167    }
168
169    #[test]
170    fn div_rem_inplace_zero_divisor() {
171        let mut dividend = Polynomial::try_from(&[1u8, 2, 3][..]).expect("valid polynomial");
172
173        let divisor = Polynomial::default();
174        let mut quotient = Polynomial::default();
175
176        let result = dividend.div_rem_inplace(&divisor, &mut quotient);
177
178        assert!(result.is_err());
179    }
180}