## Abstract

Disorder and many body interactions are known to impact transport and thermalization in competing ways, with the dominance of one or the other giving rise to fundamentally different dynamical phases. Here we investigate the spin diffusion dynamics of ^{13}C in diamond, which we dynamically polarize at room temperature via optical spin pumping of engineered color centers. We focus on low-abundance, strongly hyperfine-coupled nuclei, whose role in the polarization transport we expose through the integrated impact of variable radio-frequency excitation on the observable bulk ^{13}C magnetic resonance signal. Unexpectedly, we find good thermal contact throughout the nuclear spin bath, virtually independent of the hyperfine coupling strength, which we attribute to effective carbon-carbon interactions mediated by the electronic spin ensemble. In particular, observations across the full range of hyperfine couplings indicate the nuclear spin diffusion constant takes values up to two orders of magnitude greater than that expected from homo-nuclear spin couplings.

## INTRODUCTION

Although the quest to understand the roles of disorder and couplings in the out-of-equilibrium dynamics of many-body systems goes back several decades (*1*), the field is presently witnessing a resurgence, in part, because of its intrinsic connection to the development of novel quantum technologies. Progress has been made largely possible by captivating experiments in cold gases where the coupling to outer reservoirs can be virtually suppressed and the evolution of each of the atoms in the interacting ensemble is probed individually (*2*–*4*). An example of recently observed phenomena is many-body localization (MBL) (*5*, *6*), a process where, despite the interactions between its inner units, the system fails to thermalize, i.e., its long-term properties cannot be captured by conventional equilibrium statistical mechanics (*7*, *8*). Unlike Anderson localization (*9*), interparticle couplings lead to dephasing of individual, initially localized states (*7*). However, the absence of exchange between different MBL modes endows these systems with a long-term memory, which makes them potentially useful platforms to store and retrieve quantum information.

Interacting spins in diamond provide an intriguing platform to investigate the interplay between localization and thermalization because electrons and nuclei feature species-specific interactions and concentrations that can be tuned and dynamically controlled. Hyperfine couplings with paramagnetic centers can take extreme values (exceeding hundreds of megahertz for first shell carbons), while the low gyromagnetic ratio and natural abundance of ^{13}C spins make homonuclear couplings orders of magnitude weaker (~100 Hz). Given our understanding of thermalization as a spin diffusion process, the large frequency mismatch between hyperfine-coupled and bulk nuclei immediately raises questions on the system’s ability to reach equilibrium. This problem, paramount to interpreting nuclear spin-lattice relaxation (*10*) but equally relevant to carrier transport (*9*, *11*), has been traditionally explained through the notion of a “spin diffusion barrier,” i.e., a virtual line in the space around a paramagnetic center separating “a frozen core” of nuclei unable to communicate (i.e., “flip-flop”) with bulk spins (*12*–*14*) .

Here, we combine optical excitation and nuclear magnetic resonance (NMR) at low magnetic fields to investigate the generation and transport of nuclear magnetization in a diamond crystal hosting nitrogen-vacancy (NV) centers. Formed by a substitutional nitrogen immediately adjacent to a vacancy, these spin-1 point defects polarize efficiently under green illumination, which can be exploited to dynamically polarize the ^{13}C nuclei in the crystal. Working under “energy matching” conditions, where NVs cross-relax with surrounding spin-1/2 nitrogen impurities or “P1 centers,” we find that strongly hyperfine-coupled carbons can efficiently exchange polarization with bulk nuclei; this process is made possible by many-body interactions involving electron and nuclear spins through mechanisms that we formally capture via a nuclear spin–only effective Hamiltonian. Furthermore, we measure nuclear spin diffusion constants across a range of hyperfine couplings orders of magnitude greater than the nuclear Larmor frequency and find values ~100-fold bigger than those possible via homonuclear couplings, a phenomenon that we interpret in terms of electron-mediated interactions between distant carbons.

## RESULTS

^{13}C hyperfine spectroscopy at low magnetic fields

