# Concurrent multi-peak Bragg coherent x-ray diffraction imaging of 3D nanocrystal lattice displacement via global optimization

Siddharth Maddali, Travis Frazer, Nazar Delegan, Katherine Harmon, Sean Sullivan, Marc Allain, Wonsuk Cha, Alan Dibos, Ishwor Poudyal, Saugat Kandel, et al.

► To cite this version:

Siddharth Maddali, Travis Frazer, Nazar Delegan, Katherine Harmon, Sean Sullivan, et al. Concurrent multi-peak Bragg coherent x-ray diffraction imaging of 3D nanocrystal lattice displacement via global optimization. npj Computational Materials, 2023, 9 (1), pp.77. 10.1038/s41524-023-01022-7. hal-04412853

HAL Id: hal-04412853

https://hal.science/hal-04412853v1

Submitted on 19 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
ARTICLE

OPEN

<figure><img src="image_2.png" /></figure>

Check for updates

# Concurrent multi-peak Bragg coherent x-ray diffraction imaging of 3D nanocrystal lattice displacement via global optimization

Siddharth Maddali<fnref n="1" />, Travis D. Frazer<fnref n="1" />, Nazar Delegan<fnref n="1" />,², Katherine J. Harmon<fnref n="1" />, Sean E. Sullivan<fnref n="1" />,², Marc Allain³, Wonsuk Cha<fnref n="1" />,⁴, Alan Dibos⁵, Ishwor Poudyal<fnref n="1" />, Saugat Kandel<fnref n="1" />,⁴, Youssef S. G. Nashed⁶,⁸, F. Joseph Heremans<fnref n="1" />,²,⁷, Hoydoo You<fnref n="1" />, Yue Cao<fnref n="1" /> and Stephan O. Hruszkewycz<fnref n="1" />✉

In this paper we demonstrated a method to reconstruct vector-valued lattice distortion fields within nanoscale crystals by optimization of a forward model of multi-reflection Bragg coherent diffraction imaging (MR-BCDI) data. The method flexibly accounts for geometric factors that arise when making BCDI measurements, is amenable to efficient inversion with modern optimization toolkits, and allows for globally constraining a single image reconstruction to multiple Bragg peak measurements. This is enabled by a forward model that emulates the multiple Bragg peaks of a MR-BCDI experiment from a single estimate of the 3D crystal sample. We present this forward model, we implement it within the stochastic gradient descent optimization framework, and we demonstrate it with simulated and experimental data of nanocrystals with inhomogeneous internal lattice displacement. We find that utilizing a global optimization approach to MR-BCDI affords a reliable path to convergence of data which is otherwise challenging to reconstruct.

npj Computational Materials (2023)9:77; https://doi.org/10.1038/s41524-023-01022-7

## INTRODUCTION

Bragg coherent diffraction imaging (BCDI) is a synchrotron-based characterization technique that utilizes coherent x-ray illumination and three-dimensional (3D) phase retrieval algorithms to inter-rogate the internal structure of sub-micron-sized crystalline materials at spatial scales of tens of nanometers<fnref n="1" />,². This is possible because Bragg coherent diffraction is sensitive to the shape of the crystal and its internal atomic lattice displacement field projected along the direction of the scattering vector. This information is encoded in the amplitude and phase of the complex-valued wave that creates a fringe pattern about the diffracted Bragg peak. However, as with other x-ray scattering techniques which face the crystallographic "phase problem," BCDI measurements record only the intensity of the diffracted wave, providing no direct information as to the phase. Thus phase retrieval approaches are needed to gain access to the reconstructed 3D image of the crystal and its internal lattice distortion field. The typical approach to solving this phase problem in BCDI is via iterative methods³⁻⁵ that update an estimate of the 3D image based on the error between its diffraction (emulated via a 'forward model') and the experimentally observed diffraction of a Bragg peak from the sample. BCDI of this sort uses a family of relatively straightforward iterative phase retrieval algorithms to reconstruct the spatial distribution of a single scalar component of the lattice distortion field of a nanocrystal, with many examples in the literature from diverse materials systems⁶⁻⁸. BCDI methods have also been developed to directly retrieve the full 3D lattice distortion field in a nanocrystal by solving the phase retrieval problems for

multiple BCDI datasets from different independent Bragg reflections of the same crystal.

The multi-reflection approach (MR-BCDI) involves a more complex measurement as compared to standard BCDI, requiring multiple Bragg peaks from a single nanoparticle to be measured⁹. This further complicates the image reconstruction process by imposing additional constraints, which then become difficult to reconcile simultaneously with the target image of interest. Various strategies of MR-BCDI image reconstruction have been developed that either utilize phase retrieval methods adapted from standard single-Bragg-peak approaches to reconstruct experimental data sets¹⁰⁻¹², or that explore more sophisticated multi-Bragg peak reconstruction concepts within the realm of simplified simulated data¹³,¹⁴ that are more difficult to translate to as-measured MR-BCDI data sets.

In this work, we introduce an MR-BCDI reconstruction approach that can reconstruct experimental data directly and that is compatible with modern optimization methods. As compared to other MR-BCDI methods, our approach enables global fitting to the full set of intensity patterns from multiple Bragg peaks in order to determine the scattering volume of a crystal and its internal lattice distortion field, affording a much more reliable path to convergence. In this paper, we present this forward model, we implement it within the stochastic gradient descent optimization framework common in machine learning model training, and we demonstrate the method with simulated and experimental data. We focus in particular on enabling accurate reconstructions of crystals that contain discontinuities in their lattice distortion fields.

<footnote n="1">¹Materials Science Division, Argonne National Laboratory, Lemont, IL 60439, USA. ²Center for Molecular Engineering, Argonne National Laboratory, Lemont, IL 60439, USA. ³Aix-Marseille Univ, CNRS, Centrale Marseille, Institut Fresnel, Marseille, France. ⁴X-ray Science Division, Argonne National Laboratory, Lemont, IL 60439, USA. ⁵Center for Nanoscale Materials, Argonne National Laboratory, Lemont, IL 60439, USA. ⁶Mathematics & Computer Science Division, Argonne National Laboratory, Lemont, IL 60439, USA. ⁷Pritzker School of Molecular Engineering, University of Chicago, 5640 South Ellis Avenue, Chicago, IL 60637, USA. ⁸Present address: SLAC National Accelerator Laboratory, Menlo Park, CA, USA. ✉email: shrus@anl.gov</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 1 BCDI measurement schematic.** BCDI setup at Beamline 34-ID-C of the Advanced Photon Source showing the laboratory and detector frames $[\hat{s}_1 \hat{s}_2 \hat{s}_3]$ and $[\hat{k}_1 \hat{k}_2 \hat{k}_3]$ for a single Bragg condition. The diffractometer rotations in this figure can be redefined as needed to generate forward models of MR-BCDI measurements at other synchrotron beam lines with differing geometries. $k_i$ and $k_f$ are drawn to satisfy the current diffraction condition with the area detector shown.</figcaption></figure>

# RESULTS

## The multi-reflection BCDI forward model

A comprehensive forward model for MR-BCDI must tackle several prominent challenges inherent to the geometry of the measurement. First, each Bragg peak measurement is sensitive to a specific projection of the displacement field **u** within the crystal sample along the reciprocal lattice vector **G** of the measured reflection. Second, the crystal and detector must be physically rotated to different orientations to satisfy each different Bragg reflection, resulting in different orientations (or views) of the sample with respect to the detector. Third, the reciprocal space sampling basis of each Bragg peak is generally different and non-orthogonal. This comes about because parallel 2D slices of the 3D diffraction intensity patterns about each Bragg peak are obtained via a fine angular scan (rocking curve) of the crystal by one of the sample rotation stages of the diffractometer after the sample and detector are oriented for diffraction. Thus, flexible multi-axis rotation and adaptable and accurate re-sampling of the sample object space to the different bases dictated by the measured data are needed.

A final requirement, as discussed in the next section, is that in order for this forward model to be suitable for inversion using modern numerical optimization toolkits available as part of machine learning packages, it should be composable of differentiable elementary functions. This allows global optimization of a loss function based on rapid evaluation of gradients through the numerical process of automatic differentiation, which has recently gained traction in x-ray coherent diffraction imaging¹⁵. As we describe here, our MR-BCDI forward model is designed to accommodate these requirements.

