Hi all,
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.
For example, on an Intel i7-4710MQ, running with 8 threads this code
solves approximately 1.2e7 quartic polynomial equations per second.
(With OMP_NUM_THREADS=4 the solver gives approximately 8.9e6 solutions
per second, so it seems HyperThreading in this case helps about 30%.)
The algorithms for cubics and quadratics come from Numerical Recipes
(3rd ed.), and the quartic problem is internally reduced to a cubic and
two quadratics, using well-known standard tricks.
Since to my understanding, for solving polynomial equations NumPy only
provides roots(), which works one problem at a time, I'm writing to
inquire whether there would be interest to include this functionality to
NumPy, if I contribute the code (and clean it up for integration)?
I have some previous open source contribution experience. I have
authored the IVTC and Phosphor deinterlacers for VLC, and a modular
postprocessing filter framework for the Panda3D game engine (posted at
the forums on panda3d.org, projected for version 1.10).
Currently the polynomial solver is available in a git repository hosted
by our university, the relevant files being:
https://yousource.it.jyu.fi/jjrandom2/miniprojects/blobs/master/misc/tworods/polysolve2.pyx
https://yousource.it.jyu.fi/jjrandom2/miniprojects/blobs/master/misc/tworods/polysolve.py
https://yousource.it.jyu.fi/jjrandom2/miniprojects/blobs/master/misc/tworods/polysolve2example.py
https://yousource.it.jyu.fi/jjrandom2/miniprojects/blobs/master/misc/tworods/setup.py
https://yousource.it.jyu.fi/jjrandom2/miniprojects/blobs/master/misc/tworods/compile.sh
I'm working on a research grant, which allows me to hold the copyright;
license is BSD.
Polysolve2 is the fast Cython+OpenMP version, while polysolve is an
earlier (slower) experimental version using only NumPy array operations.
The example module contains a usage example, and setup (in the standard
manner) instructs distutils and Cython as to how to build polysolve2
(including setting the relevant math flags and enabling OpenMP).
(The setup script includes also some stuff specific to my solver for the
Ziegler problem; basically, the "tworods" module can be ignored. I
apologize for the inconvenience.)
Thoughts?
Best regards,
-J
-------------------------------------------------
Juha Jeronen, Ph.D.
juha.jero...@jyu.fi
University of Jyväskylä
Department of Mathematical Information Technology
-------------------------------------------------
_______________________________________________
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
https://mail.scipy.org/mailman/listinfo/numpy-discussion
