Cowley–Moodie Multislice Solution

• Cowley–Moodie multislice solution is a computational framework that models dynamical electron scattering in crystalline materials using a thin-slice decomposition approach.
• It employs alternating transmission and Fresnel propagation operators to simulate phase modulation and diffraction effects in TEM.
• The transmission-matrix formalism facilitates rigorous eigenstructure analysis and precise extraction of physical parameters like the mean inner potential.

The Cowley–Moodie multislice solution is a foundational theoretical and computational framework for modeling dynamical electron scattering in crystalline materials, particularly within transmission electron microscopy (TEM). It systematically decomposes the three-dimensional crystal potential into consecutive thin slices, capturing both phase modulation by the local projected potential and Fresnel propagation between slices. Recent reformulations cast the multislice algorithm into a “transmission-matrix” formalism, enabling rigorous eigenstructure analysis and direct comparison with the Bloch-wave scattering matrix. This establishes the conditions under which both formulations are formally equivalent (modulo $2\pi$ phase ambiguities), and supports physically meaningful parameter extraction, such as the mean inner potential, directly from transmission matrix determinants (Bangun et al., 2024).

1. Multislice Decomposition and Transmission Operators

The Cowley–Moodie approach partitions the specimen’s three-dimensional electrostatic potential $V(x,y,z)$ into $M$ slices of uniform thickness $\Delta z$. The projected potential for slice $n$ is

$V_n(x,y) = \int_{z_n}^{z_n+\Delta z} V(x,y,z)\,dz.$

Each slice imparts a phase shift to the traversing electron wave, described by the transmission operator

$$T_n(x,y) = \exp[i\,\sigma\,V_n(x,y)], \qquad n=1,\ldots,M,$$

where the interaction constant is

$\sigma = \frac{2\pi m e}{h^2 k_0},$

with $m$ as the relativistic electron mass, $e$ the elementary charge, $h$ Planck’s constant, and $k_0=2\pi/\lambda$ the relativistically-corrected electron wavenumber.

Upon traversing slice $n$, the wave immediately after is

$\Psi_n^+(x,y) = T_n(x,y)\,\Psi_{n-1}(x,y).$

2. Propagation Between Slices: Fresnel Operator

Electron propagation between adjacent slices over distance $\Delta z$ is described as Fresnel (paraxial) propagation. In Fourier space, this is a multiplication: $$\Psi^-(q_x,q_y) = \exp\Big[-i\pi\lambda\Delta z(q_x^2 + q_y^2)\Big]\,\Psi^+(q_x,q_y).$$ The equivalent real-space operation is implemented efficiently via fast Fourier transforms (FFT) as

$$\Psi^-(x,y) = \mathcal{F}^{-1}_{2D} \left[ P(q_x,q_y)\,\mathcal{F}_{2D}[\Psi^+(x,y)] \right],$$

with $P(q_x,q_y)$ as the Fresnel propagator.

3. Multislice Transmission Matrix Construction

The complete multislice algorithm applies alternating transmission and propagation operators for $M$ slices, such that the exit wave is

$$\Psi_M^-(x,y) = \Big[P\,T_M\,P\,T_{M-1}\cdots P\,T_1\Big]\Psi_{\text{in}}(x,y).$$

Upon discretizing the lateral plane to an $N \times N$ grid, each transmission operator $T_n$ becomes a diagonal $N^2 \times N^2$ matrix $O_n$, while propagation between slices is mediated by the matrix $G = F_{2D}^{-1} D F_{2D}$, where $D$ is a diagonal matrix of Fresnel phase factors.

The full real-space multislice operator is

$$A = \prod_{n=1}^{M} G\,O_n \in \mathbb{C}^{N^2 \times N^2}.$$

For reciprocal-space (diffraction) analysis, the transmission matrix is conjugated by the $2$D Fourier transform: $\widehat{S} = F_{2D} A F_{2D}^{-1}.$

4. Bloch-Wave Scattering Matrix and Eigenstructure

The Bloch-wave formulation constructs the scattering matrix by exponentiating the crystal “structure matrix” $B$, assembled from Fourier components of $V$: $$S = \exp\left(i2\pi T B\right) = C\,\Lambda\,C^{-1},$$ with $T = M\Delta z$ the total thickness, $\{\gamma_j\}$ the eigenvalues of $B$, and $C$ the eigenvector matrix of Bloch-wave Fourier coefficients. The diagonal matrix $\Lambda$ contains eigenphases $\exp(2\pi i \gamma_j T)$. Both $S$ and $\widehat{S}$ are unitary.

Diagonalization of the multislice transmission matrix yields

$$\widehat{S} = W V W^{-1}, \quad V_{jj} = e^{i2\pi\theta_j},$$

with $|V_{jj}| = 1$.

5. Formal Equivalence: Spectral and Structural Analysis

Quantitative comparison of $S$ and $\widehat{S}$ uses the Frobenius norm,

$$\|S-\widehat S\|_F^2 = 2N^2 - 2\sum_{j=1}^{N^2} \cos\left[2\pi(\gamma_j T - \theta_j)\right].$$

The matrices are equivalent if the following hold:

• The eigenvector bases coincide up to a $2$D Fourier transform: $C = F_{2D} W$.
• The eigenphase differences are integer multiples of $2\pi$: $\gamma_j T - \theta_j = 2\pi n_j$, $n_j\in\mathbb{Z}$.

Thus, the spectral representations of multislice and Bloch-wave formulations are congruent modulo $2\pi$ phase factors under these conditions (Bangun et al., 2024).

6. Physical Parameter Extraction: Mean Inner Potential

The determinant of the transmission matrix encodes the total projected potential: $$\det(\widehat{S}) = \prod_{n=1}^{M} \det(G\,O_n) = \exp\Big[-i\big(\pi\lambda Q T - \sigma V_{\text{total}}\big)\Big],$$ where $Q = \sum_{q_x,q_y} (q_x^2 + q_y^2)$ and $V_{\text{total}} = \sum_n \sum_{x,y} V_n(x,y)$. The mean inner potential (MIP), a quantity of central importance in electron scattering, is then extracted as

$$V_{\text{total}} = \frac{\pi\lambda Q T - \arg\det(\widehat S)}{\sigma} + 2\pi k, \qquad k \in \mathbb{Z},$$

and the normalized MIP for the unit-cell volume $\Omega$ as $V_0 = V_{\text{total}}/\Omega$.

7. Computational Considerations and Underlying Assumptions

The Cowley–Moodie multislice approach relies on several critical approximations:

• The thin-slice (phase-grating) approximation: $T_n\approx\exp[i\sigma V_n]$ assumes linearizable intra-slice multiple scattering.
• Paraxial propagation validity: $\lambda\Delta z(q_x^2 + q_y^2) \ll 1$ and small electron tilt angles are required.
• Efficient implementation: Each slice update via $G\,O_n$ exploits FFTs for computational scaling of $O(N^2\log N)$, which is tractable for large grids. By contrast, direct Bloch-wave eigendecomposition is $O(N^6)$, making the matrix multislice formulation preferable for practical simulation at scale.

The formalism thus not only establishes a rigorous equivalence between widely used dynamical electron scattering models but also provides an efficient, physically insightful, and numerically robust method for TEM simulation and analysis (Bangun et al., 2024).