The forward model for MR-BCDI is built upon a reference frame in which relevant conventions and operations can be defined. To describe the nanocrystal containing a distortion field we wish to reconstruct, we adopt a right-handed ortho-normal laboratory frame $[\hat{s}_1 \hat{s}_2 \hat{s}_3]$ (Fig. 1). This frame is defined such that the incident x-ray beam path $\hat{k}_i$ is coincident with $\hat{s}_3$, and $\hat{s}_2$ is vertically upward. In this frame, we define the following quantities needed for the forward model:

* The three-dimensional scalar field $\mathcal{A}(\mathbf{x})$ and the 3D vector field $\mathbf{u}(\mathbf{x})$, respectively denoting the spatial distribution of scattering amplitude and the spatial distribution of the relative lattice displacement vector within the crystal at locations $\mathbf{x} \in \mathbb{R}^3$. We assume $\mathcal{A}(\mathbf{x}) \in [0, 1]$ and define $\mathcal{A}: \mathbb{R}^3 \to \mathbb{R}$ and $\mathbf{u}: \mathbb{R}^3 \to \mathbb{R}^3$.

* The bounding box $\mathcal{V}$ that should be larger than the expected size of the crystal. Outside of this volume, the condition $\mathcal{A}(\mathbf{x}) = 0$ and $\mathbf{u}(\mathbf{x}) = 0$ is enforced. $\mathcal{A}$ and $\mathbf{u}$ are discretized on a given grid within $\mathcal{V}$ with a voxel size $s_0$ along each of the Cartesian axes $[\hat{s}_1 \hat{s}_2 \hat{s}_3]$. Typical voxel sizes $s_0$ for BCDI experiments are ~10-nm due to the resolution limits of the method, a length scale that spans many lattice unit cells in the physical crystal.

* The reciprocal lattice vectors corresponding to the multiple Bragg reflections ($M$ in number) measured within a MR-BCDI data set are denoted by the set $\{\mathbf{G}_i\}_{i=1}^M$. The physical units of the $\mathbf{G}_i$ are chosen to be the inverse of those chosen for $\mathbf{u}(\mathbf{x})$, i.e., $\|\mathbf{G}_i\| = 1/d_i$, where $d_i$ is the corresponding spacing of atomic planes diffracting to a given peak.

* The $3 \times 3$ rotation matrices $\{\mathcal{R}_i\}_{i=1}^M$ that act upon the nanocrystal in the lab frame to bring its crystal lattice into each of the $\mathbf{G}_i$ Bragg conditions. This matrix represents a composite rotation of the object to orient the set of sample rotation stages and the detector positioning stages used in the MR-BCDI measurement, providing the orientation of the sample as viewed from the detector. We note that $\mathcal{R}_i$ represents an active rotation of the object that can be implemented by applying a complementary passive rotation $\mathcal{R}_i^{-1}$ to the coordinates of that object, as is done in this work.

* The $3 \times 3$ matrices $\{\mathbf{B}_{\text{real}}^{(i)}\}_{i=1}^M$, whose columns form a basis of 3D sampling vectors in real space in the frame $\hat{k}_1, \hat{k}_2, \hat{k}_3$ of the detector (see Fig. 1). The elements of $\{\mathbf{B}_{\text{real}}^{(i)}\}_{i=1}^M$ are determined by the corresponding reciprocal basis matrix $\{\mathbf{B}_{\text{recip}}^{(i)}\}_{i=1}^M$ which depends on the detector pixel size, detector orientation, the angular increments of the rocking curve, sample rotation axis, and the signal-space array size chosen to encompass the Bragg peak. It is important to note that $\{\mathbf{B}_{\text{real}}^{(i)}\}_{i=1}^M$ for every Bragg peak scan is unique.

* Dataset-dependent scaling factors $\{\chi_i | \chi_i \in \mathbb{R}\}_{i=1}^M$ for the scattering amplitude $\mathcal{A}$ that rescales the sample exit wave amplitude of the model to correspond to the relative intensities of the set of 3D Bragg peak measurements.

These elements can be used to make a MR-BCDI forward model that addresses the above-mentioned challenges and makes the inversion of the forward model tractable. The projection of $\mathbf{u}(\mathbf{x})$ along $\mathbf{G}_i$ allows the 3D complex-valued amplitude and phase of the crystal exit wave field for a given Bragg reflection to be determined as has been established in the BCDI literature:
$\psi^{(i)}(\mathbf{x}) = A(\mathbf{x}) \exp[i2\pi\mathbf{G}^T\mathbf{u}(\mathbf{x})]$. The complex-valued object consistent with the i'th sample orientation in the reference frame of the detector is given by $\psi^{(i)}(\mathbf{R}_i^{-1}\mathbf{x}) = A(\mathbf{R}_i^{-1}\mathbf{x}) \exp[i2\pi\mathbf{G}^T\mathbf{u}(\mathbf{R}_i^{-1}\mathbf{x})]$ and represents the 3D object exit wave to be propagated to the far field for a particular Bragg reflection.

With the object in the detector coordinate frame, the issue of non-orthogonal reciprocal space sampling in BCDI and how to account for its impact in the reconstruction image space also needs to be accounted for. Here, we use the construction of $\mathbf{B}_{\text{real}}$ as derived in Reference$^{16,17}$ as a means to account for this effect. $\mathbf{B}_{\text{real}}^{(i)}$ is defined with respect to the i'th detector frame and originates from the reciprocal space sampling basis $\mathbf{B}_{\text{recip}}^{(i)}$. Using the real space basis $\mathbf{B}_{\text{real}}^{(i)}$, the i'th object exit wave can be appropriately discretized by applying it to the 3D pixel array index space of the Bragg peak measurement. This can be implemented by defining $\mathbf{m}$, a $3 \times N_{\text{vox}}$ matrix of integer array coordinates with $N_{\text{vox}}$ being the total number of voxels in the data of the BCDI scan under consideration. Thus, we can substitute $\mathbf{x}$ with the discretized $\mathbf{B}_{\text{real}}^{(i)}\mathbf{m}$:

$$ \psi_{\mathbf{m}}^{(i)} = \psi^{(i)}(\mathbf{R}_i^{-1}\mathbf{B}_{\text{real}}^{(i)}\mathbf{m}) = A(\mathbf{R}_i^{-1}\mathbf{B}_{\text{real}}^{(i)}\mathbf{m}) \exp[i2\pi\mathbf{G}^T\mathbf{u}(\mathbf{R}_i^{-1}\mathbf{B}_{\text{real}}^{(i)}\mathbf{m})] \quad (1) $$

The quantity $\mathbf{B}_{\text{real}}\mathbf{m}$ produces a grid of sampling points in real space along a non-orthogonal basis. The discrete Fourier transform (DFT) of $\psi_{\mathbf{m}}^{(i)}$ from Equation (1) yields the far field diffraction pattern $\psi^{(i)}$ sampled on the conjugate grid $\mathbf{B}_{\text{recip}}^{(i)}\mathbf{n}$ of the measurement data space, where the integer vector set $\mathbf{n}$ is equivalent to $\mathbf{m}$. Both $\mathbf{m}$ and $\mathbf{n}$ define the equivalent 3D array index space of the diffraction pattern ($\mathbf{n}$) and of the sample ($\mathbf{m}$). We distinguish these indices in a manner consistent with the convention of DFT notation. This construction links the discrete sampling of both sample space and signal space via the DFT and affords significant flexibility to apply different rotation and coordinate transformations of the sample corresponding to different Bragg peak measurements:

$$ \psi_{\mathbf{m}}^{(i)} = \psi^{(i)}(\mathbf{R}_i^{-1}\mathbf{B}_{\text{real}}^{(i)}\mathbf{m}) \stackrel{\mathcal{F}}{\rightleftharpoons} \psi^{(i)}(\mathbf{R}_i^{-1}\mathbf{B}_{\text{recip}}^{(i)}\mathbf{n}) = \psi_{\mathbf{n}}^{(i)} \quad (2) $$

