Lab Sheet 11

from header import *

Exercise 1.1

def mandelbrot_test(c, maxiter=100):
    """
    Test if a point c is in the Mandelbrot set.  Return 1 if it is, 0 otherwise.
    The function maps efficiently over numpy arrays.
    """
    z = 0 * c  # array of zeros of the same shape as c
    m = z + 1  # array of ones of the same shape as c
    for k in range(maxiter):
        # Throughout this loop, m is an array of 1s and 0s and z is an array of complex numbers
        # An entry in m is 1 if the iteration starting with the corresponding c has 
        # not yet escaped from the circle of radius 2, and 0 otherwise.
        # The corresponding entry in z is the current value of the iteration if m = 1,
        # and 10 if m = 0.
        z = m * (z**2 + c) + 10 * (1 - m)
        m = (abs(z) <= 2)*1
    return m

n = 1000
xs = np.linspace(-2.0, 0.4, n+1) # real parts
ys = np.linspace(-1.2, 1.2, n+1) # imaginary parts
cs = xs.reshape(1, n+1) + 1j*ys.reshape(n+1, 1) # broadcast sum gives a matrix of complex numbers
ms = mandelbrot_test(cs, 100)

img = np.expand_dims(ms, axis=2) * np.array([-255, -180, -255], dtype=int) + np.array([255, 255, 255], dtype=int)
fig = go.Figure(data=[go.Image(x0=xs.min(), y0=ys.min(), dx=xs[1]-xs[0], dy=ys[1]-ys[0], z=img)])
ts = np.linspace(0, 2*np.pi, 200)
fig.add_trace(go.Scatter(x=np.cos(ts)/2-np.cos(2*ts)/4, y=np.sin(ts)/2-np.sin(2*ts)/4, mode='lines', line=dict(color='red')))
fig.add_trace(go.Scatter(x=np.cos(ts)/4-1, y=np.sin(ts)/4, mode='lines', line=dict(color='red')))
fig.add_trace(go.Scatter(x=-0.1225+0.0945*np.cos(ts), y=0.7449+0.0945*np.sin(ts), mode='lines', line=dict(color='red')))
fig.add_trace(go.Scatter(x=-0.1225+0.0945*np.cos(ts), y=-0.7449+0.0945*np.sin(ts), mode='lines', line=dict(color='red')))

fig.update_layout(width=600, height=600, showlegend=False, xaxis_visible=False, yaxis_visible=False)
fig.show()

Exercise 1.2

def show_curve(n, fig, col, c):
    xs, ys = np.meshgrid(np.linspace(-3,3,200), np.linspace(-3,3,200))
    zs = (np.abs(xs)**n + np.abs(ys)**n)/n**(n/2)
    fig.add_trace(go.Contour(x=xs[0], y=ys[:,0], z=zs, contours=dict(type='constraint', value=1), line_color=c), row=1, col=col)

fig = plotly.subplots.make_subplots(rows=1, cols=3, subplot_titles=('n=6','n=2','n=3'))
show_curve(6, fig, 1, 'red')
show_curve(2, fig, 2, 'red')
show_curve(3, fig, 3, 'red')
fig.update_layout(width=1200, height=400, showlegend=False)
fig.show()

fig = plotly.subplots.make_subplots(rows=1, cols=1)
cols = ['red', 'green', 'blue', 'cyan', 'magenta', 'yellow', 'black', 'orange', 'purple']
for n in range(1,10):
    show_curve(n, fig, 1, cols[n-1])
fig.update_layout(width=400, height=400, showlegend=False)
fig.show()

Exercise 1.3

ns = np.arange(101)
ys = [1]
for n in range(1,101):
    ys.append(ys[-1]*50/n)
zs = ys[50] * np.exp(-(ns-50)**2/100)
fig = go.Figure()
fig.update_layout(width=600, height=400, showlegend=False)
fig.add_trace(go.Scatter(x=ns, y=zs, mode='lines', line=dict(color='blue')))
fig.add_trace(go.Scatter(x=ns, y=ys, mode='markers', marker=dict(color='red')))
fig.show()

Exercise 1.4

zeta = sp.functions.special.zeta_functions.zeta
f = lambda t: sp.re(zeta(0.5+1j*t))
g = lambda t: sp.im(zeta(0.5+1j*t))

ts0 = np.linspace(0, 15, 1000)
zs0 = np.array([zeta(0.5 + 1j*t).evalf() for t in ts0]).astype(complex)
ts1 = np.linspace(0, 50, 1000)
zs1 = np.array([zeta(0.5 + 1j*t).evalf() for t in ts1]).astype(complex)
fig = plotly.subplots.make_subplots(rows=1, cols=2)
fig.update_layout(width=1000, height=400, showlegend=False)
fig.add_trace(go.Scatter(x=ts0, y=np.real(zs0), mode='lines', line=dict(color='red')), row=1, col=1)
fig.add_trace(go.Scatter(x=ts0, y=np.imag(zs0), mode='lines', line=dict(color='blue')), row=1, col=1)
fig.add_trace(go.Scatter(x=np.real(zs1), y=np.imag(zs1), mode='lines', line=dict(color='green')), row=1, col=2)
fig.show()

