[SciPy-user] Orthogonal polynomial evaluation

Pauli Virtanen pav at iki.fi
Mon Apr 13 20:21:15 EDT 2009


Mon, 13 Apr 2009 20:09:58 +0000, Pauli Virtanen wrote:
[clip]
> I wonder if simply making orthopoly1d evaluate the polynomials using the
> Newton scheme would fix stability issues... Needs some thinking.

This seems to actually help quite a lot: see the branch at

	http://github.com/pv/scipy-work/tree/ticket/921-orthogonal

For example,

Before (Horner):

>>> import scipy.special as sc
>>> sc.chebyt(300)(0.999) - np.cos(300*np.arccos(0.999))
-2.9807351165313966e+117

After (using roots):

>>> import scipy.special as sc
>>> sc.chebyt(300)(0.999) - np.cos(300*np.arccos(0.999))
1.9129142714291447e-13

Specialized evaluation routines could of course be faster; finding the 
zeros of the polynomial only for evaluating it is overkill...

-- 
Pauli Virtanen




More information about the SciPy-User mailing list