Figure 1 (A to C) summarizes the conditions in our experiments. We study a diamond sample with a large NV and P1 content [~10 and ~50 parts per million (ppm), respectively] produced via high-energy electron irradiation and annealing. We operate in the regime of “cross-relaxation” where the separation between the ∣*m*_{S} = 0⟩ and ∣*m*_{S} = −1⟩ energy levels of the NV spin approximately matches the P1 Zeeman splitting in an external magnetic field *B*^{(0)}, whose exact value depends on the angle θ with the NV axis (*15*, *16*). Optical pumping of the NV induces dynamic nuclear polarization (DNP) of bulk ^{13}C, which we subsequently detect using a field cycling protocol (Fig. 1C). Figure 1D shows the amplitude of the observed ^{13}C NMR signal as a function of the optical pumping magnetic field *B*. The DNP generation can be simplistically understood through an energy-conserving NV-P1-^{13}C process where nuclear spins polarize positively or negatively depending on the sign of the difference between the NV and P1 transitions above or below *B*^{(0)}. On the other hand, the fact that sizable DNP can be observed for a field mismatch as large as ~0.4 mT (corresponding to hyperfine couplings of order ~10 MHz) immediately points to nontrivial channels of polarization transfer from nuclear spins strongly coupled to defects.

To measure the ^{13}C spectrum at a given optical pumping magnetic field, we apply a radio frequency (RF) pulse immediately after laser illumination (before sample shuttling; Fig. 1E) within a range around the ^{13}C Zeeman frequency. The pulse duration (1 s) is chosen so as to make the up/down ^{13}C spin populations equal when on resonance, hence leading to a “dip” in the observed signal amplitudes plotted as a function of the RF frequency. Figure 1F shows the results for variable RF power. In the limit of weak RF excitation (−25 dBm), the NMR linewidth amounts to ~1.5 kHz, coincident with that observed at high field (see section S1). Stronger RF power results in broader dips, a consequence of the greater excitation bandwidth; in the experiments below, we use an RF power of −8 dBm, which confines the effect to a ~13-kHz band around resonance.

The ability to manipulate ^{13}C spins gives us the opportunity to probe the transport of spin magnetization from paramagnetic centers to “bulk” (i.e., very weakly coupled) carbons as it cascades down across nuclear spins with different hyperfine couplings under NV-P1 cross-relaxation. Given the multispin nature of the dynamics at work, this process is better visualized in frequency space as a sequence of jumps along a chain formed by groups of carbons with varying hyperfine coupling (and, hence, different resonance frequencies; Fig. 2A). Nuclear spins proximal to paramagnetic centers (NVs or P1s) are normally invisible in the standard NMR signal because of their comparatively low abundance and extreme hyperfine-induced gradients. Nevertheless, their ability to mediate the transfer of magnetization to bulk spins can be selectively exposed through the accumulated effect of RF excitation on the polarization buildup during optical spin pumping.

Initial evidence revealing the nontrivial role of hyperfine-coupled carbons is shown in Fig. 2C, where we compare the NMR signal amplitudes following simultaneous RF and laser excitation (Fig. 2B). Accompanying the expected dip near the ^{13}C Zeeman transition (ω_{I}~2π × 560 kHz), we observe (partial) NMR signal reduction over a wide frequency range (green squares), far exceeding the excitation bandwidth (faint, 13-kHz-broad Gaussian in the back here serving as a reference). We find that this effect persists at even higher frequencies, where intercarbon flip-flops should be strongly suppressed. This is further shown in Fig. 2D, where we measure the equivalent of a hyperfine-resolved spectrum over a 160-MHz range, selectively sensitive to nuclear spins participating in the magnetization transport. We identify several high-frequency regions where RF excitation has a substantial impact on the observed NMR signal, suggesting that localization, the regime naively anticipated for hyperfine-coupled carbons in a dilute nuclear spin system such as diamond, cannot capture the dynamics at play. Very much on the contrary, we show next that most nuclear spins communicate efficiently with each other despite their relatively large frequency mismatch.

To shed light on the underlying mechanisms, we start with a comparison between the RF absorption spectrum in Fig. 2D and the set of hyperfine couplings to NVs and P1s (respectively, colored bands in the background and vertical bars in Fig. 3A). We find a moderate correlation between the two. For example, the dip at ~40 MHz, associated with a second shell carbon around the P1 center (*17*), suggests that substitutional nitrogen plays an important role in enabling spin exchange between near-defect and bulk nuclei. The dip disappears if one shifts the magnetic field from *B*^{( + )} to *B*^{( − )}, a change of only ~0.2 mT (see Fig. 2E), suggesting that spin diffusion emerges from a multispin process requiring precise alignment between the NV, P1, and ^{13}C energy levels. This notion is consistent with the very premise of DNP near ~51 mT, arising from nuclear spin–assisted NV-P1 cross-relaxation at these fields (*15*, *18*). On the other hand, one cannot rule out spin-lattice relaxation effects, as the bulk carbon *T*_{1} time is also seen to moderately change, from ~5 to ~7 s when transitioning from *B*^{( + )} to *B*^{( − )}. Last, the ~97.5-MHz resonance, which we could not match to any reported ^{13}C site near the NV or P1, may instead correspond to polarization pathways involving the nuclear spin of the ^{14}N host at the P1 [known to participate in the polarization transfer (*15*, *16*, *18*)]; additional work, however, will be needed to clarify its origin.