A necessary condition for implementing our forward model is that, prior to rotations into the detector frame, $A(\mathbf{x})$ and $\mathbf{u}(\mathbf{x})$ need to be discretized in the lab frame along an orthonormal grid with spacing $s_0$ that is similar to the magnitude of the values of $\mathbf{B}_{\text{real}}$ for the series of Bragg reflection data sets. As a consequence, the elements of all $\mathbf{B}_{\text{real}}^{(i)}$ should be narrowly distributed about a mean value, which can be used to designate $s_0$. Additionally, the degree of orthogonality of each $\mathbf{B}_{\text{real}}^{(i)}$ basis should be as high as possible. Both of these factors inform the design of angular rocking curve scans with suitable angular motors and angular increments at each Bragg peak. (More detail on this topic is found in the Supplementary Information.)

Another important aspect is the size of the object bounding box $V$ relative to the total size of the buffered array in which it is nested. For each combined re-sampling and rotation transformation $T_i \equiv \mathbf{R}_i^{-1}\mathbf{B}_{\text{real}}^{(i)}$, the object will be sheared, rotated, and rescaled, as will the bounding box. Since our implementation of the MR-BCDI forward problem uses the fast Fourier transform, it must be ensured that the bounding box $V$, after transformation into each of the $M$ detector frames, spans less than half of the number of voxels along each array axis. The implications of this requirement in terms of the measurement is that the intensity fringe pattern of each Bragg peak should be oversampled. This mirrors the well-known sampling condition for single-peak BCDI measurements.

A remaining requirement of the forward model, provided that our aim is to use modern optimization for its inversion, is that Equations (1) and (2) be implemented with readily differentiable operators. In order to achieve this efficiently, we used Fourier interpolation based operations akin to those in the image

processing community$^{18-20}$ to implement the re-sampling and rotation $T_i$. More details regarding the operations we used are given in the Supplementary Information.

## The MR-BCDI inverse problem

Having established a forward model for multi-reflection BCDI, establishing an inverse problem becomes a matter of expressing a loss function and implementing an automatic-differentiation-based optimization algorithm. The loss function we utilize for MR-BCDI, $\mathcal{L}_{\text{multi}}$, represents the mean error in signal amplitude per BCDI data set, aggregated over the number of data sets:

$$ \mathcal{L}_{\text{multi}}[A, \mathbf{u}] = \sum_{i=1}^{M} \frac{1}{N_{\text{vox}}} \sum_{\mathbf{n}} \left[ \left\| \mathcal{F} \left\{ \chi_i A(T_i \mathbf{m}) e^{i2\pi \mathbf{G}^T \mathbf{u}(T, \mathbf{m})} \right\} \right\|_{\mathbf{n}} - \sqrt{l_{\mathbf{n}}^{(i)}} \right]^2 \quad (3) $$

The quantities $A$, $\mathbf{u}$, and the global scaling factors $\chi_i$ are treated as parameters to be optimized. In this expression, $\chi$ represents the set of $M$ scaling factors, and $l_{\mathbf{n}}^{(i)}$ is the measured diffraction data corresponding to the i'th Bragg reflection. The scaling factors $\chi$ are needed to reconcile the different diffraction signal strengths encountered at different Bragg peaks with the constraint that $A \in [0, 1]$. $N_{\text{vox}}$ is the number of voxels in a single BCDI data set, is $128^3$ in all of the examples in this paper. We note that by using $\sqrt{1}$, the MR-BCDI loss function mirrors what is typically used in coherent diffraction phase retrieval approaches$^{21,22}$.

It is important to consider the degree to which the above optimization problem is overdetermined. This can be established by taking $V$ to have an edge size of $n$ voxels within a buffered array of size $N$ (with $n < N/2$). Combining $\chi$, $A$, and the three scalar components of $\mathbf{u}$ yields $M + 4n^3$ unknowns. The number of measurements (total pixels in the data constraint space) is given by $MN^3$, for $M$ BCDI scans. The MR-BCDI problem is therefore overdetermined by a factor of $MN^3/(M + 4n^3)$. In the case of our first numerical test presented below, $M=5$, $N=128$ and $n=46$. This gives an overdetermination factor of ~27 which is significantly higher than the minimum required 1. We note further that we require $M \ge 3$ Bragg peak measurements with non-coplanar $\mathbf{G}_i$ reciprocal lattice vectors to ensure that the optimization is overdetermined in $\mathbf{u}(\mathbf{x})$, as has been established in prior MR-BCDI literature$^{11,13,14,23}$.

It is convenient for optimization to implement the constraint $A \in [0, 1]$ analytically by defining an activation field $\alpha$ in order to reconcile the fact that the relative scaling of the expected values of $A$ and the values of $\mathbf{u}$ may take on very different values and variances, which can complicate optimization. Thus, the scalar field $\alpha(\mathbf{x})$ is defined such that:

$$ A(\mathbf{x}) = \frac{1}{2} \left[ 1 + \tanh \left( \frac{\alpha(\mathbf{x})}{a_0} \right) \right] \quad (4) $$

where the constant $a_0$ scales the rate of gradient descent of $A$. This allows $\alpha$ to take on any value but will result in $A$ being bracketed within the interval $[0, 1]$. With this definition, we see that $A \approx 0$ for $\alpha(\mathbf{x}) \ll 0$ and $A \approx 1$ for $\alpha(\mathbf{x}) \gg 1$. The transition between 0 and 1 in the neighborhood of $\alpha(\mathbf{x}) = 0$ is controlled by the positive hyperparameter $a_0$. Qualitatively, a larger value of $a_0$ results in a smaller gradient descent step for $A$ at $\mathbf{x}$ (as is seen by differentiating Equation (4) with respect to $\alpha$), and a more gradual and controlled approach of $A$ to its intended value. Smaller $a_0$ results in $A$ converging through sporadic switching between values close to 0 and 1, which is much less controlled. We found that $a_0 = 1$ best served our global convergence rate, and have used this for all reconstructions in this paper. Given this analytical framing of the constraint, the optimization problem to be solved in practice is:

$$ \begin{gathered} a^*, \mathbf{u}^*, \mathbf{X}^* = \arg \min_{\mathbf{a}, \mathbf{u}, \mathbf{X}} \mathcal{L}_{\text{multi}}[A[\mathbf{a}], \mathbf{u}, \mathbf{X}] \quad (5) \\ A^* = A[a^*] \end{gathered} $$
where the * indicates optimal values that produce minimum error.

The loss minimization problem framed above can be solved by implementing an optimization algorithm such as stochastic gradient descent (SCG)²⁴,²⁵. Modern machine learning software toolkits offer the opportunity to do so efficiently and with relative ease. In each iteration of a gradient descent approach, the gradient of the loss function is evaluated with respect to the quantities of interest, and the quantities are updated via an increment proportional to the negative of the gradient. Thus, the total error is expected to be lower for the next iteration, and new gradients and updates are calculated and applied iteratively until improvements in error plateau. Stochastic gradient descent, which is applied in this work, applies gradient descent to different changing subsets of the constraint set. This strategy is valuable in cases where the starting guesses of the quantities of interest are likely far from the eventual solutions, and coarse updates based on gradients with respect to subsets of the full constraint space prove effective in reducing error in early iterations. In our work, $a^*, u^*, x^*$ are updated first by constraining to two of the available Bragg peaks data sets. Then a gradually stricter global constraint to the data is imposed as the reconstruction progresses by constraining to an increasing number of Bragg peaks until the complete data set is used as the basis of the update.

Owing to the complexity of our MR-BCDI model, we implement SCG using automatic differentiation (AD) using the PyTorch toolkit, which simplifies evaluation of the loss function and its gradient. AD-based optimization evaluates the gradient automatically and exactly via the chain rule applied repeatedly to the constituent elementary functions. This feature allows increasingly complicated forward models of coherent diffraction to be developed corresponding to more sophisticated measurements, enabling image inversion without the need to derive loss function gradients analytically. Given the complexity of the forward model presented here, MR-BCDI leverages this ability. The use of AD optimization does require that the mathematical operations composing the forward model be differentiable, as we have taken care to enforce here. In this work, we implement the “Adam” optimization algorithm, which is an extension of SCG that is commonly used in machine learning, in a Python-based MR-BCDI reconstruction software²⁶ for CPU and GPU hardware using PyTorch²⁷.

