Convergence order in electromagnetic simulation: why mesh refinement compounds

Convergence order is one of the most important properties of a numerical method, but it is often given less attention than mesh size in production simulation. A method has order $p$ if its error satisfies $\|u - u_h\| \le C\,h^p$ for sufficiently small mesh spacing $h$, where $u$ is the exact solution and $u_h$ the numerical approximation. The order $p$ is a fixed property of the discretization and governs how rapidly error decreases with refinement.

For high-frequency electromagnetic scattering the natural mesh-size parameter is not $h$ alone but the dimensionless product $k_0 h$, where $k_0$ is the wavenumber of the incident field. Equivalently, $k_0 h = 2\pi h / \lambda$ — mesh spacing measured in fractions of a wavelength. A solver that satisfies $\|u - u_h\| \le C\,(k_0 h)^p$ holds its error constant across any electrical size as long as the mesh density per wavelength is held constant. Planck's solver achieves this with $p \ge 6$.

The cost of running a simulation is set by the number of degrees of freedom (DOFs) — the count of unknowns the solver must compute (one per basis function or grid point). In three dimensions DOFs scale as $N \sim h^{-3}$, so reducing error by factor $F$ requires growing the DOF count by $F^{3/p}$ — a formula that depends only on $p$ and the dimension. Going from second-order to sixth-order reduces per-decade DOF cost by a factor of ten, and the gap compounds: at four decades of additional accuracy, $p = 2$ requires $10^6$ more DOFs while $p = 6$ requires $10^2$.

One decade of additional accuracy (factor $F = 10$):

Most production electromagnetic simulators are second-order. FDTD (Yee scheme [Yee 1966]) is structurally fixed at $p = 2$. FEM defaults to second-order curvilinear elements; higher-order options exist [Demkowicz 2006] but $p > 3$ is rare in production due to conditioning and basis-complexity costs. The method of moments uses the Rao–Wilton–Glisson basis [Rao 1982] with $p \approx 1.5$ by default; higher-order basis function (HOBF) implementations exist [Graglia et al. 1997; Notaroš 2008] but production runs typically remain at $p \le 3$ because singular-kernel quadrature at high order is delicate, and insufficient quadrature degrades the effective order regardless of basis-function choice. To our knowledge, no commercial solver supports $p \ge 6$ in production on general 3D engineering geometry with material inhomogeneity and sharp interfaces.

Planck Labs is developing a high-order electromagnetic solver with $p \ge 6$ convergence on arbitrary 3D geometry.

References

• Demkowicz, L. (2006). Computing with hp-Adaptive Finite Elements, Volume 1. Chapman & Hall/CRC.
• Graglia, R. D., Wilton, D. R., & Peterson, A. F. (1997). Higher order interpolatory vector bases for computational electromagnetics. IEEE Trans. Antennas Propag. 45, 329–342.
• Notaroš, B. M. (2008). Higher Order Frequency-Domain Computational Electromagnetics. IEEE Trans. Antennas Propag. 56(8), 2251–2276.
• Rao, S. M., Wilton, D. R., & Glisson, A. W. (1982). Electromagnetic scattering by surfaces of arbitrary shape. IEEE Trans. Antennas Propag. 30, 409–418.
• Yee, K. S. (1966). Numerical solution of initial boundary value problems involving Maxwell's equations in isotropic media. IEEE Trans. Antennas Propag. 14, 302–307.