The absence of RF absorption is also an important indicator. For example, the flat response in Fig. 2D near ν_{RF} ~ 130 MHz, coincident with the hyperfine splitting of first shell carbons around the NV (*19*), indicates that these sites do not partake in the polarization transfer process, hence suggesting that select nuclear spins, featuring exceedingly strong hyperfine interactions, fail to thermalize with the rest (see below). By the same token, no RF dips are observable between ~50 and ~90 MHz (omitted in Fig. 3A for simplicity), a range with no hyperfine-coupled carbons (*17*, *20*, *21*).

More generally, the amplitude of the RF absorption dip reflects on the number of diffusion channels available to the system near a given excitation frequency ν_{RF} (Fig. 3B). A complete transport blockade, manifesting in the form of a full-contrast dip, is possible only when the nuclear spins resonant with the applied RF intervene in every polarization transfer event. As the number of alternative channels increases, the RF-induced contrast diminishes because most spin diffusion pathways do not involve resonant nuclei. The latter, of course, depends on the granularity of the frequency jumps δν_{d}(ν_{RF}), characterizing the multispin configurational change during spin diffusion; greater RF absorption can be regained as δν_{d}(ν_{RF}) becomes comparable to (or smaller than) the RF bandwidth δν_{b} (~13 kHz in the present experiments). We believe that this interplay is responsible for the DNP signal response below ~10 MHz, where the number of carbon sites with comparable hyperfine couplings (and, thus, the number of spin diffusion pathways) increases rapidly, while the nuclear spin energy difference δν_{d}(ν_{RF}) in each jump gradually fades away. On a related note, a close inspection of Fig. 2D shows a slight offset relative to the signal amplitude observed in the absence of RF (faint horizontal lines). We presently ignore its origin but hypothesize that it could stem from weak RF absorption between many-body electron spin states (i.e., weakly allowed “zero-quantum” transitions), which subsequently causes nuclear spin relaxation. Additional experiments, however, are mandatory to clarify this point.

### Effective Hamiltonian and spin diffusion dynamics

Deriving a Hamiltonian that correctly reproduces the behavior of interacting electron and nuclear spin ensembles from first principles, a problem at the center of ongoing efforts (*22*), remains a challenging task. Here, we capture the dynamics at play by considering a pair of carbons, each interacting with one of two P1 centers, which, in turn, couple dipolarly to each other (Fig. 3C). Focusing first on the “hyperfine-dominated” regime (where ‖*A*_{1}‖ ~ ‖*A*_{2}‖ > ℐ_{d} > ω_{I}), we find that the polarization can flow from one carbon to the other with an effective rate _{d} is the interelectronic dipolar coupling constant and *J*_{eff} can reach sizable values when the electron spin concentration is sufficiently high. As an illustration, for an electron spin dipolar coupling 〈ℐ_{d}/2π〉~1 MHz [corresponding to a nitrogen concentration of ~10 ppm (*23*)], we obtain *J*_{eff}/2π ~ 1 kHz for

While the above effective coupling allows most hyperfine-shifted nuclei to communicate, we also find that transport can be suppressed if the hyperfine shift difference δ*A* = ‖*A*_{1}‖ − ‖*A*_{2}‖ between the two carbons is large. More formally, we express the condition for delocalization as*A* (i.e., when ‖*A*_{1}‖ ≫ ‖*A*_{2}‖; see section S3). This is likely the scenario for first shell carbons (*A*~130 MHz), separated from the rest by a large spectral gap (see Fig. 3A). For completeness, it is worth mentioning that in the “dipolar-dominated” regime (where ℐ_{d} > *A*, ω_{I}), the effective nuclear spin coupling takes the form *J*_{eff} ~ *A*^{2}/(4ℐ_{d}). This expression shows, as expected, vanishing interaction for nuclei decoupled from paramagnetic defects (ℐ_{d} > ω_{I} > *A*), but it also suggests that *J*_{eff} can be quite strong, potentially exceeding 10 kHz in the narrow window where ℐ_{d} > *A* > ω_{I} (sections S2 and S3).