Lastly, we note that the optimal solution in Equation (5) does not preclude the universal twin solution, akin to the twin image degeneracy in single-reflection BCDI²⁸. In other words, the crystal defined by the transformations $A(x) \to A(-x)$ and $u(x) \to -u(-x)$ in the same lab frame is also a solution to the MR-BCDI problem, in that it results in indistinguishable diffraction patterns for the same diffraction geometry. Thus, no matter its implementation, MR-BCDI is still susceptible to a reconstruction degeneracy, but our global optimization will converge to one twin solution or the other without requiring intervention to align individual phase retrieval instances¹². Typical ambiguities associated with real-space coordinate origin and phase offset present in single-Bragg-peak BCDI are mitigated in our approach because a single object is the subject of iterative update from all Bragg peaks. This forces consistency of the spatial position of the object. Also, because our reconstruction problem aims to image lattice distortion rather than the phase of the exit wavefield, phase offsets are not an issue.

## Method demonstration

In this section we present three reconstruction results with the MR-BCDI technique described above. Two of these use numerically synthesized coherent Bragg diffraction intensity data sets from digital nanocrystals. In the numerical studies, one crystal has a slowly varying continuous inhomogeneous distortion field and the other has a discontinuous, winding displacement field, which emulates two orthogonal screw dislocations within the crystal. The

third reconstruction is from experimental diffraction data consisting of six Bragg reflections of a silicon carbide (SiC) nanocrystal measured using BCDI at Beamline 34-ID-C of the Advanced Photon Source. In each case, we present the MR-BCDI reconstruction of the scattering amplitude $A(x)$ and the lattice deformation $u(x)$. All reconstructions were performed using a single Nvidia Tesla P100 GPU with 16 GB of RAM, and typical reconstruction times were 2-4 hours.

We have adopted the matrix coordinate convention for all the cross-section plots in this paper. As an example, a cross-section image labeled as $X-Y$ implies that the $X$-axis is directed from top to bottom along the image, and the $Y$-axis from left to right. The third axis (in this example, the $Z$-axis) emerges out of the plane of the figure in order to maintain right-handedness. The $X$, $Y$, and $Z$ axes correspond to the laboratory frame of reference (i.e., the $\hat{s}_1$, $\hat{s}_2$ and $\hat{s}_3$ directions in Fig. 1, respectively). With this convention, we ensure that the reference frame of the reconstructed crystal is consistent with the programmed order of the array axes in the multi-dimensional FFT routines in Python and PyTorch.

**Simulated crystal without dislocations.** A synthetic crystal with arbitrary facets was generated on a Cartesian grid of size $128 \times 128 \times 128$ voxels, with a voxel size of $s_0 = 12$ nm along each axis in the laboratory frame ($\hat{s}_1 \hat{s}_2 \hat{s}_3$) in Fig. 1). To emulate realistic BCDI scans at a diverse set of sample and detector orientations, a face-centered cubic (FCC) gold lattice with a lattice constant of $a_0 = 4.078$ Å was assumed for the simulated particle. The orientation of the FCC unit cell with respect to the lab frame was assigned arbitrarily. A slowly-varying internal field $u(x)$ was created by first generating uniform random samples for the three components of $u(x)$ ($-0.1a_0 \le u(x) \le 0.1a_0$) and then retaining the long-period variations by convolving each $u_i(x)$ with a low-pass filter. The crystal spanned (39, 39, 40) voxels along the laboratory frame axes, and a cubic bounding box $V$ of size $46 \times 46 \times 46$ voxels was chosen for the reconstruction.

Numerical MR-BCDI diffraction data were generated in a manner consistent with the measurement of a real FCC-structured crystal at the BCDI diffractometer at the APS Sector 34-ID-C beamline. The angular rotations available at that diffractometer (Fig. 1) were emulated as coordinate rotations about cardinal lab frame axes as follows: the sample rotation stage $\theta$ and the detector rotation stage $\delta$ were coordinate transformations about $\hat{s}_2$, and the sample rotation stage $\phi$ and the detector rotation stage $\gamma$ are coordinate transformations about $\hat{s}_1$. With this convention, the orientations of $\theta$, $\delta$, and $\gamma$ needed to satisfy the [111], [111], [220], [202] and [022] Bragg reflections with an x-ray energy of 9 keV were determined (with $\phi$ held at 0) and used to make up the composite rotation $R_i$ as discussed in detail in²⁶. Since either the $\theta$ or $\phi$ sample rotation axes could be used to emulate the angular rocking curve for each Bragg peak, the rotation axis that provided a higher degree of orthogonality of the basis set [q,q,q,k] was chosen (details in the Supplementary Information). Following the forward model described above, $\psi_m^{(i)}$ and $\psi_n^{(i)}$ were calculated for each Bragg peak. To create a simulated MR-BCDI data set, the $\psi_n^{(i)}$ were squared to generate the far-field intensity distributions, the peak intensities for each data set were scaled to $10^5$ photon counts in the highest-intensity pixel (consistent with high-quality experimental BCDI data), and Poissonian fluctuations consistent with the scaled intensity fields were added to emulate counting statistics.

For the MR-BCDI reconstruction, the crystal was initialized to a constant-amplitude cube ($A = 1$) occupying the entire $46 \times 46 \times 46$ bounding box with no interior lattice distortion ($u = 0$). This was achieved by setting $a(x)$ from Equation (4) to 2 and $a_0 = 1$ for all pixels within $V$. Each of the global scaling factors $\chi_i$ was initialized to match the total energy $\sum_n ||\Psi_n||^2$ of the corresponding simulated Bragg intensity pattern. An Adam optimizer²⁹ with an initial learning rate of 0.005 was employed.
<table><caption>**Table 1.** Sequential optimization plan for dislocation-free crystal.</caption><thead><tr><th>Epoch</th><th>Number of minibatches</th><th>Scans in minibatch</th><th>Iterations per minibatch</th></tr></thead><tbody><tr><td>1</td><td>400</td><td>2</td><td>6</td></tr><tr><td>2</td><td>320</td><td>3</td><td>12</td></tr><tr><td>3</td><td>240</td><td>3</td><td>25</td></tr><tr><td>4</td><td>160</td><td>3</td><td>50</td></tr><tr><td>5</td><td>80</td><td>3</td><td>100</td></tr><tr><td>6</td><td>40</td><td>4</td><td>200</td></tr><tr><td>7</td><td>20</td><td>4</td><td>400</td></tr><tr><td>8</td><td>1</td><td>5 (full dataset)</td><td>2000</td></tr></tbody></table>

Note: The BCDI scans within each minibatch were randomly selected.

A total variation (TV) regularizer was applied to $a(x)$ with a coupling constant of $10^{-5}$. During optimization, $u(x)$ was constrained so as not to extend beyond the magnitude of the lattice plane separations $\pm a_0/2$ in each of the [100], [010] and [001] crystallographic directions of the model. This constraint may be employed without loss of generality and is discussed further in the Supplementary Information. The optimization was carried out sequentially over increasingly larger 'minibatches' of the BCDI data, as summarized in Table 1 of Methods. After each optimization epoch, a median filter with a kernel size of $3 \times 3 \times 3$ voxels was applied to the reconstructed lattice distortion field $u(x)$ in order to remove spurious isolated discontinuities induced by phase wraps. In the final epoch the full set of 5 Bragg peaks was used, and for the last 1000 iterations of this epoch, the

<figure><img src="image_4.png" /><figcaption>**Fig. 2** MR-BCDI results from simulated crystal. a Ground truth (top) and reconstructed (bottom) cross-sections of the electron density $A(x)$ of the synthetic strained crystal with no internal dislocations. The array size represents the original bounding box chosen for reconstruction, i.e., 46 pixels along each lab frame axis. The spatial resolution was estimated by fitting an error function to the reconstructed electron density along the dashed line (see Fig. 3c). b Orthogonal cross sections of ground truth lattice distortion components $u_1(x)$, $u_2(x)$ and $u_3(x)$. The length scale is the same as in a. c Corresponding sections of the reconstructed lattice distortion. The scales shown for $u$ in this and subsequent Figures are in units of fractional lattice constant [$-0.05a_0$, $0.05a_0$]. For comparison, this scale translates to a scale of [$-0.1$, $0.1$] radians when considering the $G^T u$ scalar field projection corresponding to the [002] Bragg peak.</figcaption></figure>
voxels in the crystal for which $A > 0.2$ were optimized with the other voxels held constant to refine the solution within the crystal interior.

