Orthogonal Polynomials¶
This module wraps some of the orthogonal/special functions in the Maxima package “orthopoly”. This package was written by Barton Willis of the University of Nebraska at Kearney. It is released under the terms of the General Public License (GPL). Send Maxima-related bug reports and comments on this module to willisb@unk.edu. In your report, please include Maxima and specfun version information.
The Chebyshev polynomial of the first kind arises as a solution to the differential equation
and those of the second kind as a solution to
The Chebyshev polynomials of the first kind are defined by the recurrence relation
The Chebyshev polynomials of the second kind are defined by the recurrence relation
For integers
, they satisfy the orthogonality relations
and
They are named after Pafnuty Chebyshev (alternative transliterations: Tchebyshef or Tschebyscheff).
The Hermite polynomials are defined either by
(the “probabilists’ Hermite polynomials”), or by
(the “physicists’ Hermite polynomials”). Sage (via Maxima) implements the latter flavor. These satisfy the orthogonality relation
They are named in honor of Charles Hermite.
Each Legendre polynomial
is an
-th degree polynomial. It may be expressed using Rodrigues’ formula:
These are solutions to Legendre’s differential equation:
and satisfy the orthogonality relation
The Legendre function of the second kind
is another (linearly independent) solution to the Legendre differential equation. It is not an “orthogonal polynomial” however.
The associated Legendre functions of the first kind
can be given in terms of the “usual” Legendre polynomials by
Assuming
, they satisfy the orthogonality relation:
where
is the Kronecker delta.
The associated Legendre functions of the second kind
can be given in terms of the “usual” Legendre polynomials by
They are named after Adrien-Marie Legendre.
Laguerre polynomials may be defined by the Rodrigues formula
They are solutions of Laguerre’s equation:
and satisfy the orthogonality relation
The generalized Laguerre polynomials may be defined by the Rodrigues formula:
(These are also sometimes called the associated Laguerre polynomials.) The simple Laguerre polynomials are recovered from the generalized polynomials by setting
.
They are named after Edmond Laguerre.
Jacobi polynomials are a class of orthogonal polynomials. They are obtained from hypergeometric series in cases where the series is in fact finite:
where
is Pochhammer’s symbol (for the rising factorial), (Abramowitz and Stegun p561.) and thus have the explicit expression
They are named after Carl Jacobi.
Ultraspherical or Gegenbauer polynomials are given in terms of the Jacobi polynomials
with
by
They satisfy the orthogonality relation
for
. They are obtained from hypergeometric series in cases where the series is in fact finite:
where
is the falling factorial. (See Abramowitz and Stegun p561)
They are named for Leopold Gegenbauer (1849-1903).
For completeness, the Pochhammer symbol, introduced by Leo August
Pochhammer, , is used in the theory of special
functions to represent the “rising factorial” or “upper factorial”
On the other hand, the “falling factorial” or “lower factorial” is
in the notation of Ronald L. Graham, Donald E. Knuth and Oren Patashnik in their book Concrete Mathematics.
Note
The first call of any of these will usually cost a bit extra (it loads “specfun”, but I’m not sure if that is the real reason). The next call is usually faster but not always.
Todo
Implement associated Legendre polynomials and Zernike polynomials. (Neither is in Maxima.) Wikipedia article Associated_Legendre_polynomials Wikipedia article Zernike_polynomials
REFERENCES:
[ASHandbook] | (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11) Abramowitz and Stegun: Handbook of Mathematical Functions, http://www.math.sfu.ca/ cbm/aands/ |
[EffCheby] | (1, 2) Wolfram Koepf: Effcient Computation of Chebyshev Polynomials in Computer Algebra Computer Algebra Systems: A Practical Guide. John Wiley, Chichester (1999): 79-99. |
AUTHORS:
- David Joyner (2006-06)
- Stefan Reiterer (2010-)
-
class
sage.functions.orthogonal_polys.
ChebyshevPolynomial
(name, nargs=2, latex_name=None, conversions={})¶ Bases:
sage.functions.orthogonal_polys.OrthogonalPolynomial
Abstract base class for Chebyshev polynomials of the first and second kind.
EXAMPLES:
sage: chebyshev_T(3,x) 4*x^3 - 3*x
-
class
sage.functions.orthogonal_polys.
Func_chebyshev_T
¶ Bases:
sage.functions.orthogonal_polys.ChebyshevPolynomial
Chebyshev polynomials of the first kind.
REFERENCE:
- [ASHandbook] 22.5.31 page 778 and 6.1.22 page 256.
EXAMPLES:
sage: chebyshev_T(5,x) 16*x^5 - 20*x^3 + 5*x sage: var('k') k sage: test = chebyshev_T(k,x) sage: test chebyshev_T(k, x)
-
eval_algebraic
(n, x)¶ Evaluate
chebyshev_T
as polynomial, using a recursive formula.INPUT:
n
– an integerx
– a value to evaluate the polynomial at (this can be any ring element)
EXAMPLES:
sage: chebyshev_T.eval_algebraic(5, x) 2*(2*(2*x^2 - 1)*x - x)*(2*x^2 - 1) - x sage: chebyshev_T(-7, x) - chebyshev_T(7,x) 0 sage: R.<t> = ZZ[] sage: chebyshev_T.eval_algebraic(-1, t) t sage: chebyshev_T.eval_algebraic(0, t) 1 sage: chebyshev_T.eval_algebraic(1, t) t sage: chebyshev_T(7^100, 1/2) 1/2 sage: chebyshev_T(7^100, Mod(2,3)) 2 sage: n = 97; x = RIF(pi/2/n) sage: chebyshev_T(n, cos(x)).contains_zero() True sage: R.<t> = Zp(2, 8, 'capped-abs')[] sage: chebyshev_T(10^6+1, t) (2^7 + O(2^8))*t^5 + (O(2^8))*t^4 + (2^6 + O(2^8))*t^3 + (O(2^8))*t^2 + (1 + 2^6 + O(2^8))*t + (O(2^8))
-
eval_formula
(n, x)¶ Evaluate
chebyshev_T
using an explicit formula. See [ASHandbook] 227 (p. 782) for details for the recurions. See also [EffCheby] for fast evaluation techniques.INPUT:
n
– an integerx
– a value to evaluate the polynomial at (this can be any ring element)
EXAMPLES:
sage: chebyshev_T.eval_formula(-1,x) x sage: chebyshev_T.eval_formula(0,x) 1 sage: chebyshev_T.eval_formula(1,x) x sage: chebyshev_T.eval_formula(2,0.1) == chebyshev_T._evalf_(2,0.1) True sage: chebyshev_T.eval_formula(10,x) 512*x^10 - 1280*x^8 + 1120*x^6 - 400*x^4 + 50*x^2 - 1 sage: chebyshev_T.eval_algebraic(10,x).expand() 512*x^10 - 1280*x^8 + 1120*x^6 - 400*x^4 + 50*x^2 - 1
-
class
sage.functions.orthogonal_polys.
Func_chebyshev_U
¶ Bases:
sage.functions.orthogonal_polys.ChebyshevPolynomial
Class for the Chebyshev polynomial of the second kind.
REFERENCE:
- [ASHandbook] 22.8.3 page 783 and 6.1.22 page 256.
EXAMPLES:
sage: R.<t> = QQ[] sage: chebyshev_U(2,t) 4*t^2 - 1 sage: chebyshev_U(3,t) 8*t^3 - 4*t
-
eval_algebraic
(n, x)¶ Evaluate
chebyshev_U
as polynomial, using a recursive formula.INPUT:
n
– an integerx
– a value to evaluate the polynomial at (this can be any ring element)
EXAMPLES:
sage: chebyshev_U.eval_algebraic(5,x) -2*((2*x + 1)*(2*x - 1)*x - 4*(2*x^2 - 1)*x)*(2*x + 1)*(2*x - 1) sage: parent(chebyshev_U(3, Mod(8,9))) Ring of integers modulo 9 sage: parent(chebyshev_U(3, Mod(1,9))) Ring of integers modulo 9 sage: chebyshev_U(-3,x) + chebyshev_U(1,x) 0 sage: chebyshev_U(-1,Mod(5,8)) 0 sage: parent(chebyshev_U(-1,Mod(5,8))) Ring of integers modulo 8 sage: R.<t> = ZZ[] sage: chebyshev_U.eval_algebraic(-2, t) -1 sage: chebyshev_U.eval_algebraic(-1, t) 0 sage: chebyshev_U.eval_algebraic(0, t) 1 sage: chebyshev_U.eval_algebraic(1, t) 2*t sage: n = 97; x = RIF(pi/n) sage: chebyshev_U(n-1, cos(x)).contains_zero() True sage: R.<t> = Zp(2, 6, 'capped-abs')[] sage: chebyshev_U(10^6+1, t) (2 + O(2^6))*t + (O(2^6))
-
eval_formula
(n, x)¶ Evaluate
chebyshev_U
using an explicit formula. See [ASHandbook] 227 (p. 782) for details on the recurions. See also [EffCheby] for the recursion formulas.INPUT:
n
– an integerx
– a value to evaluate the polynomial at (this can be any ring element)
EXAMPLES:
sage: chebyshev_U.eval_formula(10, x) 1024*x^10 - 2304*x^8 + 1792*x^6 - 560*x^4 + 60*x^2 - 1 sage: chebyshev_U.eval_formula(-2, x) -1 sage: chebyshev_U.eval_formula(-1, x) 0 sage: chebyshev_U.eval_formula(0, x) 1 sage: chebyshev_U.eval_formula(1, x) 2*x sage: chebyshev_U.eval_formula(2,0.1) == chebyshev_U._evalf_(2,0.1) True
-
class
sage.functions.orthogonal_polys.
OrthogonalPolynomial
(name, nargs=2, latex_name=None, conversions={})¶ Bases:
sage.symbolic.function.BuiltinFunction
Base class for orthogonal polynomials.
This class is an abstract base class for all orthogonal polynomials since they share similar properties. The evaluation as a polynomial is either done via maxima, or with pynac.
Convention: The first argument is always the order of the polynomial, the others are other values or parameters where the polynomial is evaluated.
-
eval_formula
(*args)¶ Evaluate this polynomial using an explicit formula.
EXAMPLES:
sage: from sage.functions.orthogonal_polys import OrthogonalPolynomial sage: P = OrthogonalPolynomial('testo_P') sage: P.eval_formula(1,2.0) Traceback (most recent call last): ... NotImplementedError: no explicit calculation of values implemented
-
-
sage.functions.orthogonal_polys.
gegenbauer
(n, a, x)¶ Returns the ultraspherical (or Gegenbauer) polynomial for integers
.
Computed using Maxima.
REFERENCE:
- [ASHandbook] 22.5.27
EXAMPLES:
sage: x = PolynomialRing(QQ, 'x').gen() sage: ultraspherical(2,3/2,x) 15/2*x^2 - 3/2 sage: ultraspherical(2,1/2,x) 3/2*x^2 - 1/2 sage: ultraspherical(1,1,x) 2*x sage: t = PolynomialRing(RationalField(),"t").gen() sage: gegenbauer(3,2,t) 32*t^3 - 12*t
- Check that trac ticket #17192 is fixed::
sage: x = PolynomialRing(QQ, ‘x’).gen() sage: ultraspherical(0,1,x) 1
sage: ultraspherical(-1,1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -1
sage: ultraspherical(-7,1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -7
-
sage.functions.orthogonal_polys.
gen_laguerre
(n, a, x)¶ Returns the generalized Laguerre polynomial for integers
. Typically,
or
.
REFERENCES:
- Table on page 789 in [ASHandbook].
EXAMPLES:
sage: x = PolynomialRing(QQ, 'x').gen() sage: gen_laguerre(2,1,x) 1/2*x^2 - 3*x + 3 sage: gen_laguerre(2,1/2,x) 1/2*x^2 - 5/2*x + 15/8 sage: gen_laguerre(2,-1/2,x) 1/2*x^2 - 3/2*x + 3/8 sage: gen_laguerre(2,0,x) 1/2*x^2 - 2*x + 1 sage: gen_laguerre(3,0,x) -1/6*x^3 + 3/2*x^2 - 3*x + 1
- Check that trac ticket #17192 is fixed::
sage: x = PolynomialRing(QQ, ‘x’).gen() sage: gen_laguerre(0,1,x) 1
sage: gen_laguerre(-1,1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -1
sage: gen_laguerre(-7,1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -7
-
sage.functions.orthogonal_polys.
gen_legendre_P
(n, m, x)¶ Returns the generalized (or associated) Legendre function of the first kind.
The awkward code for when m is odd and 1 results from the fact that Maxima is happy with, for example,
, but Sage is not. For these cases the function is computed from the (m-1)-case using one of the recursions satisfied by the Legendre functions.
REFERENCE:
- Gradshteyn and Ryzhik 8.706 page 1000.
EXAMPLES:
sage: P.<t> = QQ[] sage: gen_legendre_P(2, 0, t) 3/2*t^2 - 1/2 sage: gen_legendre_P(2, 0, t) == legendre_P(2, t) True sage: gen_legendre_P(3, 1, t) -3/2*(5*t^2 - 1)*sqrt(-t^2 + 1) sage: gen_legendre_P(4, 3, t) 105*(t^2 - 1)*sqrt(-t^2 + 1)*t sage: gen_legendre_P(7, 3, I).expand() -16695*sqrt(2) sage: gen_legendre_P(4, 1, 2.5) -583.562373654533*I
-
sage.functions.orthogonal_polys.
gen_legendre_Q
(n, m, x)¶ Returns the generalized (or associated) Legendre function of the second kind.
Maxima restricts m = n. Hence the cases m n are computed using the same recursion used for gen_legendre_P(n,m,x) when m is odd and 1.
EXAMPLES:
sage: P.<t> = QQ[] sage: gen_legendre_Q(2,0,t) 3/4*t^2*log(-(t + 1)/(t - 1)) - 3/2*t - 1/4*log(-(t + 1)/(t - 1)) sage: gen_legendre_Q(2,0,t) - legendre_Q(2, t) 0 sage: gen_legendre_Q(3,1,0.5) 2.49185259170895 sage: gen_legendre_Q(0, 1, x) -1/sqrt(-x^2 + 1) sage: gen_legendre_Q(2, 4, x).factor() 48*x/((x + 1)^2*(x - 1)^2)
-
sage.functions.orthogonal_polys.
hermite
(n, x)¶ Returns the Hermite polynomial for integers
.
REFERENCE:
- [ASHandbook] 22.5.40 and 22.5.41, page 779.
EXAMPLES:
sage: x = PolynomialRing(QQ, 'x').gen() sage: hermite(2,x) 4*x^2 - 2 sage: hermite(3,x) 8*x^3 - 12*x sage: hermite(3,2) 40 sage: S.<y> = PolynomialRing(RR) sage: hermite(3,y) 8.00000000000000*y^3 - 12.0000000000000*y sage: R.<x,y> = QQ[] sage: hermite(3,y^2) 8*y^6 - 12*y^2 sage: w = var('w') sage: hermite(3,2*w) 8*(8*w^2 - 3)*w
- Check that trac ticket #17192 is fixed::
sage: x = PolynomialRing(QQ, ‘x’).gen() sage: hermite(0,x) 1
sage: hermite(-1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -1
sage: hermite(-7,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -7
-
sage.functions.orthogonal_polys.
jacobi_P
(n, a, b, x)¶ Returns the Jacobi polynomial
for integers
and a and b symbolic or
and
. The Jacobi polynomials are actually defined for all a and b. However, the Jacobi polynomial weight
isn’t integrable for
or
.
REFERENCE:
- Table on page 789 in [ASHandbook].
EXAMPLES:
sage: x = PolynomialRing(QQ, 'x').gen() sage: jacobi_P(2,0,0,x) 3/2*x^2 - 1/2 sage: jacobi_P(2,1,2,1.2) # random output of low order bits 5.009999999999998
- Check that trac ticket #17192 is fixed::
sage: x = PolynomialRing(QQ, ‘x’).gen() sage: jacobi_P(0,0,0,x) 1
sage: jacobi_P(-1,0,0,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -1
sage: jacobi_P(-7,0,0,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -7
-
sage.functions.orthogonal_polys.
laguerre
(n, x)¶ Return the Laguerre polynomial for integers
.
REFERENCE:
- [ASHandbook] 22.5.16, page 778 and page 789.
EXAMPLES:
sage: x = PolynomialRing(QQ, 'x').gen() sage: laguerre(2,x) 1/2*x^2 - 2*x + 1 sage: laguerre(3,x) -1/6*x^3 + 3/2*x^2 - 3*x + 1 sage: laguerre(2,2) -1
- Check that trac ticket #17192 is fixed::
sage: x = PolynomialRing(QQ, ‘x’).gen() sage: laguerre(0,x) 1
sage: laguerre(-1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -1
sage: laguerre(-7,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -7
-
sage.functions.orthogonal_polys.
legendre_P
(n, x)¶ Returns the Legendre polynomial of the first kind.
REFERENCE:
- [ASHandbook] 22.5.35 page 779.
EXAMPLES:
sage: P.<t> = QQ[] sage: legendre_P(2,t) 3/2*t^2 - 1/2 sage: legendre_P(3, 1.1) 1.67750000000000 sage: legendre_P(3, 1 + I) 7/2*I - 13/2 sage: legendre_P(3, MatrixSpace(ZZ, 2)([1, 2, -4, 7])) [-179 242] [-484 547] sage: legendre_P(3, GF(11)(5)) 8
-
sage.functions.orthogonal_polys.
legendre_Q
(n, x)¶ Returns the Legendre function of the second kind.
Computed using Maxima.
EXAMPLES:
sage: P.<t> = QQ[] sage: legendre_Q(2, t) 3/4*t^2*log(-(t + 1)/(t - 1)) - 3/2*t - 1/4*log(-(t + 1)/(t - 1)) sage: legendre_Q(3, 0.5) -0.198654771479482 sage: legendre_Q(4, 2) 443/16*I*pi + 443/16*log(3) - 365/12 sage: legendre_Q(4, 2.0) 0.00116107583162324 + 86.9828465962674*I
-
sage.functions.orthogonal_polys.
ultraspherical
(n, a, x)¶ Returns the ultraspherical (or Gegenbauer) polynomial for integers
.
Computed using Maxima.
REFERENCE:
- [ASHandbook] 22.5.27
EXAMPLES:
sage: x = PolynomialRing(QQ, 'x').gen() sage: ultraspherical(2,3/2,x) 15/2*x^2 - 3/2 sage: ultraspherical(2,1/2,x) 3/2*x^2 - 1/2 sage: ultraspherical(1,1,x) 2*x sage: t = PolynomialRing(RationalField(),"t").gen() sage: gegenbauer(3,2,t) 32*t^3 - 12*t
- Check that trac ticket #17192 is fixed::
sage: x = PolynomialRing(QQ, ‘x’).gen() sage: ultraspherical(0,1,x) 1
sage: ultraspherical(-1,1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -1
sage: ultraspherical(-7,1,x) Traceback (most recent call last): ... ValueError: n must be greater than -1, got n = -7