Struct Poly 

pub struct Poly<const N: usize> {
    pub(crate) coeffs: [f64; N],
}

A polynomial whose degree is known at compile-time.

Although this supports polynomials of arbitrary degree, it is intended for low-degree polynomials. For example, the coefficients are stored in an array, and so they will be stack-allocated (unless you Box the Poly, of course) tend to be copied around.

Polynomial multiplication is not yet implemented, because doing it “nicely” would require const generic expressions: ideally we’d do something like

impl<N, M> Mul<Poly<M>> for Poly<N> {
    type Output = Poly<{M + N - 1}>;
}

It’s possible to work around this with macros, but there are lots of possibilities and I didn’t feel like it was worth the trouble (and the hit to compilation time).

Fields

coeffs: [f64; N]

Implementations

impl Poly<4>

fn critical_points(&self) -> Option<(f64, f64)>

Computes the critical points of this cubic, as long as the discriminant of the derivative is positive. The return values are in increasing order.

Some corner cases worth noting:

• If the discriminant is zero, returns nothing. That is, we don’t find double-roots of the derivative. These aren’t needed anyway, as these critical points are used to construct bracketing intervals for the roots. Double-roots of the derivative are inflection points of the cubic, and they aren’t needed to bracket roots.
• If the derivative is linear or close to it, we might return +/- infinity as one of the roots.
• Unless some input is NaN, we don’t return NaN.

fn rescaled_critical_points(&self) -> Option<(f64, f64)>

fn one_root( &self, lower: f64, upper: f64, lower_val: f64, upper_val: f64, x_error: f64, ) -> f64

fn first_root(self, lower: f64, upper: f64, x_error: f64) -> Option<f64>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 3>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. In fact, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to optimize a quartic, because double-roots of the derivative aren’t local extrema.

pub(crate) fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, _scratch: &mut ArrayVec<f64, M>, )

impl<const N: usize> Poly<N>

pub const fn new(coeffs: [f64; N]) -> Poly<N>

Creates a new polynomial with the provided coefficients.

The constant coefficient comes first, then the linear coefficient, and so on. So if you pass [c, b, a] you’ll get the polynomial a x^2 + b x + c.

pub fn coeffs(&self) -> &[f64; N]

The coefficients of this polynomial.

In the returned array, the coefficient of x^i is at index i.

pub fn eval(&self, x: f64) -> f64

Evaluates this polynomial at a point.

pub fn max_abs_coefficient(&self) -> f64

Returns the largest absolute value of any coefficient.

Always returns a non-negative number, or NaN if some coefficient is NaN.

pub fn is_finite(&self) -> bool

Are all the coefficients finite?

impl Poly<3>

pub fn deriv(&self) -> Poly<2>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<2>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<4>

pub fn deriv(&self) -> Poly<3>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<3>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<5>

pub fn deriv(&self) -> Poly<4>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<4>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<6>

pub fn deriv(&self) -> Poly<5>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<5>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<7>

pub fn deriv(&self) -> Poly<6>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<6>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<8>

pub fn deriv(&self) -> Poly<7>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<7>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<9>

pub fn deriv(&self) -> Poly<8>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<8>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<10>

pub fn deriv(&self) -> Poly<9>

Compute the derivative of this polynomial, as a polynomial with one less coefficient.

pub fn deflate(&self, root: f64) -> Poly<9>

Divide this polynomial by the polynomial x - root, returning the quotient (as a polynomial with one less coefficient) and ignoring the remainder.

If root is actually a root of self (as the name suggests it should be, but this is not actually required), the remainder will be zero. In general, the remainder will be self.eval(root).

impl Poly<5>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 4>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. For example, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to find critical points for optimizing a polynomial, because roots that don’t come with a sign change aren’t local extrema.

fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, scratch: &mut ArrayVec<f64, M>, )

impl Poly<6>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 5>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. For example, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to find critical points for optimizing a polynomial, because roots that don’t come with a sign change aren’t local extrema.

fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, scratch: &mut ArrayVec<f64, M>, )

impl Poly<7>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 6>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. For example, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to find critical points for optimizing a polynomial, because roots that don’t come with a sign change aren’t local extrema.

fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, scratch: &mut ArrayVec<f64, M>, )

impl Poly<8>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 7>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. For example, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to find critical points for optimizing a polynomial, because roots that don’t come with a sign change aren’t local extrema.

fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, scratch: &mut ArrayVec<f64, M>, )

impl Poly<9>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 8>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. For example, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to find critical points for optimizing a polynomial, because roots that don’t come with a sign change aren’t local extrema.

fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, scratch: &mut ArrayVec<f64, M>, )

impl Poly<10>

pub fn roots_between( self, lower: f64, upper: f64, x_error: f64, ) -> ArrayVec<f64, 9>

Computes all roots between lower and upper, to the desired accuracy.

We make no guarantees about multiplicity. For example, if there’s a double-root that isn’t a triple-root (and therefore has no sign change nearby) then there’s a good chance we miss it altogether. This is fine if you’re using this root-finding to find critical points for optimizing a polynomial, because roots that don’t come with a sign change aren’t local extrema.