Figure 2 shows a comparison of the ground truth and the reconstructed quantities $A(\mathbf{x})$ and $\mathbf{u}(\mathbf{x})$ on the lab frame grid after the optimization schedule was completed. Orthogonal cross-sections of the three components of $\mathbf{u}(\mathbf{x})$ are shown from the reconstruction alongside those from the known model (ground truth). The degree of similarity between the ground truth and reconstructed $\mathbf{u}(\mathbf{x})$ quantities is evidenced by comparing the histograms of the voxel-by-voxel components of $\mathbf{u}(\mathbf{x})$

(Fig. 3b). Point-to-point residuals were also calculated for all voxels with $A > 0.5$, and the histogram of these residuals is sharply peaked at $u = 0$, indicating good agreement between reconstruction and ground truth. The second peak in residuals originates from imperfect spatial alignment of the estimated object and ground truth, which results in higher residuals near the object edges. The spatial resolution of the reconstruction was estimated by fitting an error function to the density profile $A(\mathbf{x})$ along the dashed line shown in Fig. 3c. This error function has a characteristic width of $\sigma = 0.788$ pixels, which gives a spatial resolution estimate of $\sqrt{8 \ln 2\sigma} \approx 1.86$ pixels. Figure 4 shows the

<figure><img src="image_5.png" /><figcaption>**Fig. 3** Figures of merit for simulated crystal reconstruction. a Trend in multi-reflection loss function given by Equation (3) for the dislocation-free synthetic crystal. The units of this loss function are chosen not to mirror a physical quantity, but as result of a noise model robust to Poissonian noise that is routinely used in BCDI phase retrieval approaches. The red lines denote the beginning of a new optimization epoch. b Histogram of the simulated and reconstructed components of the vector $\mathbf{u}(\mathbf{x})$, along with point-to-point residuals of each vector component. c Error function fit to the reconstructed profile along the dashed line in Fig. 2a (X-Y slice). The spatial resolution was estimated from this fit to be ≈ 1.86 pixels, or ≈ 23 nm.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 4** Comparison of simulated and reconstructed diffraction patterns. Cross-sections of the noisy simulated (top) and reconstructed (bottom) coherent diffraction patterns from the five Bragg peaks chosen for the dislocation-free synthetic crystal. The reconstructed diffraction patterns are shown here without the scaling effect of the $\chi_\mu$ but rather scaled to match the corresponding simulated intensity peaks. To fully match the color scales, each reconstructed diffraction was clipped below to the smallest nonzero photon count in the corresponding signal (i.e., 1 photon).</figcaption></figure>
<table><caption>**Table 2.** Optimization plan for crystal with orthogonal screw dislocations along [111] and [220].</caption><thead><tr><th>Epoch</th><th>Number of minibatches</th><th>Scans in minibatch</th><th>Iterations per minibatch</th></tr></thead><tbody><tr><td>1</td><td>800</td><td>2</td><td>3</td></tr><tr><td>2</td><td>400</td><td>2</td><td>6</td></tr><tr><td>3</td><td>320</td><td>2</td><td>12</td></tr><tr><td>4</td><td>240</td><td>3</td><td>25</td></tr><tr><td>5</td><td>160</td><td>3</td><td>50</td></tr><tr><td>6</td><td>80</td><td>3</td><td>100</td></tr><tr><td>7</td><td>40</td><td>3</td><td>200</td></tr><tr><td>8</td><td>20</td><td>3</td><td>400</td></tr><tr><td>9</td><td>1</td><td>4 (full dataset)</td><td>2000</td></tr></tbody></table>

Note: The BCDI scans within each minibatch were randomly selected.

comparison of two-dimensional cuts of the reconstructed diffraction patterns at the end of the last epoch alongside the ground truth simulated diffraction used as constraints, showing good agreement.

The progression of the loss function $\mathcal{L}_{\text{multi}}[\mathcal{A}, \mathbf{u}]$ over the entire optimization process is shown in Fig. 3a and has several noteworthy features. The dashed lines in Fig. 3a indicate the beginning of each new optimization epoch. Within each epoch, the rapid oscillations are attributed to the abrupt change in the loss function landscape due to each new randomized minibatch. For example, in the first epoch, the loss function is successively optimized over 400 sets of 2 randomly selected BCDI scans for 6 iterations each. Whenever a new randomized set of scans is used, the loss function landscape abruptly changes from its previous state, and the gradient descent evaluates a different error.

We also performed single-Bragg-peak phase retrieval of the 5 BCDI simulated data sets using standard methods$^{3,4}$ in order to demonstrate the advantage of our MR-BCDI method, and the results are shown in Supplementary Fig. 2. It is apparent that the MR-BCDI reconstructed electron density from Fig. 2a is more uniform than what was obtained from phase retrieval of individual peaks, indicating that applying the global set of constraints to a single $\mathcal{A}(\mathbf{x})$ mitigates spurious fluctuations that are difficult to avoid in single-peak BCDI reconstructions. A deterioration of spatial resolution was also observed with the single-peak reconstructions (3.06 pixels resolution estimated by error function fitting), again stemming from the fact that single peak reconstructions do not benefit from a coupling across related data.

**Simulated crystal with screw dislocations.** A second synthetic crystal with arbitrary facets was generated on the same lab-frame grid as above ($s_0 = 12$ nm), and as before, the crystallographic lattice structure of gold was adopted with a new initial lab frame crystal orientation. Using this crystal lattice orientation frame, a $\mathbf{u}(\mathbf{x})$ field was calculated to simulate two spatially separated screw dislocations within the crystal volume, with dislocation lines along the orthogonal [111] and [220] crystallographic directions. For this simulation, the respective Burgers vector magnitudes $||\mathbf{b}||$ for the screw dislocations were set to: $||\mathbf{b}||_{111} = a_0/\sqrt{3}$ and $||\mathbf{b}||_{220} = a_0/\sqrt{8}$. The $\mathbf{u}(\mathbf{x})$ fields for each dislocation were calculated with the continuum model of lattice distortion for screw dislocations, as in other BCDI work$^{30,31}$, and the net displacement field was modeled as the vector sum of the displacement fields of the pair of screw dislocations. The simulation procedure of calculating the BCDI data sets, including the selection criterion for Bragg reflections, selection of rocking

directions, and introduction of Poisson noise, followed the same procedure as for the dislocation-free crystal. Four Bragg peaks with high reciprocal-space mutual orthogonality were chosen for modeling BCDI intensity patterns: [200], [002], [202] and [220]. The simulation box $V$ was chosen to be $40 \times 40 \times 40$ voxels in size. The Adam optimizer was initialized with a learning rate of 0.02, the object initialization was done as in the previous example, and a similar optimization schedule was used (see Table 2 in Methods). The final optimization stretch for the interior voxels of the object was carried out for 5000 iterations constrained by all four Bragg peak data sets. A median filter of size $7 \times 7 \times 7$ voxels was applied to the components of $\mathbf{u}(\mathbf{x})$ after each optimization epoch.

Figure 5 compares the simulated and reconstructed electron densities $\mathcal{A}(\mathbf{x})$ and $\mathbf{u}(\mathbf{x})$, and Fig. 6 compares Bragg peak intensity distributions. (Reconstruction metrics are shown in Supplementary Fig. 4). As in the dislocation-free case, there is good agreement between the simulation and reconstructed images as well as in the diffraction intensity patterns. One noteworthy feature of this reconstruction is the presence of regions of low electron density where the discontinuous lattice distortion in the $\mathbf{u}(\mathbf{x})$ field intersects the plane of the figure (regions indicated with blue arrows) and the surface of the crystal (regions circled in red). This effect is attributed to the fact that the highest spatial resolution components of the signal (the “high-q” regions of the detector) are suppressed with the introduction of Poisson counting statistics into the simulated BCDI data sets, effectively imposing a low-pass spatial filter on the reconstructed image.

