X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Finterpolation.py;h=a5ad78ddb32589137b72658afa4fe6401f17b202;hb=6e68bda35776924bee44e934b862540543335731;hp=cbba68e03a4e2bd3ae5601ec51844aada3cf435d;hpb=e15316c191095ad432b6e28148bee01df6ffebc5;p=sage.d.git diff --git a/mjo/interpolation.py b/mjo/interpolation.py index cbba68e..a5ad78d 100644 --- a/mjo/interpolation.py +++ b/mjo/interpolation.py @@ -18,7 +18,7 @@ def lagrange_denominator(k, xs): The product of all xs[j] with j != k. """ - return product([xs[k] - xs[j] for j in range(0, len(xs)) if j != k]) + return product([xs[k] - xs[j] for j in xrange(len(xs)) if j != k]) def lagrange_coefficient(k, x, xs): @@ -44,11 +44,15 @@ def lagrange_coefficient(k, x, xs): A symbolic expression of one variable. + SETUP:: + + sage: from mjo.interpolation import lagrange_coefficient + TESTS:: sage: xs = [ -pi/2, -pi/6, 0, pi/6, pi/2 ] sage: lagrange_coefficient(0, x, xs) - 1/8*(pi - 6*x)*(pi - 2*x)*(pi + 6*x)*x/pi^4 + 1/8*(pi + 6*x)*(pi - 2*x)*(pi - 6*x)*x/pi^4 """ numerator = lagrange_psi(x, xs)/(x - xs[k]) @@ -60,7 +64,7 @@ def lagrange_coefficient(k, x, xs): def lagrange_polynomial(x, xs, ys): """ - Return the Lagrange form of the interpolation polynomial in `x` of + Return the Lagrange form of the interpolating polynomial in `x` at the points (xs[k], ys[k]). INPUT: @@ -75,6 +79,10 @@ def lagrange_polynomial(x, xs, ys): A symbolic expression (polynomial) interpolating each (xs[k], ys[k]). + SETUP:: + + sage: from mjo.interpolation import lagrange_polynomial + TESTS:: sage: xs = [ -pi/2, -pi/6, 0, pi/6, pi/2 ] @@ -88,8 +96,8 @@ def lagrange_polynomial(x, xs, ys): True """ - ls = [ lagrange_coefficient(k, x, xs) for k in range(0, len(xs)) ] - sigma = sum([ ys[k] * ls[k] for k in range(0, len(xs)) ]) + ls = [ lagrange_coefficient(k, x, xs) for k in xrange(len(xs)) ] + sigma = sum([ ys[k] * ls[k] for k in xrange(len(xs)) ]) return sigma @@ -111,6 +119,10 @@ def lagrange_interpolate(f, x, xs): A polynomial in ``x`` which interpolates ``f`` at ``xs``. + SETUP:: + + sage: from mjo.interpolation import lagrange_interpolate + EXAMPLES: We're exact on polynomials of degree `n` if we use `n+1` points:: @@ -135,7 +147,11 @@ def divided_difference_coefficients(xs): Assuming some function `f`, compute the coefficients of the divided difference f[xs[0], ..., xs[n]]. - TESTS: + SETUP:: + + sage: from mjo.interpolation import divided_difference_coefficients + + TESTS:: sage: divided_difference_coefficients([0]) [1] @@ -145,7 +161,7 @@ def divided_difference_coefficients(xs): [1/2/pi^2, -1/pi^2, 1/2/pi^2] """ - coeffs = [ QQ(1)/lagrange_denominator(k, xs) for k in range(0, len(xs)) ] + coeffs = [ QQ(1)/lagrange_denominator(k, xs) for k in xrange(len(xs)) ] return coeffs @@ -166,6 +182,10 @@ def divided_difference(xs, ys): The (possibly symbolic) divided difference function. + SETUP:: + + sage: from mjo.interpolation import divided_difference + TESTS:: sage: xs = [0] @@ -183,7 +203,7 @@ def divided_difference(xs, ys): We try something entirely symbolic:: - sage: f = function('f', x) + sage: f = function('f')(x) sage: divided_difference([x], [f(x=x)]) f(x) sage: x1,x2 = SR.var('x1,x2') @@ -215,7 +235,11 @@ def newton_polynomial(x, xs, ys): A symbolic expression. - TESTS: + SETUP:: + + sage: from mjo.interpolation import lagrange_polynomial, newton_polynomial + + TESTS:: sage: xs = [ -pi/2, -pi/6, 0, pi/6, pi/2 ] sage: ys = map(sin, xs) @@ -229,7 +253,7 @@ def newton_polynomial(x, xs, ys): N = SR(0) - for k in range(0, degree+1): + for k in xrange(degree+1): term = divided_difference(xs[:k+1], ys[:k+1]) term *= lagrange_psi(x, xs[:k]) N += term @@ -304,7 +328,11 @@ def hermite_interpolant(x, xs, ys, y_primes): A symbolic expression. - TESTS: + SETUP:: + + sage: from mjo.interpolation import hermite_interpolant + + TESTS:: sage: xs = [ 0, pi/6, pi/2 ] sage: ys = map(sin, xs) @@ -319,10 +347,10 @@ def hermite_interpolant(x, xs, ys, y_primes): """ s1 = sum([ ys[k] * hermite_coefficient(k, x, xs) - for k in range(0, len(xs)) ]) + for k in xrange(len(xs)) ]) s2 = sum([ y_primes[k] * hermite_deriv_coefficient(k, x, xs) - for k in range(0, len(xs)) ]) + for k in xrange(len(xs)) ]) return (s1 + s2)