ℐ_{d}-induced state mixing activates transitions at frequencies other than those expected for pure nuclear spin flips. This is shown in Fig. 3 (C and D), where we plot the calculated nuclear spin polarization in a ^{13}C-P1-P1-^{13}C chain under continuous RF excitation, assuming that both carbons start from a polarized state (see also section S4). When ℐ_{d} ~ 0, the system absorbs selectively at the single nuclear spin hyperfine transitions. As ℐ_{d} increases, however, new dips corresponding to simultaneous nuclear and electron spin flips emerge. Given the range of possible spatial configurations in disordered spin ensembles, RF excitation should therefore yield broad bands of less-than-optimal DNP crudely centered around the hyperfine transitions, in qualitative agreement with our observations.

From the above considerations, we surmise that the ensemble of paramagnetic defects can be thought of as an underlying network providing the couplings required for nuclear spins to thermalize (*24*); correspondingly, the spin Hamiltonian for a group of (otherwise noninteracting) *N*_{I} carbon spins **I*** _{i}* takes the form (see section S3)

*i*th nuclear spin site and

*J*

_{eff, zz}and

*J*

_{eff, xy}represent effective electron spin–mediated internuclear couplings, which, in general, must be seen as functions of the applied magnetic field and electron spin concentration.

While chain-like systems are often integrable, added spatial dimensions break any underlying symmetry and typically render the dynamics chaotic. A realistic simulation of the system at hand requires, therefore, the use of multidimensional spin arrays, an increasingly challenging task as the number of nuclei grows. Here, we qualitatively test the dynamics of the Hamiltonian in Eq. 2 using a model nuclear spin set of 22 carbons in a Cayley tree geometry, assuming that only the central spin is initially polarized (Fig. 3E). The effective couplings between nuclear spins in different rings grow to the outside of the tree, so as to emulate the transition between the hyperfine-dominated (i.e., *A* > ℐ_{d} > ω_{I}) and the dipolar-dominated (i.e., ℐ_{d} > *A* > ω_{I}) regimes. It is known that the interplay between the terms linear and bilinear in *I ^{z}* [respectively corresponding to the local potential and interaction terms of the Hubbard Hamiltonian in a carrier transport picture (

*7*) ] may lead to MBL. To make the numerical problem tractable, we assume below that the flip-flop terms are dominant, and thus, the system is in an ergodic phase.

To compute the many-body spin dynamics, we use a Trotter-Suzuki (TS) decomposition assisted by quantum parallelism (*25*). Unlike other, more common approaches (*26*, *27*), this technique does not require truncation of the Hilbert space and is thus applicable to long times (see Materials and Methods). As shown in Fig. 3F, we observe a diffusive (i.e., recurrence-free) evolution, pointing to the onset of quantum chaos (*28*–*30*). Chaoticity arises in the Cayley geometry as a consequence of the system branching, effectively enlarging the size of the accessible Hilbert space as the polarization moves from inner to outer rings. Note that despite the growing internuclear couplings, the characteristic time constant (of order ~2 ms) is uniform across the tree structure, dictated by the higher-order (and hence weaker) effective electronic couplings communicating the central spin with nuclei in the first ring.

Experimentally, we probe the time scale of spin diffusion in our sample via the protocol in Fig. 4A where we evenly distribute RF pulses of fixed duration throughout the illumination interval; the pulse length is chosen so as to ensure several ^{13}C Rabi cycles (section S1). The upper half of Fig. 4B shows an example plot corresponding to RF at 10 MHz. For interpulse intervals τ ≳ 10 ms, we find that the effect of RF pulses on the hyperpolarization amplitude is negligible, an early indication that spin diffusion takes place on a time scale faster than that deriving from direct internuclear dipolar couplings (averaging ~100 Hz in nonenriched diamond). Overall, our data can be reasonably described via a stretched exponential dependence of the form *S* = *S*_{0} − *S*_{1} exp(− (τ/τ_{d})^{ε}), where τ is the interpulse separation and all other variables are fitting parameters, with τ_{d} representing the characteristic nuclear spin diffusion time. We find that ε < 1, typically indicative of heterogeneity in the underlying physical process (*31*, *32*). This idea is consistent with the multichannel nature of the transport dynamics at play, here expressed via the probability distribution ℒ(μ, ε), satisfying _{d} (lower half plot in Fig. 4B). In particular, we identify a broad set of fast transport processes whose rates extend beyond ~1 ms^{−1} (shadowed tail in the plot).