As before, a comparison was also made to BCDI images reconstructed from the individual Bragg peak data sets from the dislocated crystal lattice using convention phase retrieval approaches. The recipe used was the same for the dislocation-free crystal. A representative amplitude cross-section of a single-peak BCDI image from the [202] Bragg peak is shown in Fig. 7. From this image, it is clear by comparing to Fig. 5 that the conventional phase retrieval approach from a single peak struggles to converge to the uniform amplitude distribution expected in the ground truth.

**Silicon carbide nanocrystal.** Here, we describe a MR-BCDI reconstruction from experimental data acquired from a silicon carbide (SiC) nanoparticle fabricated from a bulk substrate in a manner similar to fabrication methods of SiC quantum sensors$^{32-37}$. The sample was fabricated using lithography and wet etching methods from a bulk single crystal of single-polymorph SiC with a hexagonal 4H structure (lattice parameters $a = 3.073$ Å and $c = 10.053$ Å). (The fabrication process is described in the Methods Section.) The shape of the sample was intentionally chosen to be an asymmetric “D”-shaped column with tapering edges fabricated such that the macroscopic facets approximately aligned to low-index lattice planes in the underlying crystal, as shown schematically in Fig. 8a. The sample that was studied is shown in the scanning electron micrograph in Fig. 8b, oriented with the edge of the “D” shape laying flat on the silicon substrate that was used to support the crystal during BCDI measurement. As an intentional result of the fabrication process, only the region of the sample inside the red circle retained a highly crystalline structure and contributed to the measured BCDI peaks.

The beamline at 34-ID-C of the Advanced Photon Source was used to measure BCDI data from six Bragg reflections from this crystal: [1011], [1010], [0112], [0110], [0111] and [1013]. At each Bragg condition, rocking curves were performed via fine angular steps of the $\theta$ motor in order to record 3D coherent diffraction intensity data sets with a pixelated area detector, resulting in raw BCDI data sets of $256 \times 256$ images with 80 steps along the rocking curve. In terms of signal strength, the maximum pixel intensities in each of the 6 measured Bragg peaks were respectively: 40,980, 16,072, 32,661, 14,577, 34,123 and 7805 peak photons.
<figure><img src="image_7.png" /><figcaption>**Fig. 5** MR-BCDI results from simulated crystal with dislocations. **a** Ground truth (top) and reconstructed (bottom) electron density cross-sections of the synthetic strained crystal with two orthogonal screw dislocations. The array size represents the original bounding box chosen for reconstruction, i.e., 40 pixels along each lab frame axis. The arrows indicate the points of intersection of the dislocation cores with the figure plane. **b** Orthogonal cross sections of ground truth lattice distortion components $u_1(x)$, $u_2(x)$ and $u_3(x)$. **c** Corresponding cross section for the reconstructed lattice distortion. The length scale is the same as **a**. The regions of discontinuity in the lattice distortion match with the drop in electron density in **a**. The arrows show the location of the modeled dislocation cores. The circled regions show imperfect reconstruction at locations where the lattice discontinuity meets the crystal surface.</figcaption></figure>

Several steps were taken to process this data for MR-BCDI reconstruction. The data sets were clipped to a size of 128 × 128 pixels in the detector plane, and the angular rocking curve dimensions of the data sets were zero-padded to a size of 128, resulting in a cubic data array space of 128 × 128 × 128 voxels for each Bragg peak. This uniform array size ensured that the $T_i$ transformation operations act as intended. The lab frame crystal orientation of the SiC particle was determined by at least squares optimization of the crystal orientation (expressed as a 3-parameter rotation vector) over the $\delta, \gamma$ and $\theta$ motor orientations of all 6 BCDI scans with the Nelder-Mead optimizer³⁸. From this crystal orientation, the set of appropriately oriented reciprocal space lattice vectors $G_i$ that are needed for MR-BCDI reconstruction were generated.

The MR-BCDI reconstruction was performed on the data from the SiC particle in much the same way as for the simulated data. Figure 8c shows the measured and inferred diffraction patterns, showing close agreement, as in the previous cases. Figure 9a depicts the progression of the loss function for the optimization scheme described in Table 3 of the Supplementary Information. Figure 9b shows the electron density after optimization of the interior voxels with $A > 0.1$. Also shown is the contour at $A = 0.65$, depicting the approximate surface of the crystal. Figure 9c shows the reconstructed components of $u(x)$ within this contour. Figure 10 shows the corresponding 3D isosurface plots of the components of $u$, with the color scale denoting the spatially varying lattice distortion in nanometers. We see from the absence of discontinuities in the $u$ field that the particle contains a
<figure><img src="image_8.png" /><figcaption>**Fig. 6** Comparison of simulated and reconstructed diffraction patterns. Cross-sections of the noisy simulated (top) and reconstructed (bottom) coherent diffraction patterns from the five Bragg peaks chosen for the gold crystal with two orthogonal screw dislocations.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 7** Reconstructed amplitude from a single Bragg peak of the simulated data set. a An individual amplitude reconstruction $||\psi^{(0)}||$ from phase retrieval applied to the BCDI virtual data set of the [202] Bragg peak of the crystal with screw dislocations, using the recipe in Section B of the Supplementary Information. The cross section shown displays prominent fluctuations in electron density that are not present in the MR-BCDI reconstructions by virtue of the coupling effect of multiple data sets in a single optimization of the object. Note that the aspect ratio of the central slice of the object is not as it appears in Fig. 5a because the output of standard single-peak BCDI phase retrieval is not shear corrected.</figcaption></figure>

smoothly varying lattice displacement field and no dislocations. Further, morphological features of the reconstructed image correspond as expected with the features of the SiC crystal, including the "D"-shaped cross-section of the nanoparticle, the particle aspect ratio, and the tapering edges.

## DISCUSSION

In this paper we demonstrated a method by which to reconstruct vector-valued lattice distortion fields within nanocrystals by optimization of a comprehensive forward model of multi-Bragg-peak BCDI data. This forward model flexibly accounts for important geometric factors that arise when making BCDI measurements, is amenable to efficient inversion with modern optimization toolkits, and allows for globally constraining a single image reconstruction to multiple Bragg peak measurements. A key feature of our method is the fact that the forward model we developed to reconcile multiple Bragg peak measurements from a single originating crystal does so by emulating the native measurement space associated with BCDI rocking curve

data sets. We also note that our reconstruction approach obviates the need for trial-and-error fixed point projection recipes and reduces the need for intermittent support updates via the shrink-wrap algorithm⁴. Our formulation of MR-BCDI offers another potential advantage in that it permits the analytical interrogation of the high-dimensional solution space in the neighborhood of the true solution, opening the door to rigorous uncertainty analysis where fixed-point-iteration-based methods fall short.

In implementing our MR-BCDI method, we have demonstrated successful reconstructions of the 3D electron density and displacement field for both numerical and experimental data. Crucially, we have shown that crystals with both smoothly varying displacement fields as well as ones with discontinuities in lattice displacement due to, for example, the presence of dislocations, are reconstructable with our approach. In addition, this work brings coherent diffraction methods further into the scope of global optimization techniques¹⁵,³⁹⁻⁴², thereby enabling the use of highly optimized software packages capable of handling large data volumes and running on seemingly ever-improving high-performance computing hardware. The use of AD within a machine learning framework also makes the development and testing of more sophisticated cost functions that make use of regularizers based on a-priori information for potentially improved MR-BCDI reconstructions in certain situations.

The demonstration of lattice displacement reconstruction within an individual SiC nanoparticle motivates the use of MR-BCDI in science domains where structural inhomogeneities at nanometric length scales impact materials performance and properties. For example, the case of SiC demonstrated in this work is pertinent to the field of quantum sensing. In a SiC nanoparticle quantum sensor, the sensitivity of near-surface optically active point defects to changes in temperature, magnetic field, and mechanical stress is exploited in order to measure these quantities in media where the nanoparticles can be dispersed³²⁻³⁷. Latent inhomogeneous lattice distortion fields within these crystal sensors detract from the precision of such measurements. Thus, quantum sensor fabrication methods aim to minimize such lattice distortion fields. As we demonstrated here, MR-BCDI provides a means by which to assess this key metric.

## METHODS

### Sequential optimization plans

