"""Test full QFT propagation (Klein‑Gordon and Dirac)."""
import matplotlib.pyplot as plt
import numpy as np
from mercurial.simulation.engine import SimulationEngine
engine = SimulationEngine(dim=10)
# Test 1: Klein‑Gordon with Gaussian initial condition
[docs]
def gaussian_phi(nx, ny, x0, y0, sigma):
x = np.linspace(-5, 5, nx)
y = np.linspace(-5, 5, ny)
X, Y = np.meshgrid(x, y)
return np.exp(-((X - x0) ** 2 + (Y - y0) ** 2) / (2 * sigma**2))
nx, ny = 64, 64
initial_phi = gaussian_phi(nx, ny, 0, 0, 1.0)
res_kg = engine.run_klein_gordon(
(0.0, 5.0), nx=nx, ny=ny, dx=0.2, m=0.5, initial_phi=initial_phi, save_every=20
)
# Plot final field
plt.figure(figsize=(12, 5))
plt.subplot(1, 2, 1)
plt.imshow(initial_phi, cmap="RdBu", origin="lower")
plt.colorbar()
plt.title("Initial φ")
plt.subplot(1, 2, 2)
plt.imshow(res_kg["final_phi"], cmap="RdBu", origin="lower")
plt.colorbar()
plt.title("Final φ after 5s")
plt.tight_layout()
plt.show()
# Test 2: Dirac 1D
res_dirac = engine.run_dirac((0.0, 2.0), nx=256, dx=0.05, m=0.5, save_every=20)
# Plot probability density evolution (space-time)
rho_stack = np.array(res_dirac["rho"])
plt.figure(figsize=(12, 5))
plt.imshow(rho_stack.T, aspect="auto", origin="lower", extent=[0, 2, -6.4, 6.4])
plt.xlabel("Time (s)")
plt.ylabel("Position")
plt.colorbar(label="Probability density")
plt.title("Dirac Field Evolution")
plt.show()
print(f"Klein‑Gordon final energy: {np.sum(res_kg['energy'][-1]):.3f}")