fn roots_between_with_buffer<const M: usize>( self, lower: f64, upper: f64, x_error: f64, out: &mut ArrayVec<f64, M>, scratch: &mut ArrayVec<f64, M>, )

impl Poly<3>

pub fn roots(&self) -> ArrayVec<f64, 2>

Returns the roots of this quadratic, in increasing order.

Double-roots are only counted once.

pub fn positive_discriminant_roots(&self) -> Option<(f64, f64)>

Returns the two distinct roots of this quadratic, but only if the discriminant is positive.

fn positive_discriminant_roots_scaled(&self) -> Option<(f64, f64)>

Trait Implementations

impl<const N: usize> Add<&Poly<N>> for Poly<N>

type Output = Poly<N>

The resulting type after applying the + operator.

fn add(self, rhs: &Poly<N>) -> Poly<N>

Performs the + operation.

impl<const N: usize> Add<Poly<N>> for &Poly<N>

type Output = Poly<N>

The resulting type after applying the + operator.

fn add(self, rhs: Poly<N>) -> Poly<N>

Performs the + operation.

impl<const N: usize> Add for Poly<N>

type Output = Poly<N>

The resulting type after applying the + operator.

fn add(self, rhs: Poly<N>) -> Poly<N>

Performs the + operation.

impl<const N: usize> AddAssign<&Poly<N>> for Poly<N>

fn add_assign(&mut self, rhs: &Poly<N>)

Performs the += operation.

impl<const N: usize> AddAssign for Poly<N>

fn add_assign(&mut self, rhs: Poly<N>)

Performs the += operation.

impl<const N: usize> Clone for Poly<N>

fn clone(&self) -> Poly<N>

Returns a duplicate of the value.

fn clone_from(&mut self, source: &Self)

Performs copy-assignment from source.

impl<const N: usize> Debug for Poly<N>

fn fmt(&self, f: &mut Formatter<'_>) -> Result

Formats the value using the given formatter.

impl<const N: usize> Div<f64> for &Poly<N>

type Output = Poly<N>

The resulting type after applying the / operator.

fn div(self, scale: f64) -> Poly<N>

Performs the / operation.

impl<const N: usize> Div<f64> for Poly<N>

type Output = Poly<N>

The resulting type after applying the / operator.

fn div(self, scale: f64) -> Poly<N>

Performs the / operation.

impl<const N: usize> DivAssign<f64> for Poly<N>

fn div_assign(&mut self, scale: f64)

Performs the /= operation.

impl<const N: usize> Mul<f64> for &Poly<N>

type Output = Poly<N>

The resulting type after applying the * operator.

fn mul(self, scale: f64) -> Poly<N>

Performs the * operation.

impl<const N: usize> Mul<f64> for Poly<N>

type Output = Poly<N>

The resulting type after applying the * operator.

fn mul(self, scale: f64) -> Poly<N>

Performs the * operation.

impl<const N: usize> MulAssign<f64> for Poly<N>

fn mul_assign(&mut self, scale: f64)

Performs the *= operation.

impl<const N: usize> PartialEq for Poly<N>

fn eq(&self, other: &Poly<N>) -> bool

Tests for self and other values to be equal, and is used by ==.

fn ne(&self, other: &Rhs) -> bool

Tests for !=. The default implementation is almost always sufficient, and should not be overridden without very good reason.

impl<const N: usize> PartialOrd for Poly<N>

fn partial_cmp(&self, other: &Poly<N>) -> Option<Ordering>

This method returns an ordering between self and other values if one exists.

fn lt(&self, other: &Rhs) -> bool

Tests less than (for self and other) and is used by the < operator.

fn le(&self, other: &Rhs) -> bool

Tests less than or equal to (for self and other) and is used by the <= operator.

fn gt(&self, other: &Rhs) -> bool

Tests greater than (for self and other) and is used by the > operator.

fn ge(&self, other: &Rhs) -> bool

Tests greater than or equal to (for self and other) and is used by the >= operator.

impl<const N: usize> Sub<&Poly<N>> for Poly<N>

type Output = Poly<N>

The resulting type after applying the - operator.

fn sub(self, rhs: &Poly<N>) -> Poly<N>

Performs the - operation.

impl<const N: usize> Sub<Poly<N>> for &Poly<N>

type Output = Poly<N>

The resulting type after applying the - operator.

fn sub(self, rhs: Poly<N>) -> Poly<N>

Performs the - operation.

impl<const N: usize> Sub for Poly<N>

type Output = Poly<N>

The resulting type after applying the - operator.

fn sub(self, rhs: Poly<N>) -> Poly<N>

Performs the - operation.

impl<const N: usize> SubAssign<&Poly<N>> for Poly<N>

fn sub_assign(&mut self, rhs: &Poly<N>)

Performs the -= operation.

impl<const N: usize> SubAssign for Poly<N>

fn sub_assign(&mut self, rhs: Poly<N>)

Performs the -= operation.

impl<const N: usize> Copy for Poly<N>

impl<const N: usize> StructuralPartialEq for Poly<N>