In Tables 1–3 we specify the optimization plans used in the three reconstructions corresponding to the simulated crystal without
<figure><img src="image_10.png" /><figcaption>**Fig. 8** Silicon carbide nanoparticle used for experimental demonstration of MR-BCDI. **a** Schematic of the design of the lithographically fabricated SiC nanoparticle with facets aligned closely with crystallographic directions of the parent 4H single crystal substrate. **b** A scanning electron micrograph of the nanoparticle used in the study. The portion of the particle that retained a high degree of crystallinity from which BCDI scans were done is circled in red. **c** A central slice through the six Bragg peaks measured from this crystal (top row), alongside this same slice after MR-BCDI optimization and image reconstruction (bottom row).</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 9** Loss function and interior lattice displacement of SiC particle reconstruction. **a** Loss function trend for the SiC crystal reconstruction, for the optimization plan shown in Table 3 of the Supplementary Information. **b** Final estimated SiC nanocrystal electron density $A$, from a highly permissive mask threshold of 0.1. The red line is the contour at a threshold of $A = 0.65$, representing the approximate crystal surface. **c** Reconstructed components of the lattice distortion within the SiC crystal (centered at zero volume-averaged distortion in each dimension), within the estimated crystal surface at $A = 0.65$.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**Fig. 10** Surface lattice displacement of SiC. a-c Isosurface plots of the reconstructed $u_1(x)$, $u_2(x)$ and $u_3(x)$ respectively in the SiC crystal (color scale in nanometers), with an isovalue of 0.65 in electron density.</figcaption></figure>

<table><caption>**Table 3.** Optimization plan for the silicon carbide nanocrystal.</caption><thead><tr><th>Epoch</th><th>Number of minibatches</th><th>Scans in minibatch</th><th>Iterations per minibatch</th></tr></thead><tbody><tr><td>1</td><td>800</td><td>4</td><td>6</td></tr><tr><td>2</td><td>640</td><td>4</td><td>12</td></tr><tr><td>3</td><td>480</td><td>4</td><td>25</td></tr><tr><td>4</td><td>320</td><td>5</td><td>50</td></tr><tr><td>5</td><td>160</td><td>5</td><td>100</td></tr><tr><td>6</td><td>80</td><td>5</td><td>200</td></tr><tr><td>7</td><td>40</td><td>5</td><td>400</td></tr><tr><td>8</td><td>1</td><td>6 (full dataset)</td><td>2000</td></tr></tbody></table>

Note: The BCDI scans within each minibatch were randomly selected.

dislocations (Table 1), the simulated crystal with dislocations (Table 2), and the experimental data set from a nanoparticle of SiC (Table 3).

## Silicon carbide nanoparticle fabrication

Here, we specify the method of fabrication of the SiC crystal used in the experimental demonstration of our MR-BCDI approach. The SiC nanoparticles were fabricated from a single 4H-SiC wafer with an *i-p-n* doping structure. The wafer had a 400nm-thick intrinsic SiC layer with an intrinsic $<10^{15}$atoms/cm$^3$ defect content, a 2µm-thick *p*-type layer with $<10^{19}$atoms/cm$^3$ of aluminum dopant atoms, and a 0.5mm-thick *n*-type layer with $<10^{18}$atoms/cm$^3$ of nitrogen dopant atoms. The wafer was diced (5 × 5 mm), electron-beam lithographically patterned (bi-layer PMMA 495K-A3/950K-A6; 1000µC/cm$^2$ dose, 3-minute development in 1:3 MIBK:IPA at room temperature), and metallized (10 nm of Ti and 60 nm of Cr via e-beam evaporation) to transfer the D-shaped array pattern onto the SiC. The D shape flat of the lithographic pattern was aligned to the crystallographic [1120] direction (see Fig. 8 in the main text). To fabricate particles that adopt the profile of the surface pattern, the SiC was etched (SF6-Ar ICP/RIE dry-etching) with the Ti/Cr acting as a hard mask. The 1.2 µm etch depth exposed the *p*-type layer for subsequent photo-electrochemical (PEC) etching of the nano-particles. The etch-defined nanopillars were then PEC etched (0.2M KOH, -0.3V biased etch under 365nm UV illumination at < 500 mW) to partially and selectively remove the *p*-type SiC to allow for nanoparticle formation and subsequent detachment. For this particular experimental run, the *p*-type SiC was underetched, resulting in a weakly scattering, porous *p*-type tether remaining which did not contribute to the BCDI signal (as shown in Fig. 8b in the main text). These nanoparticle arrays were stamp-transferred to a PMMA-coated silicon substrate. The wafers were baked at the PMMA glass transition temperature, which

locked the nanoparticles on the Si substrate. The nanoparticles were arranged on the Si wafer in 100 µm pitched arrays bordered by macroscopic fiducials consisting of tightly packed nanoparticles. Finally, the PMMA was O$_2$ etched away, leaving pristine SiC nanoparticle arrays on the Si substrate that were then conformally covered with alumina (22 nm thick) via ALD to help adhere them to the substrate.

## DATA AVAILABILITY

The simulated and experimental BCDI data sets, along with the reconstruction results, are publicly available in HDF5 format in the *siddharth-maddali/mrbcdi* repository in Github, through the Git Large File System (LFS).

## CODE AVAILABILITY

The MR-BCDI Python module, along with the Jupyter notebooks for reconstructions, are publicly available in the *siddharth-maddali/mrbcdi* repository in Github. The software for conventional phase retrieval is publicly available in the *siddharth-maddali/Phaser* repository in Github.

Received: 2 August 2022; Accepted: 15 April 2023; Published online: 23 May 2023

## REFERENCES

1. Robinson, I. K., Vartanyants, I. A., Williams, G. J., Pfeifer, M. A. & Pitney, J. A. Reconstruction of the shapes of gold nanocrystals using coherent x-ray diffraction. *Phys. Rev. Lett.* **87**, 195505 (2001).

2. Miao, J., Ishikawa, T., Robinson, I. K. & Murnane, M. M. Beyond crystallography: Diffractive imaging using coherent x-ray light sources. *Science* **348**, 530–535 (2015).

3. Fienup, J. R. Phase retrieval algorithms: a comparison. *Appl. Opt.* **21**, 2758–2769 (1982).

4. Marchesini, S. et al. X-ray image reconstruction from a diffraction pattern alone. *Phys. Rev. B* **68**, 140101 (2003).

5. Marchesini, S. Invited article: A unified evaluation of iterative projection algorithms for phase retrieval. *Rev. Sci. Instrum.* **78**, 011301 (2007).

6. Kawaguchi, T. et al. Gas-induced segregation in pt-rh alloy nanoparticles observed by in situ bragg coherent diffraction imaging. *Phys. Rev. Lett.* **123**, 246001 (2019).

7. Singer, A. et al. Nucleation of dislocations and their dynamics in layered oxide cathode materials during battery charging. *Nat. Energy* **3**, 641–647 (2018).

8. Hruszkewycz, S. O. et al. Strain annealing of sic nanoparticles revealed through bragg coherent diffraction imaging for quantum technologies. *Phys. Rev. Mater.* **2**, 086001 (2018).

9. Pateras, A. et al. Combining Laue diffraction with Bragg coherent diffraction imaging at 34-ID-C. *J. Synchrotron Radiat.* **27**, 1430 (2020).

10. Newton, M. C., Leake, S. J., Harder, R. & Robinson, I. K. Three-dimensional imaging of strain in a single ZnO nanorod. *Nat. Mater.* **9**, 120 (2009).

11. Hofmann, F. et al. Nanoscale imaging of the full strain tensor of specific dislocations extracted from a bulk sample. *Phys. Rev. Mater.* **4**, 013801 (2020).

12. Wilkin, M. J. et al. Experimental demonstration of coupled multi-peak bragg coherent diffraction imaging with genetic algorithms. *Phys. Rev. B* **103**, 214103 (2021).
13. Newton, M. C. Concurrent phase retrieval for imaging strain in nanocrystals. *Phys. Rev. B* **102**, 014104 (2020).

14. Gao, Y., Huang, X., Yan, H. & Williams, G. J. Bragg coherent diffraction imaging by simultaneous reconstruction of multiple diffraction peaks. *Phys. Rev. B* **103**, 014102 (2021).

15. Kandel, S. et al. Using automatic differentiation as a general framework forptychographic reconstruction. *Opt. Express* **27**, 18653–18672 (2019).

