Figure 2.7:

Function f(x)=H(x) and truncated Legendre-Fourier series approximations with n=10,50,100.

Code for Figure 2.7

Text of the GNU GPL.

main.py


 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import gamma, lpmv

nxs = 301
x = np.linspace(-1, 1, nxs)
nterms = [10, 50, 100]
u = np.zeros((len(nterms), nxs))
ntermsmax = nterms[-1]
phi = np.zeros((ntermsmax+1, nxs))

for m in range(ntermsmax+1):
    P = lpmv(0, m, x)
    phi[m, :] = P * np.sqrt((2 * m + 1) / 2)

flist = (np.sign(x) + 1) / 2

for i in range(len(nterms)):
    k = 0
    un0 = np.sqrt(np.pi) / (2 * gamma(1 - k / 2) * gamma((3 + k) / 2)) * np.sqrt((2 * k + 1) / 2)
    u[i, :] += un0 * phi[k, :]
    for k in range(1, nterms[i] + 1, 2):
        un0 = np.sqrt(np.pi) / (2 * gamma(1 - k / 2) * gamma((3 + k) / 2)) * np.sqrt((2 * k + 1) / 2)
        u[i, :] += un0 * phi[k, :]

plt.plot(x, flist, '-', x, u[0, :], '--', x, u[1, :], ':', x, u[2, :], '-.')
plt.legend(['H(x)', 'n=10', 'n=50', 'n=100'])
plt.show(block=False)

table = np.column_stack( (x, flist, u.T) )
with open("LegendreStep.dat", "w") as f:
    np.savetxt(f, table, fmt='%f', header="table")