From ef767a54fae8977cb15d7e0849a8cb626c4b0f2c Mon Sep 17 00:00:00 2001 From: Charles Harris Date: Mon, 30 Jan 2012 20:10:52 -0700 Subject: ENH: Improve the computation of polynomials from roots. The original method was overly sensitive to roundoff. Of the two approaches considered, gauss integration or binary subdivision of the roots, the latter is more compatible with using other number representations such as mpmath. No method is going to be suitable for large numbers of arbitrary zeros but the current method is a significant improvement. --- numpy/polynomial/polynomial.py | 15 +++++++++++---- 1 file changed, 11 insertions(+), 4 deletions(-) (limited to 'numpy/polynomial/polynomial.py') diff --git a/numpy/polynomial/polynomial.py b/numpy/polynomial/polynomial.py index 207511216..c8894c68e 100644 --- a/numpy/polynomial/polynomial.py +++ b/numpy/polynomial/polynomial.py @@ -186,10 +186,17 @@ def polyfromroots(roots) : return np.ones(1) else : [roots] = pu.as_series([roots], trim=False) - prd = np.array([1], dtype=roots.dtype) - for r in roots: - prd = polysub(polymulx(prd), r*prd) - return prd + roots.sort() + n = len(roots) + p = [polyline(-r, 1) for r in roots] + while n > 1: + m = n//2 + tmp = [polymul(p[i], p[i+m]) for i in range(m)] + if n%2: + tmp[0] = polymul(tmp[0], p[-1]) + p = tmp + n = m + return p[0] def polyadd(c1, c2): -- cgit v1.2.1