6 files changed, 66 insertions, 24 deletions
diff --git a/numpy/polynomial/chebyshev.py b/numpy/polynomial/chebyshev.py
index 022dccfa7..ab4c96a54 100644
--- a/numpy/polynomial/chebyshev.py
+++ b/numpy/polynomial/chebyshev.py
@@ -540,10 +540,17 @@ def chebfromroots(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 = chebsub(chebmulx(prd), r*prd)
-        return prd
+        roots.sort()
+        n = len(roots)
+        p = [chebline(-r, 1) for r in roots]
+        while n > 1:
+            m = n//2
+            tmp = [chebmul(p[i], p[i+m]) for i in range(m)]
+            if n%2:
+                tmp[0] = chebmul(tmp[0], p[-1])
+            p = tmp
+            n = m
+        return p[0]
 
 
 def chebadd(c1, c2):
diff --git a/numpy/polynomial/hermite.py b/numpy/polynomial/hermite.py
index 627232d7c..32b50a866 100644
--- a/numpy/polynomial/hermite.py
+++ b/numpy/polynomial/hermite.py
@@ -287,10 +287,17 @@ def hermfromroots(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 = hermsub(hermmulx(prd), r*prd)
-        return prd
+        roots.sort()
+        n = len(roots)
+        p = [hermline(-r, 1) for r in roots]
+        while n > 1:
+            m = n//2
+            tmp = [hermmul(p[i], p[i+m]) for i in range(m)]
+            if n%2:
+                tmp[0] = hermmul(tmp[0], p[-1])
+            p = tmp
+            n = m
+        return p[0]
 
 
 def hermadd(c1, c2):
diff --git a/numpy/polynomial/hermite_e.py b/numpy/polynomial/hermite_e.py
index daed05008..0511307a4 100644
--- a/numpy/polynomial/hermite_e.py
+++ b/numpy/polynomial/hermite_e.py
@@ -287,10 +287,17 @@ def hermefromroots(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 = hermesub(hermemulx(prd), r*prd)
-        return prd
+        roots.sort()
+        n = len(roots)
+        p = [hermeline(-r, 1) for r in roots]
+        while n > 1:
+            m = n//2
+            tmp = [hermemul(p[i], p[i+m]) for i in range(m)]
+            if n%2:
+                tmp[0] = hermemul(tmp[0], p[-1])
+            p = tmp
+            n = m
+        return p[0]
 
 
 def hermeadd(c1, c2):
diff --git a/numpy/polynomial/laguerre.py b/numpy/polynomial/laguerre.py
index fc9afec00..866aab5c4 100644
--- a/numpy/polynomial/laguerre.py
+++ b/numpy/polynomial/laguerre.py
@@ -283,10 +283,17 @@ def lagfromroots(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 = lagsub(lagmulx(prd), r*prd)
-        return prd
+        roots.sort()
+        n = len(roots)
+        p = [lagline(-r, 1) for r in roots]
+        while n > 1:
+            m = n//2
+            tmp = [lagmul(p[i], p[i+m]) for i in range(m)]
+            if n%2:
+                tmp[0] = lagmul(tmp[0], p[-1])
+            p = tmp
+            n = m
+        return p[0]
 
 
 def lagadd(c1, c2):
diff --git a/numpy/polynomial/legendre.py b/numpy/polynomial/legendre.py
index 8fd64985b..b76e678c2 100644
--- a/numpy/polynomial/legendre.py
+++ b/numpy/polynomial/legendre.py
@@ -314,10 +314,17 @@ def legfromroots(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 = legsub(legmulx(prd), r*prd)
-        return prd
+        roots.sort()
+        n = len(roots)
+        p = [legline(-r, 1) for r in roots]
+        while n > 1:
+            m = n//2
+            tmp = [legmul(p[i], p[i+m]) for i in range(m)]
+            if n%2:
+                tmp[0] = legmul(tmp[0], p[-1])
+            p = tmp
+            n = m
+        return p[0]
 
 
 def legadd(c1, c2):
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):