Electromagnetic pulse propagation in passive media by path integral methods

TL;DR

The paper develops a time-domain solver for Maxwell's equations in passive media by casting the problem into a Hamiltonian (path-integral) framework and solving with a split-operator, pseudospectral approach. By augmenting electromagnetic fields with matter degrees of freedom, it achieves dispersion-free evolution at the Nyquist limit, exact Gauss-law enforcement, and (un)conditional stability results for multi-resonant Lorentz media, with quadratic convergence in time step and potential higher-order extensions. The method provides explicit operator exponentials for Lorentz models, analyzes energy and norm conservation, and offers robust boundary treatment strategies, making it applicable to photonic media, scattering from dispersive targets, and time-varying media. While FFT-based solvers face aliasing and boundary challenges, the approach promises high accuracy with favorable scaling ($O(N\log N)$ per step) and strong theoretical guarantees on stability and convergence.

Abstract

A novel time domain solver of Maxwell's equations in passive (dispersive and absorbing) media is proposed. The method is based on the path integral formalism of quantum theory and entails the use of ({\it i}) the Hamiltonian formalism and ({\it ii}) pseudospectral methods (the fast Fourier transform, in particular) of solving differential equations. In contrast to finite differencing schemes, the path integral based algorithm has no artificial numerical dispersion (dispersive errors), operates at the Nyquist limit (two grid points per shortest wavelength in the wavepacket) and exhibits an exponential convergence as the grid size increases, which, in turn, should lead to a higher accuracy. The Gauss law holds exactly with no extra computational cost. Each time step requires $O(N\log_2 N)$ elementary operations where $N$ is the grid size. It can also be applied to simulations of electromagnetic waves in passive media whose properties are time dependent when conventional stationary (scattering matrix) methods are inapplicable. The stability and accuracy of the algorithm are investigated in detail.