ts = np.linspace(14.1, 14.15, 1000)
zs = np.array([zeta(0.5 + 1j*t).evalf() for t in ts]).astype(complex)
fig = go.Figure()
fig.update_layout(width=600, height=400, showlegend=False)
fig.add_trace(go.Scatter(x=ts, y=np.real(zs), mode='lines', line=dict(color='red')))
fig.add_trace(go.Scatter(x=ts, y=np.imag(zs), mode='lines', line=dict(color='blue')))
fig.show()

Exercise 1.4

x = sp.symbols('x')
y = sp.exp(-1/x)

display(sp.simplify(sp.diff(y,x)))
display(sp.simplify(sp.diff(y,x,2)))
display(sp.simplify(sp.diff(y,x,3)))
display(sp.simplify(sp.diff(y,x,4)))
display(sp.simplify(sp.diff(y,x,5)))

$\displaystyle \frac{e^{- \frac{1}{x}}}{x^{2}}$

$\displaystyle \frac{\left(1 - 2 x\right) e^{- \frac{1}{x}}}{x^{4}}$

$\displaystyle \frac{\left(6 x^{2} - 6 x + 1\right) e^{- \frac{1}{x}}}{x^{6}}$

$\displaystyle \frac{\left(- 24 x^{3} + 36 x^{2} - 12 x + 1\right) e^{- \frac{1}{x}}}{x^{8}}$

$\displaystyle \frac{\left(120 x^{4} - 240 x^{3} + 120 x^{2} - 20 x + 1\right) e^{- \frac{1}{x}}}{x^{10}}$

p = [sp.simplify(sp.diff(y,x,k) * x**(2*k)/y) for k in range(10)]
for k in range(10):
    display(Latex("$p_{" + str(k) + "}(x) = " + sp.latex(p[k]) + "$"))

$p_{0}(x) = 1$

$p_{1}(x) = 1$

$p_{2}(x) = 1 - 2 x$

$p_{3}(x) = 6 x^{2} - 6 x + 1$

$p_{4}(x) = - 24 x^{3} + 36 x^{2} - 12 x + 1$

$p_{5}(x) = 120 x^{4} - 240 x^{3} + 120 x^{2} - 20 x + 1$

$p_{6}(x) = - 720 x^{5} + 1800 x^{4} - 1200 x^{3} + 300 x^{2} - 30 x + 1$

$p_{7}(x) = 5040 x^{6} - 15120 x^{5} + 12600 x^{4} - 4200 x^{3} + 630 x^{2} - 42 x + 1$

$p_{8}(x) = - 40320 x^{7} + 141120 x^{6} - 141120 x^{5} + 58800 x^{4} - 11760 x^{3} + 1176 x^{2} - 56 x + 1$

$p_{9}(x) = 362880 x^{8} - 1451520 x^{7} + 1693440 x^{6} - 846720 x^{5} + 211680 x^{4} - 28224 x^{3} + 2016 x^{2} - 72 x + 1$

xs = np.linspace(1e-6, 2, 1000)
ys = sp.lambdify(x, y, 'numpy')(xs)
fig = go.Figure()
fig.update_layout(width=1200, height=400, showlegend=False)
fig.add_trace(go.Scatter(x=xs, y=ys, mode='lines', line=dict(color='blue')))
for k in range(1,5):
    yk = sp.diff(y,x,k)
    yk_fun = sp.lambdify(x, yk, 'numpy')
    ys = yk_fun(xs)
    ys = ys / np.max(ys)
    fig.add_trace(go.Scatter(x=xs, y=ys, mode='lines'))
fig.show()

fig = go.Figure()
fig.update_layout(width=1200, height=400, showlegend=False, yaxis_range=[0,0.002])
xs = np.linspace(1e-6, 0.1, 1000)
for k in range(5):
    yk = sp.diff(y,x,k)
    yk_fun = sp.lambdify(x, yk, 'numpy')
    ys = yk_fun(xs)
    fig.add_trace(go.Scatter(x=xs, y=ys, mode='lines'))
fig.show()

del x, y, p

Exercise 1.5

COSH = lambda x: (sp.exp(x) + sp.exp(-x))/2

x, t = sp.symbols('x t')
q = sp.sqrt(2)
p = sp.log(3 + 2*q)
r = x - 4*t
s = q * (x - 8*t)
T = 32 * COSH(2*r-p) + 16 * COSH(2*s-p) + 16
B = 4*(1+q) * COSH(r)*COSH(s) + (4*q-8)*sp.exp(r+s)
phi = [
 2*COSH(r)**(-2),
 4*COSH(s)**(-2),
 2*COSH(r-p)**(-2),
 4*COSH(s-p)**(-2),
 T/B**2
]
phi_fun = sp.lambdify((x,t),phi,'numpy')

K = lambda u: sp.diff(u,t) + sp.diff(u,x,3) + 6 * u * sp.diff(u,x)

[sp.simplify(K(phi[k])) for k in range(5)]

[0, 0, 0, 0, 0]

phi[4]

$\displaystyle \frac{16 \cosh{\left(2 \sqrt{2} \left(- 8 t + x\right) - \log{\left(2 \sqrt{2} + 3 \right)} \right)} + 32 \cosh{\left(8 t - 2 x + \log{\left(2 \sqrt{2} + 3 \right)} \right)} + 16}{\left(\left(-8 + 4 \sqrt{2}\right) e^{- 4 t + x + \sqrt{2} \left(- 8 t + x\right)} + \left(4 + 4 \sqrt{2}\right) \cosh{\left(\sqrt{2} \left(- 8 t + x\right) \right)} \cosh{\left(4 t - x \right)}\right)^{2}}$