Computational_Programs/node2contour.py at master · erikwendtmath/Computational_Programs · GitHub

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
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
import numpy as np
from numpy.polynomial import Polynomial
import matplotlib.pyplot as plt
from scipy.interpolate import BarycentricInterpolator
import make1DNodes

# a contour plot for the complex values of a lagrange polynomial with roots in [-1,1]
def poly_contour(ax,roots,n,nlevels, title=None):
    p = Polynomial.fromroots(roots)
    init_xvec = np.linspace(-1.5,1.5,num=n)
    init_yvec = np.linspace(-1.5,1.5,num=n)
    x, y = np.meshgrid(init_xvec,init_yvec)
    z = np.abs(p(x +1j*y))
    log_abs = np.log(np.maximum(np.abs(z), 1e-14))

    vmin, vmax = log_abs.min(), log_abs.max()
    levels = np.linspace(vmin, vmax, nlevels)

    cs = ax.contour(x, y, log_abs, levels=levels, cmap='viridis')
    ax.clabel(cs, inline=True, fontsize=8, fmt="%.2f")

    ax.set_aspect('equal')

    if title is not None:
        ax.set_title(title,fontsize=9)

    return cs

# create the Lebesgue function for a given set of roots
def lebesgueFnc(roots,x):
    roots = np.asarray(roots)
    y = np.zeros(len(x),dtype=float)

    # run through each node point
    for ind, root in enumerate(roots):
        root_set = roots[:ind]+roots[ind+1:]

    px = p(x)
    Lambda = np.zeros_like(x, dtype=float)

    for j, xj in enumerate(roots):
        denom = (x - xj) * dp(xj)

        # avoid division by zero at nodes
        lj = np.zeros_like(x, dtype=float)
        mask = np.abs(x - xj) > 1e-14
        lj[mask] = px[mask] / denom[mask]

        # enforce exact value 1 at node
        lj[~mask] = 1.0

        Lambda += np.abs(lj)

    return Lambda

# create the root sets and their titles
root_sets = {
    "Uniform": np.linspace(-1,1,64),
    "Chebyshev 1st Kind": make1DNodes.cheb1(64),
    "Chebyshev 2nd Kind": make1DNodes.cheb2(63),
    "Cantor Set": make1DNodes.cantorpts(level=6),
    "1/2.5 and 1/3.0 Weighter Cantor": make1DNodes.cantorpts(level=6,step=1/3.0,left_step = 1/2.5, right_step= 1/3.0),
    "1/4.0 Cantor": make1DNodes.cantorpts(level=6,step=1/4.0),
    "Continued Fractions Symmetric:": make1DNodes.cfpts(iter=3,maps=[2.01,-2.01,3,-3]),
    "Continued Fractions Odds:": make1DNodes.cfpts(iter=2,maps=[2,3,4,5,6,7,8,9])
}

fig , axes = plt.subplots(2,4, figsize = (8,8))

for ax, (name, roots) in zip(axes.flat, root_sets.items()):
    poly_contour(ax,roots,n= 500, nlevels = 6, title=name)

plt.tight_layout()

# roots = np.linspace(-1,1,4)
# n=4
# x = np.linspace(-1, 1, n)
# Lambda = lebesgueFnc(roots, x)
# if 1==1:
#     fig, ax = plt.subplots()

#     ax.plot(x, Lambda)
#     ax.set_xlabel("x")
#     ax.set_ylabel(r"$\Lambda(x)$")
#     ax.grid()

plt.show()