[Numpy-discussion] Cython-based OpenMP-accelerated quartic polynomial solver

[Juha Jeronen]

On 30.09.2015 13:21, Pauli Virtanen wrote:
> Juha Jeronen <juha.jeronen <at> jyu.fi> writes:
>> I recently developed a Cython-based, OpenMP-accelerated quartic (and
>> cubic, quadratic) polynomial solver to address a personal research need
>> for quickly solving a very large number of independent low-degree
>> polynomial equations for both real and complex coefficients.
> My 2c in this context would be to also think how this best fits
> in how collections of polynomials are represented in Numpy.
>
> AFAICS, Numpy's polynomial module supports evaluation of polynomial
> collections (cf. numpy.polynomial.polynomial.polyval),
> but the corresponding root finding routine
> (numpy.polynomial.polynomial.polyroots) only deals with one
> polynomial at a time.
>
> The present code in principle could be used to speed up the latter
> routine after it is generalized to multiple polynomials. The general
> case is probably doable just by coding up the companion matrix
> approach using low-level routines (or maybe with the new vectorized
> np.linalg routines).
>
> Some relevant code elsewhere for similar purposes can be
> found in scipy.inteprolate.PPoly/BPoly (and in the future BSpline).
>
> https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.PPoly.html
>
> However, since it's piecewise, there's purposefully support only
> for real-valued roots.

Thanks.

It sounds like numpy.polynomial.polynomial.polyroots() is a good 
candidate place for this, if as you suggest, it is first generalized to 
handle multiple polynomials.

  -J