a-cosmology-group/acg#174

$$ \frac{\partial^2u(x,t)}{\partial t^2} = \frac{\partial^2u(x,t)}{\partial x^2} - \sin(u(x,t)) $$

on a lattice from 0 to 10 with step $dx = 0.005$

and initial conditions:

• Take as initial data $u(x,t_0)$ and $u’(x,t_0)$ for $t_0$ a wave packet traveling at velocity $v$.

• Find $u(x,t)$ at $t = 150$.

• Match the velocity so that the packet travels from one edge of the grating to the other, but does not go outside the grating.

• Draw a picture of the motion of the packet.

i dont understand why my code doesnt work…

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
from scipy.sparse import diags
from matplotlib import animation

L = 10
x0 = 0
sigma = 0.5
p0 = 2
t_final = 150
v = L / t_final


def u0(x):
    return np.exp(-((x - x0) ** 2) / (2 * sigma**2)) * np.cos(p0 * x)


def u0_t(x):
    return -v * np.exp(-((x - x0) ** 2) / (2 * sigma**2)) * np.sin(p0 * x)


def du_dt(t, u):
    d2u_dx2 = D2.dot(u[:N])
    return np.hstack([u[N:], d2u_dx2 - np.sin(u[:N])])


dx = 0.005
x = np.arange(0, L + dx, dx)
N = len(x)

D2 = diags([1, -2, 1], [-1, 0, 1], shape=(N, N)) / dx**2

dt = 0.05

initial_conditions = np.hstack([u0(x), u0_t(x)][:N])

sol = solve_ivp(
    du_dt,
    t_span=[0, t_final],
    y0=initial_conditions,
    t_eval=np.arange(0, t_final, dt),
    method="RK45",
)

fig, ax = plt.subplots(figsize=(8, 6))
ax.set_xlim(0, L)
ax.set_ylim(-2, 2)
ax.set_title("Wave Packet Movement")
ax.set_xlabel("x")
ax.set_ylabel("u(x, t)")

(line,) = ax.plot([], [], lw=2)


def init():
    line.set_data([], [])
    return (line,)


def animate(i):
    line.set_data(x, sol.y[:N, i])
    return (line,)


frame_skip = 5
ani = animation.FuncAnimation(
    fig,
    animate,
    init_func=init,
    frames=range(0, len(sol.t), frame_skip),
    interval=200,
    blit=True,
)
# ani.save("wave_packe1.gif", writer="pillow", fps=60)
plt.show()