16. Maddali, S. et al. General approaches for shear-correcting coordinate transformations in Bragg coherent diffraction imaging. Part I. *J. Appl. Crystallogr.* **53**, 393–403 (2020).

17. Li, P. et al. General approaches for shear-correcting coordinate transformations in Bragg coherent diffraction imaging. Part II. *J. Appl. Crystallogr.* **53**, 404–418 (2020).

18. Unser, M., Thevenaz, P. & Yaroslavsky, L. Convolution-based interpolation for fast, high-quality rotation of images. *IEEE Trans. Image Process* **4**, 1371–1381 (1995).

19. Larkin, K. G., Oldfield, M. A. & Klemm, H. Fast fourier method for the accurate rotation of sampled images. *Opt. Commun.* **139**, 99–106 (1997).

20. Thévenaz, P., Blu, T. & Unser, M. Interpolation revisited [medical images application]. *IEEE Trans. Med. imaging* **19**, 739–758 (2000).

21. Godard, P., Allain, M., Chamard, V. & Rodenburg, J. Noise models for low counting rate coherent diffraction imaging. *Opt. Express* **20**, 25914–25934 (2012).

22. Godard, P. On the use of the scattering amplitude in coherent x-ray bragg diffraction imaging. *J. Appl. Crystallogr.* **54**, 797–802 (2021).

23. Kandel, S. *Using automatic differentiation for coherent diffractive imaging applications*. PhD thesis dissertation, (2021).

24. Saad, D. Online algorithms and stochastic approximations. *Online Learn. J.* **5**, 6–3 (1998).

25. Duchi, J., Hazan, E. & Singer, Y. Adaptive subgradient methods for online learning and stochastic optimization. *J. Mach. Learn. Res.* **12**, 2121–2159 (2011).

26. Maddali, S. *mrbcdi*: Differentiable multi-reflection bragg coherent diffraction imaging for lattice distortion fields in crystals. https://github.com/siddharth-maddali/mrbcdi.

27. Paszke, A. et al. Pytorch: An imperative style, high-performance deep learning library. In *Adv. Neural. Inf. Process. Syst.* (eds H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, & R. Garnett) 32, pages 8024–8035. (Curran Associates, Inc., 2019).

28. Fienup, J. R. & Wackerman, C. C. Phase-retrieval stagnation problems and solutions. *J. Opt. Soc. Am. A* **3**, 1897–1907 (1986).

29. Kingma, D. & Ba, J. Adam: A method for stochastic optimization. *arXiv preprint https://arxiv.org/abs/1412.6980* (2014).

30. Phillips, R. *Crystals, Defects and Microstructures: Modeling Across Scales*. (Cambridge University Press, 2001).

31. Ulvestad, A., Menickelly, M. & Wild, S. M. Accurate, rapid identification of dislocation lines in coherent diffractive imaging via a min-max optimization formulation. *AIP Adv.* **8**, 015114 (2018).

32. Fuchs, F. et al. Silicon carbide light-emitting diode as a prospective room temperature source for single photons. *Sci. Rep.* **3**, 1637 (2013).

33. Lohrmann, A. et al. Single-photon emitting diode in silicon carbide. *Nat. Commun.* **6**, 7783 (2015).

34. Cochrane, C. J., Blacksberg, J., Anders, M. A. & Lenahan, P. M. Vectorized magnetometer for space applications using electrical readout of atomic scale defects in silicon carbide. *Sci. Rep.* **6**, 37077 (2016).

35. Atatuë, M., Englund, D., Vamivakas, N., Lee, S.-Y. & Wrachtrup, J. Material platforms for spin-based photonic quantum technologies. *Nat. Rev. Mater.* **3**, 38–51 (2018).

36. Miao, K. C. et al. Electrically driven optical interferometry with spins in silicon carbide. *Sci. Adv.* **5**, eaay0527 (2019).

37. Bourassa, A. et al. Entanglement and control of single nuclear spins in isotopically engineered silicon carbide. *Nat. Mater.* **19**, 1319–1325 (2020).

38. Nelder, J. A. & Mead, R. A simplex method for function minimization. *Comput J.* **7**, 308–313 (1965).

39. Kandel, S. et al. Efficient ptychographic phase retrieval via a matrix-free levensberg-marquardt algorithm. *Opt. Express* **29**, 23019–23055 (2021).

40. Scheinker, A. & Pokharel, R. Adaptive 3d convolutional neural network-based reconstruction method for 3d coherent diffraction imaging. *J. Appl. Phys.* **128**, 184901 (2020).

41. Cherukara, M. J. et al. AI-enabled high-resolution scanning coherent diffraction imaging. *Appl. Phys. Lett.* **117**, 044103 (2020).

42. Harder, R. Deep neural networks in real-time coherent diffraction imaging. *IUCrJ* **8**, 1–3 (2021).

## ACKNOWLEDGEMENTS

The development of the MR-BCDI forward model and inversion approach, experimental demonstration, and design and fabrication of the SiC nanoparticles was supported by the U.S. Department of Energy (DOE), Office of Science, Basic Energy Sciences, Materials Science and Engineering Division. Additional support for materials preparation came from the Q-NEXT Quantum Center, a U.S. Department of Energy, Office of Science,

National Quantum Information Science Research Center, under Award Number DE-FOA-0002253. Silicon carbide deterministic nanoparticle fabrication and SEM characterization work was performed under proposals 72483 and 775514 in the Center for Nanoscale Materials clean room. Work performed at the Center for Nanoscale Materials, a U.S. Department of Energy Office of Science User Facility, was supported by the U.S. DOE, Office of Basic Energy Sciences, under Contract No. DE-AC02-06CH11357. Refinement of the geometric, computational and optimization concepts was supported by the European Research Council (European Union's Horizon H2020 research and innovation program grant agreement No. 724881). Generation of the simulated structures and the BCDI data acquisition was supported by the Laboratory Directed Research and Development (LDRD) funding from Argonne National Laboratory, provided by the Director, Office of Science, of the U.S. Department of Energy under Contract No. DE-AC02-06CH11357. This research uses the resources of the Advanced Photon Source, a U.S. DOE Office of Science User Facility operated for the DOE Office of Science by Argonne National Laboratory under contract No. DE-AC02-06CH11357. The authors gratefully acknowledge numerous valuable discussions with Drs. Anthony Rollett, Robert Suter and Matthew Wilkin (Carnegie Mellon University), Nicholas Porter and Dr. Richard Sandberg (Brigham Young University), Dr. Ross Harder (Argonne National Laboratory) and Dr. Anastasios Pateras (DESY). The authors gratefully acknowledge numerous valuable discussions and experimental guidance from Dr. David A. Czaplewski, Suzanne Miller, and Dr. Ralu Divan of the Center of Nanoscale Materials.

## AUTHOR CONTRIBUTIONS

The MR-BCDI algorithm was conceived and devised by S.M. and S.O.H. The Fourier transform based resampling method and other computational concepts were devised by S.M. in collaboration with M.A., S.K., Y.S.G.N. and S.O.H. The creation of the synthetic nanocrystals and simulated diffraction patterns was done by S.M. in collaboration with I.P. The fabrication of the SiC nanocrystal samples was done by N.D., S.E.S., A.D. and F.J.H. The preliminary characterization of the SiC samples with the laboratory diffractometer was done by T.D.F. in collaboration with H.Y. The measurement of the SiC coherent diffraction data was performed by T.D.F., K.J.H., W.C., Y.C. and S.O.H. The *mrbcdi* Python module was written and is currently maintained by S.M. The design of the optimization plans, MR-BCDI and phase retrieval reconstructions and post-reconstruction analyses were done by S.M. The manuscript writing, figure and table creation was done by S.M. and S.O.H. in collaboration with T.D.F. and N.D. All authors contributed to the refinement of the manuscript.

## COMPETING INTERESTS

The authors declare no competing interests.

## ADDITIONAL INFORMATION

**Supplementary information** The online version contains supplementary material available at https://doi.org/10.1038/s41524-023-01022-7.

**Correspondence** and requests for materials should be addressed to Stephan O. Hruszkewycz.

**Reprints and permission information** is available at http://www.nature.com/reprints

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_13.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© UChicago Argonne, LLC, Operator of Argonne National Laboratory 2023