Multislice Propagation & Scattering Model

• Multislice propagation is a computational framework that discretizes a complex medium into thin slices to simulate sequential wave scattering and phase modulation.
• The model integrates inelastic, magnetic, and anisotropic scattering to enable quantitative 3D structure recovery and high-resolution spectroscopy.
• Efficient FFT-based propagation and inverse reconstruction algorithms leverage sparse matrix operations and gradient optimization for precise defect localization.

Multislice propagation and scattering models constitute a unified computational framework for simulating wave interactions with structured media where multiple sequential forward scattering, phase modulation, and diffraction must be rigorously described. Initially established in the context of high-energy electron microscopy, the formalism now encompasses electromagnetic, acoustic, and elastic waves, and underpins both forward modeling and inverse retrieval in 3D structure determination, quantitative spectroscopy, and dynamical diffraction studies.

1. Theoretical Foundations of Multislice Propagation

The core principle of multislice propagation is the decomposition of a complex scattering medium into a stack of thin slabs or "slices" along the main propagation axis (conventionally $z$), each imparting a local phase shift on the wavefield, interleaved with free-space propagation over the inter-slice separation. For electrons, this rests on the paraxial Schrödinger equation: $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$ where $\lambda$ is the wavelength, $\sigma$ the interaction constant, and $V(\mathbf{r}, z)$ the projected potential (Zhang et al., 2022). For electromagnetic waves, the analogous vectorial Helmholtz equation with a full permittivity tensor enables modeling of birefringence and anisotropic media (Mu et al., 2022). The wavefunction or field is evolved slice-by-slice via:

• A local transmission operator encoding the projected potential or permittivity:

$$T_n(\mathbf{r}) = \exp\left[i \varphi_n(\mathbf{r})\right]$$

• A propagation operator, realized in real or reciprocal space, encapsulating the Fresnel or angular spectrum kernel:

$$\mathcal{P}(\Delta z)[f](\mathbf{k}) = \exp\left[-i\pi\lambda\Delta z|\mathbf{k}|^2\right] F[f](\mathbf{k})$$

Multiple scattering is naturally incorporated by recursive application of this sequence, supporting both strong phase and amplitude modulation regimes. The formalism generalizes to tensorial fields and fully vectorized formulations as required by anisotropic or magnetic scenarios (Mu et al., 2022, Lyon et al., 2021).

2. Extensions to Inelastic, Magnetic, and Tensorial Scattering

To address the complete spectrum of matter–wave interactions, multislice propagation extends to:

• Inelastic and incoherent scattering: Incorporation of a complex, absorptive potential allows modeling of vibrational, plasmonic, and core-loss channels, with the "frozen phonon" approximation sampling atomic displacements for TDS (Zeiger et al., 2021, Lei et al., 27 Jun 2025). In frequency-resolved regimes, thermostatted MD snapshots target specific phonon modes, enabling direct simulation of vibrational EELS spectra (Zeiger et al., 2021).
• Magnetic materials and vector potentials: The paraxial equation acquires a vector potential term, leading to a transmission function containing both scalar and vector phase:

$$T(\mathbf{r}_\perp) = \exp\left\{ i\sigma \int V_{\text{elec}} dz + (i e/\hbar)\int A_z dz \right\}$$

Robust parameterization of $A_z$ and $B$ from atomic DFT computations allows affordable large-scale multislice simulations for magnetic contrast and differential-phase-contrast experiments (Lyon et al., 2021).

• Anisotropic and birefringent scattering: In vector optics, the permittivity tensor $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$0 is retained, and each slice's transmission operator acts as a local tensor exponential. Free-space propagation applies per-polarization components via 2D FFTs and multiplication with spatially-varying $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$1 factors (Mu et al., 2022).

3. Algorithmic Structure, Computational Considerations, and Eigenanalysis

Algorithmic workflows are typically organized as follows:

• Slice-wise recursion: For each slice $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$2, the total wave is multiplied by $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$3 and then propagated via $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$4.
• Matrix formalism: For efficient inversion and theoretical analysis, vectorized (flattened) versions of the wavefield are constructed, with operations cast as products of sparse diagonal and unitary (Fourier) matrices. The multislice propagation chain becomes

$$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$5

where $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$6 and $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$7 are diagonal and Fourier-propagator matrices (Bangun et al., 2024, Bangun et al., 2022).

• Eigenstructure equivalence: For crystalline samples, the eigenvalues and eigenvectors of the multislice transmission matrix can be shown, up to trivial phase ambiguities, to match those of the Bloch-wave scattering matrix; the eigenvector sets are related by a 2D Fourier transform (Bangun et al., 2024). The determinant phase of the multislice matrix yields a direct estimate of the mean inner potential.

Typical computational implementations employ FFT-based propagation, per-slice transmission lookup tables (including for magnetic and inelastic terms), and, for inverse imaging, gradient-based optimization back-propagated through the full multislice chain (utilizing frameworks such as PyTorch for automated differentiation) (Myint et al., 2022).

4. Model Variants and Physical Regimes

Multislice propagation models are adapted to a wide range of experimental and physical scenarios, including:

• Channelling/high-angle scatter: For high-angle annular detectors (HAADF), thermal diffuse scattering randomizes the phase between slices, rendering the signal longitudinally incoherent and permitting efficient linear models such as the atomic lensing approximation for ADF-EDX cross-sections (Zhang et al., 2022).
• Reflection (grazing incidence) geometry: In x-ray and optical CSSI, the multislice approach is adapted to simulate multiple beam (dynamical) effects in surface-sensitive geometries by slicing in the propagation axis normal to the substrate, automatically accounting for evanescent field decay and interference (Myint et al., 2022).
• Vortex and OAM channel analysis: Cylindrical-coordinate multislice permits propagation and OAM analysis for electron vortex beams, with explicit expansion into azimuthal harmonics and tracking of angular momentum transfer and mode purity (Löffler et al., 2019).
• Time-domain multislice for layered systems: In acoustic and elastic contexts, a time-domain multislice propagator alternates convolutions with layer-wise Green's functions and algebraic interface scattering matrices (from multiple scattering theory) for exact simulation of wavepacket evolution including backward- and forward-moving components (Los et al., 2019).

5. Inverse Problems, Ptychography, and 3D Structure Recovery

Modern multislice models underpin advanced inversion schemes targeting 3D structural recovery at atomic and nanometer scales:

• Multislice ptychography: In 4D-STEM and ptychography, the multislice formalism generates the forward model relating a stack of 2D projected potentials/slices to measured diffraction patterns. Inverse algorithms employ layer-wise amplitude flow, sparse matrix decomposition, and gradient descent to reconstruct the specimen slices and possibly the probe (Bangun et al., 2022, Bhat et al., 2024, Lei et al., 27 Jun 2025).
• Hollow ptychography with simultaneous EELS: Dual-mode data acquisition partitions electrons between high-angle ptychographic detection and low-angle EELS spectroscopy, with the multislice algorithm ensuring faithful account of multiple elastic scattering and mask-aware inversion for full 3D reconstruction of structure and composition (Lei et al., 27 Jun 2025).
• Coherent surface scattering imaging: In CSSI, the multislice model is integrated with differentiable programming frameworks to reconstruct subsurface 3D morphology from single-shot grazing-incidence coherent diffraction (Myint et al., 2022).
• Defect sensitivity and depth localization: By combining accurate atomic slice potentials, tight control of experimental parameters (defocus, convergence angle, slice thickness), and frozen-phonon averaging, point defect signatures are resolved in depth with sub-nm accuracy by multislice ptychographic inversion (Bhat et al., 2024).

6. Quantitative Benchmarks, Physical Validity, and Model Limitations

• Computational scaling: Library-based atomic lensing and incoherent superposition models accelerate ADF-EDX signal prediction from $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$8 s per configuration (standard multislice) to sub-millisecond with negligible loss in accuracy for HAADF and EDX under column isolation (Zhang et al., 2022).
• Accuracy: Root-mean-square errors of atomic lensing predictions for ADF/EDX cross-sections are consistently within $$\frac{\partial\psi(\mathbf{r}, z)}{\partial z} = \left[\frac{i\lambda}{4\pi} \nabla_{\mathbf{r}}^2 + i\sigma V(\mathbf{r}, z)\right]\psi(\mathbf{r}, z)$$91 atom; mean inner potentials extracted from multislice determinants closely match DFT and holography (Bangun et al., 2024).
• Limitations: Approximations include the paraxial or one-way propagation assumption (neglecting backward or reflected scattering except explicitly included via interface matrices), first-Born or Dyson approximations for intra-slice multiple scattering, neglect of inter-column cross-talk, and single-mode inelastic loss channels. Correction for absorption, gain, or nonlinear index responses is in principle feasible but requires model extension (Mu et al., 2022, Myint et al., 2022).
• Physical regimes: Dynamical (multi-beam) effects in both transmission and reflection geometries are rigorously recovered, provided the slice thickness, numerical sampling, and grid resolution are consistent with physically relevant phase variation and beam divergence (Bangun et al., 2024, Myint et al., 2022).

7. Applications and Future Directions

The versatility of multislice propagation and scattering models is reflected in their use across:

• Quantitative ADF-EDX chemical and thickness mapping in heterogeneous nanostructures and alloys (Zhang et al., 2022).
• Determination of magnetic structure at atomic scale by DPC, EMCD, and magnon scattering (Lyon et al., 2021).
• 3D atomic tomography, correlative STEM-EELS, and localized defect imaging at sub-nanometer spatial and depth resolution (Lei et al., 27 Jun 2025, Bhat et al., 2024).
• Photonic device inverse design, high-precision birefringence mapping, and polarization-resolved optical tomography (Mu et al., 2022).
• Dynamical scattering correction in coherent surface structural imaging, enabling ab initio 3D pattern recovery even in the presence of strong surface and substrate dynamical effects (Myint et al., 2022).
• Time-domain pulse propagation in layered acoustic or electromagnetic stacks, including two-way coupling and causality-respecting boundary reflection via exact interface S-matrices (Los et al., 2019).

Research is ongoing into incorporating higher-order multiple scattering, nonlinearity, full tensorial magneto-optic effects, and hybrid real- and reciprocal-space algorithms for increasing physical fidelity and computational tractability across diverse modalities.