Hide code cell source
import matplotlib.pyplot as plt
from abtem.parametrizations import LobatoParametrization
from ase.io import read

import abtem

abtem.config.set({"local_diagnostics.progress_bar": False});

Potentials#

An electron beam interacts with a specimen principally through the Coulomb potential of its electrons and nuclei. Thus the total electrostatic potential of the sample is required for image simulations. Typically, a so-called independent atom model (IAM) is used, which calculates the potential as a superposition of atomic potentials, hence neglecting effects due to valence bonding.

Atomic potential parametrization#

The electron charge distribution of a single atom can be calculated from a first-principles electronic structure calculation, for example using the Hartree-Fock method or density functional theory. Given a charge distribution, the potential can be obtained via Poisson’s equation. Most multislice simulation codes include a parametrization of the atomic potentials, with a table of parameters for each element fitted to the potential calculated elsewhere from first principles.

We show the radial dependence of the electrostatic potential and scattering of five selected elements below. Note that the potentials tend to infinity at \(r=0\) due to the point-like nuclear charge.

symbols = ["C", "N", "Si", "Au", "U"]

parametrization = LobatoParametrization()

potentials = parametrization.line_profiles(symbols, cutoff=2, name="potential")
scattering_factor = parametrization.line_profiles(
    symbols, cutoff=3, name="scattering_factor"
)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(11, 4))
visualization = potentials.show(ax=ax1, legend=False)
visualization.set_ylim([-1e2, 2e3])

scattering_factor.show(legend=True, ax=ax2);
../../_images/76e3ff0455ccb8432a9568a6423997b5882f160c115a94567742e3dd9d47bf02.png

The default parametrization in abTEM, which is shown above, is created by Ivan Lobato[LD14]. We also implement the parametrization by Earl J. Kirkland[Kir10] and the parametrization by Peng[Pen99]. The differences between the parametrizations are generally negligible at low scattering angles, but the parametrization by Lobato is more accurate at higher angles[LD14].

Independent atom model#

The full specimen potential, \(V(r)\), is usually obtained as a linear superposition of atomic potentials

\[ V(r) = \sum_i V_i(r-r_i) \quad, \]

where \(V_i(r)\) is the atomic potential of the \(i\)’th atom. This is known as the independent atom model (IAM), and as a superposition of independent atoms, it neglects bonding effects. While this is in many cases an adequate approximation, there are situations where that is not the case ([MPS21]) which abTEM is particularly designed to address (as discussed later in the walkthrough.

Potential#

Here, we create a Potential object representing an IAM potential of SrTiO3 with the Lobato parametrization. The sampling denotes the spacing of the \(xy\)-samples of the potential, slice_thickness determines the spacing the slices in direction of electron beam and projection determines how those slices are calculated (as described later in this document).

We also repeat the structure to make a larger supercell for visualizations.

srtio3 = read("./data/SrTiO3.cif")
repeated_srtio3 = srtio3 * (2, 2, 6)

potential = abtem.Potential(
    repeated_srtio3,
    sampling=0.05,
    parametrization="lobato",
    slice_thickness=1,
    projection="finite",
)

The potential has 24 slices along the \(z\) propagation direction, as may be determined from getting its length.

len(potential)
24

The projected potential, i.e. the sum of all the slices multiplied by the thickness, may be shown using the show method.

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 4))

abtem.show_atoms(repeated_srtio3, ax=ax1, legend=True)

visualization = (potential.build() * 0.1).show(ax=ax2);
../../_images/ff9d059f57d1ade0139bd5cfa09a25968276866b38d8fab72a0bf7a5b6854884.png

The Potential may be indexed to return a subset of slices. Below we select the first five slices (using the Python syntax for indexing lists, [:5]) and show them all by setting project=False (the default is to show the projected potential) and explode=True (which automatically creates a multi-panel figure; in this instance the latter setting is assumed if the former is set). The titles of the panels are automatically set based on the depth of the slices, note that the slice thickness is not exactly \(1 \ \mathrm{Å}\), it is rounded down to the closest value fitting an integer number of times.

visualization = potential[:5].show(
    project=False,
    explode=True,
    figsize=(14, 5),
    common_color_scale=True,
    cbar=True,
)
../../_images/0a3073b54e520f548782ad5e37f720282b5e16061cb0eedf6c6e3cfd63112999.png

Note

Interactive visualizations You can use our interactive features to scroll through all the slices in your potential.

To enable interactivity, you first have to enable a compatible matplotlib backend. We recommend using ipympl:

%matplotlib ipympl

You can then create an interactive visualization by setting interact=True in the .show method. Make sure that explode is not set to True.

