The formula given on the page appears to take an exponential number of multiplications. That seems excessively inefficient.
In a finite field, at least, one could compute 2n syndromes (sums of the roots raised to 1st through 2n-th power) and then use the berlekamp massey to recover the polynomial. Would take O(n^2) operations overall. If it works over the reals I'm not sure how stable it would be numerically.
But I'm sure there has to be an more efficient approach than an exponential one. :)
def coeffs_from_roots(roots):
coefs=[1]
for r in roots:
coefs=[-r*coefs[0]] + [coefs[i]-r*coefs[i+1] for i in range(len(coefs)-1)] + [1]
coefs.reverse()
return coefs
That's merely quadratic, rather than taking exponentially many multiplications.
In a finite field, at least, one could compute 2n syndromes (sums of the roots raised to 1st through 2n-th power) and then use the berlekamp massey to recover the polynomial. Would take O(n^2) operations overall. If it works over the reals I'm not sure how stable it would be numerically.
But I'm sure there has to be an more efficient approach than an exponential one. :)