pyquartic

Modified Ferrari's quartic solver and modified Cardano's cubic solver for Python (4th and 3rd order polynomials).

Features

• Original algorithms are modified so they provide stable and correct solutions to polynomials
• Using modified algorithms from quarticequations.com
• Functions are optimized for fast computing (up to 20x faster than np.roots)
• Cython and pure-Python implementation

Usage

Cubic solver

import pyquartic
roots = pyquartic.solve_cubic(a, b, c, d)

Where a, b, c, d are cubic equation coefficients:
$ax^3 + bx^2 + cx + d = 0$

Quartic solver

import pyquartic
roots = pyquartic.solve_quartic(a, b, c, d, e)

Where a, b, c, d, e are quartic equation coefficients:
$ax^4 + bx^3 + cx^2 + dx + e = 0$

Pure-Python fallback

pyquartic_python.py can be used in case where cython extensions cant be compiled.
Output values are complex numbers in a tuple, where each number is one root.

Speed analysis

Speed analysis can be performed by running time_pyquartic.py.
Numpy's polyroots solver computes eigenvalues of a companion matrix, formed from n-th order polynomial coefficients. This method is stable and accurate, but very slow. Tests are performed on large number of coefficients (10000) where coefficients range from -10000 to 10000, times shown are mean and best times from those iterations.
Provided are three tests for cubic and quartic each: numpy's polyroots, python implementation, and cython implementation.

Cubic
np.root: 61.8645 us, best: 58.66 us
python : 10.3469 us, best: 9.488 us
cython : 3.147 us, best: 2.657 us

Quartic
np.root: 64.6545 us, best: 61.261 us
python : 16.0815 us, best: 14.85 us
cython : 3.4743 us, best: 2.969 us  

Modifications to the algorithms

Much more detailed versions and tutorials for this modifications are covered in this paper: quarticequations.com.
This project is only implementation of there described modified algorithms.

Cardano's method

Cardano's method has large round-off error at some specific cases (when parameter q approaches zero).
While mathematically correct, when calculated on computer, where numbers are rounded, this can create large errors, that are later carried to quartic solver.
To fix this, when there is only one real solution, solution from 'Numerical Recipes' is used.
When there are three real solutions, Viète’s trigonometric method is used.
More details here.

Ferrari's method

Ferrari's method uses one real root from cubic equation solved with Cardano's method.
As this root approaches zero, algorithm becomes computationally unstable due to round-off error.
Modified algorithm is called 'Modern generalization of Cardano’s Problem VIII'. It avoids this instability by adding one more calculation.
More details here.