项目作者: nschloe

项目描述 :
Orthogonal polynomials in all shapes and sizes.
高级语言: Python
项目地址: git://github.com/nschloe/orthopy.git
创建时间: 2017-08-20T13:59:37Z
项目社区:https://github.com/nschloe/orthopy

开源协议:GNU General Public License v3.0

下载



orthopy

All about orthogonal polynomials.


PyPi Version
PyPI pyversions
GitHub stars
Downloads

Discord
orthogonal

orthopy provides various orthogonal polynomial classes for
lines,
triangles,
disks,
spheres,
n-cubes,
the nD space with weight function exp(-r2)
and more.
All computations are done using numerically stable recurrence schemes. Furthermore, all
functions are fully vectorized and can return results in exact arithmetic.

Installation

Install orthopy from PyPI with

  1. pip install orthopy

How to get a license

Licenses for personal and academic use can be purchased
here.
You’ll receive a confirmation email with a license key.
Install the key with

  1. plm add <your-license-key>

on your machine and you’re good to go.

For commercial use, please contact support@mondaytech.com.

Basic usage

The main function of all submodules is the iterator Eval which evaluates the series of
orthogonal polynomials with increasing degree at given points using a recurrence
relation, e.g.,

  1. import orthopy
  2. x = 0.5
  3. evaluator = orthopy.c1.legendre.Eval(x, "classical")
  4. for _ in range(5):
  5. print(next(evaluator))
  1. 1.0 # P_0(0.5)
  2. 0.5 # P_1(0.5)
  3. -0.125 # P_2(0.5)
  4. -0.4375 # P_3(0.5)
  5. -0.2890625 # P_4(0.5)

Other ways of getting the first n items are

  1. evaluator = Eval(x, "normal")
  2. vals = [next(evaluator) for _ in range(n)]
  3. import itertools
  4. vals = list(itertools.islice(Eval(x, "normal"), n))

Instead of evaluating at only one point, you can provide any array for x; the
polynomials will then be evaluated for all points at once. You can also use sympy for
symbolic computation:

  1. import itertools
  2. import orthopy
  3. import sympy
  4. x = sympy.Symbol("x")
  5. evaluator = orthopy.c1.legendre.Eval(x, "classical")
  6. for val in itertools.islice(evaluator, 5):
  7. print(sympy.expand(val))
  1. 1
  2. x
  3. 3*x**2/2 - 1/2
  4. 5*x**3/2 - 3*x/2
  5. 35*x**4/8 - 15*x**2/4 + 3/8

All Eval methods have a scaling argument which can have three values:

  • "monic": The leading coefficient is 1.
  • "classical": The maximum value is 1 (or (n+alpha over n)).
  • "normal": The integral of the squared function over the domain is 1.

For univariate (“one-dimensional”) integrals, every new iteration contains one function.
For bivariate (“two-dimensional”) domains, every level will contain one function more
than the previous, and similarly for multivariate families. See the tree plots below.

Line segment (-1, +1) with weight function (1-x)α (1+x)β

Legendre Chebyshev 1 Chebyshev 2

Jacobi, Gegenbauer (α=β), Chebyshev 1 (α=β=-1/2), Chebyshev 2 (α=β=1/2), Legendre
(α=β=0) polynomials.

  1. import orthopy
  2. orthopy.c1.legendre.Eval(x, "normal")
  3. orthopy.c1.chebyshev1.Eval(x, "normal")
  4. orthopy.c1.chebyshev2.Eval(x, "normal")
  5. orthopy.c1.gegenbauer.Eval(x, "normal", lmbda)
  6. orthopy.c1.jacobi.Eval(x, "normal", alpha, beta)

The plots above are generated with

  1. import orthopy
  2. orthopy.c1.jacobi.show(5, "normal", 0.0, 0.0)
  3. # plot, savefig also exist

Recurrence coefficients can be explicitly retrieved by

  1. import orthopy
  2. rc = orthopy.c1.jacobi.RecurrenceCoefficients(
  3. "monic", # or "classical", "normal"
  4. alpha=0, beta=0, symbolic=True
  5. )
  6. print(rc.p0)
  7. for k in range(5):
  8. print(rc[k])
  1. 1
  2. (1, 0, None)
  3. (1, 0, 1/3)
  4. (1, 0, 4/15)
  5. (1, 0, 9/35)
  6. (1, 0, 16/63)

1D half-space with weight function xα exp(-r)

(Generalized) Laguerre polynomials.

  1. evaluator = orthopy.e1r.Eval(x, alpha=0, scaling="normal")

1D space with weight function exp(-r2)

Hermite polynomials come in two standardizations:

  • "physicists" (against the weight function exp(-x ** 2)
  • "probabilists" (against the weight function 1 / sqrt(2 * pi) * exp(-x ** 2 / 2)
  1. evaluator = orthopy.e1r2.Eval(
  2. x,
  3. "probabilists", # or "physicists"
  4. "normal"
  5. )

Associated Legendre “polynomials”

Not all of those are polynomials, so they should really be called associated Legendre
functions. The kth iteration contains 2k+1 functions, indexed from
-k to k. (See the color grouping in the above plot.)

  1. evaluator = orthopy.c1.associated_legendre.Eval(
  2. x, phi=None, standardization="natural", with_condon_shortley_phase=True
  3. )

Triangle (T2)

orthopy’s triangle orthogonal polynomials are evaluated in terms of barycentric
coordinates
, so the
X.shape[0] has to be 3.

  1. import orthopy
  2. bary = [0.1, 0.7, 0.2]
  3. evaluator = orthopy.t2.Eval(bary, "normal")

Disk (S2)

Xu Zernike Zernike 2

orthopy contains several families of orthogonal polynomials on the unit disk: After
Xu,
Zernike, and a simplified version
of Zernike polynomials.

  1. import orthopy
  2. x = [0.1, -0.3]
  3. evaluator = orthopy.s2.xu.Eval(x, "normal")
  4. # evaluator = orthopy.s2.zernike.Eval(x, "normal")
  5. # evaluator = orthopy.s2.zernike2.Eval(x, "normal")

Sphere (U3)

Complex-valued spherical harmonics, (black=zero, green=real positive,
pink=real negative, blue=imaginary positive, yellow=imaginary negative). The
functions in the middle are real-valued. The complex angle takes n turns on
the nth level.

  1. evaluator = orthopy.u3.EvalCartesian(
  2. x,
  3. scaling="quantum mechanic" # or "acoustic", "geodetic", "schmidt"
  4. )
  5. evaluator = orthopy.u3.EvalSpherical(
  6. theta_phi, # polar, azimuthal angles
  7. scaling="quantum mechanic" # or "acoustic", "geodetic", "schmidt"
  8. )

n-Cube (Cn)

C1 (Legendre) C2 C3

Jacobi product polynomials.
All polynomials are normalized on the n-dimensional cube. The dimensionality is
determined by X.shape[0].

  1. evaluator = orthopy.cn.Eval(X, alpha=0, beta=0)
  2. values, degrees = next(evaluator)

nD space with weight function exp(-r2) (Enr2)

E1r2 E2r2 E3r2

Hermite product polynomials.
All polynomials are normalized over the measure. The dimensionality is determined by
X.shape[0].

  1. evaluator = orthopy.enr2.Eval(
  2. x,
  3. standardization="probabilists" # or "physicists"
  4. )
  5. values, degrees = next(evaluator)

Other tools

Relevant publications