To capture these observations into a functional microscopic model, we now return to the notion of magnetization transport along a one-dimensional (spectral) chain formed by *m* spin sets {*N _{i}*},

*i*: 1…

*m*, each featuring resonance frequencies within bands Δν

*centered around effective hyperfine couplings ‖*

_{i}*A*‖ (Fig. 4C). Aiming at a qualitative comparison with experiment, this time, we model the transport problem classically using a set of coupled differential equations adapted to describe magnetization hops between boxes in the presence of optical and RF excitation, as well as nuclear spin-lattice relaxation (section S5). A train of RF pulses (resonant with nuclear spins within a band δν

_{i}_{b}along the chain) partly disrupts the transport of polarization and leads to a change δ

*M*in the magnetization stored in box

_{m}*m*(here, serving as the observable). This effect saturates in the limit where the interpulse separation τ is equal to (or shorter than) the interval required to replenish the magnetization in the depleted cell (roughly, the inflexion point in the sigmoid), hence allowing us to extract the characteristic spin diffusion time τ

_{d}at the corresponding excitation frequency ν

_{RF}(Fig. 4D).

Figure 4D summarizes numerical results from a chain of *m* = 40 spin cells, each connected to its immediate neighbors via transfer rates γ_{i, i + 1}. To establish a starting connection between τ_{d} and the underlying rates, we first investigate the case where γ_{i, i + 1} takes a constant value γ throughout the chain (Fig. 4D). As expected, we find that τ_{d} grows inversely with γ^{−1}, although the dependence is not linear, a consequence of the finite duration of the RF pulse. To investigate the impact of transport heterogeneity, we also consider the case where γ_{i, i + 1} takes on different values depending on the position across the chain, peaking at the midpoint. Imposing greater transfer rates between cells effectively amounts to fusing neighboring groups of spins into a larger cell, hence amplifying the impact of individual RF pulses resonant with the set and thus altering τ_{d}. Since the experimental response upon excitation at different frequencies does not substantially depart from our observations in Fig. 4B (see section S5), we tentatively conclude that the transfer rates across the hyperfine spectrum, or, more generally, the representative values from the transfer rate distributions connecting each cell in the chain with all others, are relatively uniform. On a related note, our numerical model exhibits only a moderate departure from a single exponential response (ε~0.9 in Fig. 4D). The latter could well be a consequence of the first neighbor coupling structure assumed for the spin chain, likely oversimplifying the system complexity by limiting the number of channels available to the transport of nuclear polarization.

Capitalizing on the above considerations to interpret our observations, we find that the characteristic diffusion rate ^{−1}. The agreement with the quantum model in Fig. 3 (E and F) should be considered rather fortuitous, as a numerical value of the diffusion time can only emerge from a suitable average over the set of possible spin configurations. On the other hand, given the mean intercarbon distance in diamond 〈*r*_{C}〉 = 0.5 nm, we conclude that the effective diffusion constant observed herein can be as large as ^{2} s^{−1}, about 100-fold greater than that derived from nuclear dipolar interactions alone (Fig. 4C) (*33*, *34*). This result reinforces the understanding of the cross-relaxing electron spin bath as a mediator to swiftly move around magnetization from otherwise many-body localized groups of nuclei. This behavior could prove advantageous to expedite the transport of polarization across the diamond surface into arbitrary nuclear spin targets (*16*, *35*).

## DISCUSSION

While the effective Hamiltonian in Eq. 1 supports the notion of a coherent, electron-mediated nuclear spin transport, a question of interest is whether spin-lattice relaxation (electronic or nuclear) affects the diffusion process itself (beyond imposing a limit on the polarization buildup). Supporting this notion, recent numerical studies suggest that incoherent dynamics can help drive the spin system away from “blockade” regimes, i.e., spin configurations that prevent the transport of spin polarization (*35*). In the present framework, these processes could, e.g., flip P1 centers that have previously been polarized upon cross-relaxation with the NVs. Note that coherent channels remain the main transport driver and interacting paramagnetic defects are still central to the process, but here, it is spin-lattice relaxation (not necessarily electron spin diffusion) that prepares the P1 for the next cycle of spin transport. This picture is consistent with the measured τ_{d}, comparable, on average, to the NV or P1 spin-lattice relaxation times (of order 1 ms in this diamond sample). Future experiments, for example, above and below room temperature or for samples with variable NV or P1 concentrations could help shed light on the role of incoherent processes.

