import numpy as np
import numpy.linalg as la
import matplotlib.pyplot as pt
sympy
¶import sympy as sym
# Enable "pretty-printing" in IPython
sym.init_printing()
Make a new Symbol
and work with it:
x = sym.Symbol("x")
myexpr = (x**2-3)**2
myexpr
myexpr = (x**2-3)**2
myexpr
myexpr.expand()
sym.integrate(myexpr, x)
sym.integrate(myexpr, (x, -1, 1))
Now write a function inner_product(f, g)
:
def inner_product(f, g):
return sym.integrate(f*g, (x, -1, 1))
Show that it works:
inner_product(1, 1)
inner_product(1, x)
Next, define a basis
consisting of a few monomials:
basis = [1, x, x**2, x**3]
#basis = [1, x, x**2, x**3, x**4, x**5]
And run Gram-Schmidt on it:
orth_basis = []
for q in basis:
for prev_q in orth_basis:
q = q - inner_product(prev_q, q)*prev_q
q = q / sym.sqrt(inner_product(q, q))
orth_basis.append(q)
orth_basis
These are called the Legendre polynomials.
What do they look like?
mesh = np.linspace(-1, 1, 100)
pt.figure(figsize=(8,8))
for f in orth_basis:
f = sym.lambdify(x, f)
pt.plot(mesh, [f(xi) for xi in mesh])
These functions are important enough to be included in scipy.special
as eval_legendre
:
!! Careful: The Scipy versions do not have norm 1.
import scipy.special as sps
for i in range(10):
pt.plot(mesh, sps.eval_legendre(i, mesh))
What can we find out about the conditioning of the generalized Vandermonde matrix for Legendre polynomials?
#keep
n = 20
xs = np.linspace(-1, 1, n)
V = np.array([
sps.eval_legendre(i, xs)
for i in range(n)
]).T
la.cond(V)