On Sat, Oct 3, 2009 at 7:12 PM, Charles R Harris
<charlesr.har...@gmail.com>wrote:
> Attached is the first rc of the chebyshev module. The module documentation
> is not yet complete and no doubt the rest of the documentation needs to be
> reviewed. The tests cover basic functionality at this point but need to be
> extended to cover the Chebyshev object. Nevertheless, the module should be
> usable.
>
> Note that the most convenient way to do the least squared fits is with the
> static method Chebyshev.fit, which will return a Chebyshev object that
> contains both the resulting Chebyshev series and its domain.
>
> Some naming questions remain. ISTM that "lstsq" or "leastsq" might be a
> better name than fit. Likewise, I have kept the poly1d names "deriv" and
> "integ", but "der" and "int" might be more appropriate.
>
> Operators behave as expected for +, -, and * but there is no truedivision
> unless both operands can be interpreted as scalars. When division hasn't
> been imported from __future__, the / and // operators are both floordivision
> and % returns the remainder. Divmod behaves as expected.
>
> Any feedback is welcome.
>
>
And an updated test to reflect changes in treatment of leading zeros.
Chuck
from __future__ import division
import numpy as np
import chebyshev as ch
from numpy.testing import *
from decimal import Decimal as Dec
from exceptions import TypeError, ValueError
def trim(x) :
return ch.chebtrim(x, tol=1e-6)
T0 = [ 1]
T1 = [ 1, 0]
T2 = [ 2, 0, -1]
T3 = [ 4, 0, -3, 0]
T4 = [ 8, 0, -8, 0, 1]
T5 = [ 16, 0, -20, 0, 5, 0]
T6 = [ 32, 0, -48, 0, 18, 0, -1]
T7 = [ 64, 0, -112, 0, 56, 0, -7, 0]
T8 = [ 128, 0, -256, 0, 160, 0, -32, 0, 1]
T9 = [ 256, 0, -576, 0, 432, 0, -120, 0, 9, 0.]
Tlist = [T0, T1, T2, T3, T4, T5, T6, T7, T8, T9]
def test_chebadd() :
for i in range(5) :
for j in range(5) :
msg = "At i=%d, j=%d" % (i,j)
tgt = np.zeros(max(i,j) + 1)
tgt[::-1][i] += 1
tgt[::-1][j] += 1
res = ch.chebadd([1]+[0]*i, [1]+[0]*j)
assert_equal(trim(res), trim(tgt), err_msg=msg)
def test_chebsub() :
for i in range(5) :
for j in range(5) :
msg = "At i=%d, j=%d" % (i,j)
tgt = np.zeros(max(i,j) + 1)
tgt[::-1][i] += 1
tgt[::-1][j] -= 1
res = ch.chebsub([1]+[0]*i, [1]+[0]*j)
assert_equal(trim(res), trim(tgt), err_msg=msg)
def test_chebmul() :
for i in range(5) :
for j in range(5) :
msg = "At i=%d, j=%d" % (i,j)
tgt = np.zeros(i + j + 1)
tgt[::-1][i + j] += .5
tgt[::-1][abs(i - j)] += .5
res = ch.chebmul([1]+[0]*i, [1]+[0]*j)
assert_equal(trim(res), trim(tgt), err_msg=msg)
def test_chebdiv() :
for i in range(5) :
for j in range(5) :
msg = "At i=%d, j=%d" % (i,j)
ci = [1] + [0]*i
cj = [1] + [0]*j
tgt = ch.chebadd(ci, cj)
quo, rem = ch.chebdiv(tgt, ci)
res = ch.chebadd(ch.chebmul(quo, ci), rem)
assert_equal(trim(res), trim(tgt), err_msg=msg)
def test_cheb2poly() :
for i in range(10) :
assert_equal(ch.cheb2poly([1] + [0]*i), Tlist[i])
def test_poly2cheb() :
for i in range(10) :
assert_equal(ch.poly2cheb(Tlist[i]), [1] + [0]*i)
def test_chebint() :
# check exceptions
assert_raises(ValueError, ch.chebint, [0], -1)
assert_raises(ValueError, ch.chebint, [0], 1, [0,0])
# check single integration
for i in range(5) :
scl = i + 1
pol = [1] + [0]*i
tgt = [1/scl] + [0]*i + [i]
chebpol = ch.poly2cheb(pol)
chebint = ch.chebint(chebpol, m=1, k=[i])
res = ch.cheb2poly(chebint)
assert_almost_equal(trim(res), trim(tgt))
# check multiple integrations with default k
for i in range(5) :
for j in range(2,5) :
pol = [1] + [0]*i
tgt = pol[:]
for k in range(j) :
tgt = ch.chebint(tgt, m=1)
res = ch.chebint(pol, m=j)
assert_almost_equal(trim(res), trim(tgt))
# check multiple integrations with defined k
for i in range(5) :
for j in range(2,5) :
pol = [1] + [0]*i
tgt = pol[:]
for k in range(j) :
tgt = ch.chebint(tgt, m=1, k=[k])
res = ch.chebint(pol, m=j, k=range(j))
assert_almost_equal(trim(res), trim(tgt))
def test_chebder() :
# check exceptions
assert_raises(ValueError, ch.chebder, [0], -1)
# check that zeroth deriviative does nothing
for i in range(5) :
tgt = [1] + [0]*i
res = ch.chebder(tgt, m=0)
assert_equal(trim(res), trim(tgt))
# check that derivation is the inverse of integration
for i in range(5) :
for j in range(2,5) :
tgt = [1] + [0]*i
res = ch.chebder(ch.chebint(tgt, m=j), m=j)
assert_almost_equal(trim(res), trim(tgt))
def test_chebval() :
#check empty input
assert_equal(ch.chebval([], 1).size, 0)
#check normal input)
for i in range(5) :
coef = [1] + [0]*i
assert_almost_equal(ch.chebval(1, coef), 1)
assert_almost_equal(ch.chebval(-1, coef), (-1)**i)
zpts = np.cos(np.linspace(-np.pi, 0, 2*i + 1)[1::2])
for z in zpts :
assert_almost_equal(ch.chebval(z, coef), 0)
#check that shape is preserved
for i in range(3) :
dims = [2]*i
x = np.zeros(dims)
assert_equal(ch.chebval(x, [1]).shape, dims)
assert_equal(ch.chebval(x, [1,0]).shape, dims)
assert_equal(ch.chebval(x, [1,0,0]).shape, dims)
def test_chebfromroots() :
res = ch.chebfromroots([])
assert_almost_equal(trim(res), [1])
for i in range(1,5) :
roots = np.cos(np.linspace(-np.pi, 0, 2*i + 1)[1::2])
tgt = [1] + [0]*i
res = ch.chebfromroots(roots)*2**(i-1)
assert_almost_equal(trim(res),trim(tgt))
def test_chebroots() :
assert_almost_equal(ch.chebroots([1]), [])
for i in range(1,5) :
tgt = np.linspace(-1, 1, i)
res = ch.chebroots(ch.chebfromroots(tgt))
assert_almost_equal(trim(res), trim(tgt))
def test_chebfit() :
def f(x) :
return x*(x - 1)*(x - 2)
# Test exceptions
assert_raises(ValueError, ch.chebfit, [1], [1], -1)
assert_raises(TypeError, ch.chebfit, [[1]], [1], 0)
assert_raises(TypeError, ch.chebfit, [], [1], 0)
assert_raises(TypeError, ch.chebfit, [1], [[[1]]], 0)
assert_raises(TypeError, ch.chebfit, [1, 2], [1], 0)
assert_raises(TypeError, ch.chebfit, [1], [1, 2], 0)
# Test fit
x = np.linspace(0,2)
y = f(x)
coef, domain = ch.chebfit(x, y, 3)
assert_equal(domain, [0, 2])
assert_equal(len(coef), 4)
assert_almost_equal(ch.chebval(x, coef, domain), y)
coef, domain = ch.chebfit(x, y, 4, domain=[-5, 5])
assert_equal(domain, [-5, 5])
assert_equal(len(coef), 5)
assert_almost_equal(ch.chebval(x, coef, domain), y)
def test_chebtrim() :
coef = [0, 1, -1, 2]
# Test exceptions
assert_raises(ValueError, ch.chebtrim, coef, -1)
# Test results
assert_equal(ch.chebtrim(coef), coef[1:])
assert_equal(ch.chebtrim(coef, 1), coef[3:])
assert_equal(ch.chebtrim(coef, 2), [0])
if __name__ == "__main__" :
test_chebadd()
test_chebsub()
test_chebmul()
test_chebdiv()
test_cheb2poly()
test_poly2cheb()
test_chebint()
test_chebder()
test_chebval()
test_chebfromroots()
test_chebroots()
test_chebfit()
test_chebtrim()
_______________________________________________
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