Extensions of the ideas introduced herein can provide additional insights on the complex spin dynamics at play. For example, the use of chirped micro-wave (MW) pulses to induce nuclear spin polarization (*36*, *37*), away from the NV-P1 cross-relaxation condition, can be exploited to separate the roles of NVs and P1s during the spin diffusion process. Along the same lines, microwave manipulation of the electron spin bath should give us the opportunity to controllably reintroduce localization in the nuclear spin system or to count the number of correlated carbons as the polarization spreads (*38*, *39*). Particularly attractive is the combined use of super-resolution microscopy (*40*–*42*) and magnetic resonance techniques to monitor the spin dynamics of small ensembles of nuclear spins communicating via NV-P1 networks, which could be relevant to quantum information processing with many-body disordered systems (*43*).

While our experiments centered on spins in diamond, we anticipate that similar techniques can be adapted to investigate the dynamics of other material systems hosting spin active nuclear and electronic spins. These include organic systems exhibiting (nonoptical) DNP, where simultaneous microwave and RF excitation could be exploited to gain information on nuclear spins proximal to radicals, normally invisible in standard DNP-enhanced NMR experiments.

## MATERIALS AND METHODS

### Experiment

Throughout our measurements, we use a chemical vapor deposition–grown, type 1b diamond, which was previously electron-irradiated and annealed to create NV centers throughout the bulk crystal at an approximate concentration of 10 ppm (*15*). DNP is carried out via 1-W laser excitation at 532 nm. We rely on a pair of coils and the stray field of a 9-T NMR magnet to adjust the magnetic field during optical illumination in the vicinity of 51 mT and use a home-made, compressed air–driven device to shuttle the NMR probe between the polarization (51 mT) and detection fields (9 T). We orient the diamond crystal so that the external magnetic field nearly coincides with one of the NV axes and use a single loop around the crystal as the source of RF excitation. All experiments are conducted at room temperature. Additional details, including a characterization of the ^{13}C spin response as a function of the RF power, are presented in section S1.

### Numerical simulations

We use a fourth-order TS method (*44*) to numerically evaluate the time dependence of large spin systems. The TS protocol avoids manipulating and diagonalizing the full Hamiltonian *H*, instead approximating the total evolution operator *U*(*t*) = exp { − *iHt*} by a suitable sequence of partial evolution operators*H _{k}*} corresponds to each of the single-spin and two-spin terms in the Hamiltonian

*H*, proportional to either linear (

*t*requires the successive application of the step-like evolutions. The approximated dynamics remains always unitary, and the accuracy of the approximation relies on the TS time step δ

*t*being sufficiently small as compared to the shortest local time scale of the original Hamiltonian. In our simulations, we have tuned the TS time step so that, for the system size considered (22 spins), relative error bounds are estimated to be 10

^{−4}. Regarding the computational implementation, an exponential speedup of our simulations is achieved by means of a massive parallelization scheme via general-purpose graphical processing units (

*45*).

## SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/6/18/eaaz6986/DC1

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is **not** for commercial advantage and provided the original work is properly cited.

## REFERENCES AND NOTES

**Acknowledgments:**

**Funding:**D.P., J.H., and C.A.M. acknowledge support from the NSF through grants NSF-1903839 and NSF-1619896 and from the Research Corporation for Science Advancement through a FRED award. They also acknowledge access to the facilities and research infrastructure of the NSF CREST IDEALS (grant number NSF-HRD-1547830). J.H. acknowledges support from CREST-PRF NSF-HRD 1827037. All authors acknowledge the CUNY High Performance Computing Center (HPCC). The CUNY HPCC is operated by the College of Staten Island and funded, in part, by grants from the City of New York, State of New York, CUNY Research Foundation, and NSF grants CNS-0958379, CNS-0855217, and ACI 1126113.

**Author contributions:**C.A.M., D.P., and P.R.Z. conceived the experiment. D.P. and J.H. conducted the experiments. P.R.Z. developed a model and carried out theoretical calculations. A.A., R.H.A., J.A.R., and A.P. advised on several aspects of theory and experiments. C.A.M. wrote the manuscript with input from all authors. R.H.A., J.A.R., A.P., and C.A.M. supervised the overall research effort.

**Competing interests:**All authors declare that they have no competing interests.

**Data and materials availability:**All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors. All correspondence and request for materials should be addressed to C.A.M. (cmeriles{at}ccny.cuny.edu).

- Copyright © 2020 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution NonCommercial License 4.0 (CC BY-NC).