import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
# Parameter simulasi
c = 1.0 # Kecepatan gelombang
Lx, Ly = 10.0, 10.0 # Dimensi domain
T = 10.0 # Waktu simulasi
dx = dy = 0.1 # Spasi grid spasial
dt = 0.01 # Interval waktu
frame_step = 10 # Ambil 1 frame setiap n langkah waktu
x = np.arange(0, Lx, dx)
y = np.arange(0, Ly, dy)
t = np.arange(0, T, dt)
Nx, Ny = len(x), len(y)
Nt = len(t)
# Cek CFL
CFL = c * dt / dx
print(f"Angka CFL: {CFL:.2f}")
if CFL >= 1:
print("Peringatan: Simulasi mungkin tidak stabil!")
# Kondisi awal (Gaussian pulse)
x0, y0 = Lx / 2, Ly / 2 # Posisi pusat gelombang
sigma = 0.5 # Lebar gelombang
X, Y = np.meshgrid(x, y, indexing='ij')
z0 = np.exp(-((X - x0)**2 + (Y - y0)**2) / (2 * sigma**2))
# Inisialisasi array solusi
z = np.zeros((Nx, Ny, Nt))
z[:, :, 0] = z0
# Mengisi langkah waktu pertama (n=1)
z[1:-1, 1:-1, 1] = z[1:-1, 1:-1, 0] + 0.5 * (CFL**2) * (
z[2:, 1:-1, 0] + z[:-2, 1:-1, 0] + z[1:-1, 2:, 0] + z[1:-1, :-2, 0] - 4 * z[1:-1, 1:-1, 0]
)
# Iterasi waktu dengan broadcasting
for n in range(1, Nt-1):
z[1:-1, 1:-1, n+1] = 2*z[1:-1, 1:-1, n] - z[1:-1, 1:-1, n-1] + (CFL**2) * (
z[2:, 1:-1, n] + z[:-2, 1:-1, n] + z[1:-1, 2:, n] + z[1:-1, :-2, n] - 4*z[1:-1, 1:-1, n]
)
# Boundary condition tetap (Dirichlet)
z[0, :, :] = z[-1, :, :] = z[:, 0, :] = z[:, -1, :] = 0
# Setup plot
fig, ax = plt.subplots(figsize=(6,6))
cmap = ax.imshow(z[:, :, 0], extent=[0, Lx, 0, Ly], origin='lower', cmap='bwr', vmin=-1, vmax=1)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_title('Simulasi Gelombang 2D')
# Fungsi animasi
def animate(frame):
cmap.set_array(z[:, :, frame])
return cmap,
# Buat animasi dengan frame yang direduksi
ani = FuncAnimation(
fig, animate, frames=range(0, Nt, frame_step), interval=20*frame_step, blit=True
)
plt.close()
# Tampilkan animasi
from IPython.display import HTML
HTML(ani.to_jshtml())