Building and saving the potential#

The Potential does not store the calculated potential slices. Hence, if a simulation, such as STEM requires multiple propagations, each slice have to be calculated multiple times. For this reason, abTEM often precalculates the potential whenever it has to be used more than once.

The potential can be precalculated manually using the build method, but you should typically let abTEM decide whether to precalculate the potential. Note that we also need to compute it, abTEM by default uses so-called lazy arrays (see our description of Dask for detail); a progress bar is displayed by default.

potential_array = potential.build().compute()

This returns an PotentialArray object, which stores each 2D potential slice in a 3D array. The first dimension is the slice index and the last two are the spatial dimensions.

potential_array.shape
(24, 158, 158)

The calculated potential can be stored in the open-source Zarr file format and conveniently read back in.

potential_array.to_zarr("data/srtio3_repeated_potential.zarr", overwrite=True)

abtem.from_zarr("data/srtio3_repeated_potential.zarr");

Slicing the potential#

The multislice algorithm underlying most modern TEM image simulations requires a mathematical discretization of the potential into slices, which is exact in the limit of infinitely thin slices. The slice thickness can be considered a convergence parameter, and since more slices increases the computational cost, an optimum providing sufficient precision at an acceptable cost can be selected.

A reasonable value for slice thickness is generally between \(0.5 \ \mathrm{Å}\) and \(2 \ \mathrm{Å}\); our default value is \(1 \ \mathrm{Å}\). abTEM also provides multiple options for evaluating the integrals required for slicing the potential that make slightly different tradeoffs in terms of precision and performance.

Finite projection integrals#

abTEM implements an accurate finite potential projection method. Numerical integration is used to calculate the integrals of the form

(1)#\[ V_{n}(x, y) = \int_{z_n}^{z_n+\Delta z} V(x,y,z) dz \quad , \]

where \(z_n\) is the \(z\)-position at the entrance of the \(n\)’th slice and \(\Delta z\) is the slice thickness.

We used the abTEM default method of handling finite projection integrals above; additional options and a description of the methods is given in an appendix.

Infinite projection integrals#

Since the finite projection method can be computationally demanding, abTEM also implements potentials using infinite projection integrals. The finite integrals are replaced by infinite integrals, which may be evaluated analytically

\[ \int_{-\infty}^{\infty} V(x,y,z) dz \approx \int_{z_i}^{z_i+\Delta z} V(x,y,z) dz \]

The infinite projection of the atomic potential for each atom is assigned to a single slice. The implementation uses the hybrid real-space/Fourier-space approach by van den Broek [VdBJK15].

Below we create the same SrTiO3 potential as above with infinite projections. The potential looks almost identical to the finite projection, but note that the slice at \(z=2.96 \ \mathrm{Å}\) has zero intensity because there is not atoms in this slice.

potential_infinite = abtem.Potential(
    repeated_srtio3,
    sampling=0.05,
    parametrization="lobato",
    slice_thickness=1,
    projection="infinite",
)

visualization = potential_infinite[:5].show(
    project=False,
    figsize=(14, 5),
    common_color_scale=True,
    cbar=True,
);
../../_images/7ed9b961c4543bc0761dca4d84ad322e967b1d4037aa757ea13fb166cb5eb74f.png

Using infinite projections can be much faster, especially for potentials with a large numbers of atoms. As shown by the time below, the potential takes less than 0.5 s to calculate, which you may compare this to the finite projection calculated above that took 4 s. The error introduced by using infinite integrals is in most, but not all, cases negligible compared to other sources of error.

potential_infinite.build().compute();

Crystal potentials#

Calculating the potential is generally not a significant cost for simulations where the same potential is used in many runs of the multislice algorithm. However, for simulations that require just one or a few wave functions, such as HRTEM and CBED, it can also be an advantage to use the a periodic crystal potential.

The CrystalPotential allows fast calculation for potentials of crystals by tiling a repeating unit of the potential. Below we create a CrystalPotential by tiling the already calculated finite projected potential.

crystal_potential = abtem.CrystalPotential(
    potential_unit=potential, repetitions=(5, 5, 10), seeds=None
)

crystal_potential.build().array
Array Chunk
Bytes 571.38 MiB 571.38 MiB
Shape (240, 790, 790) (240, 790, 790)
Dask graph 1 chunks in 5 graph layers
Data type float32 numpy.ndarray
790 790 240
crystal_potential.build().compute();

We also implement a fast frozen phonon algorithm for crystal potentials by reshuffling precalculated slices [Bar18]; see our walkthrough on frozen phonons.