Probing molecular dynamics at the nanoscale via an individual paramagnetic centre

We demonstrate a protocol using individual nitrogen-vacancy centres in diamond to observe the time evolution of proton spins from organic molecules located a few nanometres from the diamond surface. The protocol records temporal correlations among the interacting protons, and thus is sensitive to the local dynamics via its impact on the nuclear spin relaxation and interaction with the nitrogen vacancy. We gather information on the nanoscale rotational and translational diffusion dynamics by analysing the time dependence of the nuclear magnetic resonance signal. Applying this technique to liquid and solid samples, we find evidence that liquid samples form a semi-solid layer of 1.5-nm thickness on the surface of diamond, where translational diffusion is suppressed while rotational diffusion remains present. Extensions of the present technique could be exploited to highlight the chemical composition of molecules tethered to the diamond surface or to investigate thermally or chemically activated dynamical processes such as molecular folding.


Supplementary Figure 7. Simulated correlation envelopes for the solid--state samples
The Merckoglas data can be readily described by a completely static bath of molecules, i.e. by a purely dipolar NMR relaxation, yielding a single exponential decay envelope (yellow curve), in good agreement with Eq. 1 in the main manuscript. (b) The PDMS data can be readily reproduced by either assuming slow molecular rotation with no center--of--mass translation (blue curve, D R,surface = 0.001rad 2 /µs), or by assuming a completely static layer with a (~20%) reduced proton density of 40 nm --3 (green curve).  [nm]. (d) Influence of anomalous diffusion on the Initial signal decay. In (a) through (d) the data points correspond to NV 6 and the diamond coating is sample C (immersion oil). In all cases only the noted parameter is changed while the others have a fixed value coincident with that noted in the main text.

Signal reproducibility
All experiments were carried out with the same diamond sample, and thus with the same set of shallow NVs. However, because the diamond crystal must be physically removed from the microscope for surface coating, it is difficult to use the same individual NV for testing nuclear spins from different sample films. To circumvent this complication, we collected data from multiple (50 to 100) individual NVs exposed to the same sample film (see below). Among these only a few can be dynamically decoupled for a time sufficiently long to see the nuclear spin signature via an XY8--N sequence (a likely result of the NV depth dispersion and heterogeneity of the local concentration of paramagnetic impurities). We find that not every NV center exhibiting a nuclear--spin--induced dip under the XY8--N sequence also shows detectable nuclear spin correlation. The reason for this is still not fully understood but it may relate to the fortuitous overlap between the proton spin signature and the fourth harmonic of the 13 C--induced dip (see below): The latter can be easily mistaken by the former in an XY8--N sequence but not in the correlation protocol (where the nuclear spin Larmor frequency is directly probed). Also worth noting is the small relative amplitude of the correlation signal (≈ 10% of the maximum possible fluorescence contrast between the NV spin states), which leads to a correspondingly low SNR and makes detection difficult. Within these experimental limitations, Figs. 2 and 3 in the main text show representative data sets from the NVs that did display a sizeable correlation signal. Among the latter, excellent reproducibility was observed between data sets in all the samples we explored, as shown immediately below.

Sample A -Merckoglas
We examined a set of 52 NV centers, out of which 10 could be dynamically decoupled to coherence times of ≈ 100 µμs, allowing us to detect the proton signal via an XY8--10 sequence. Only 4 of these NVs show a detectable signal after an XY8--3 correlation sequence.
Supplementary Figure 1 shows three of those NV centers, the one presented in the main text (NV20) and two additional ones. We find a reproducible correlation signal, described by an exponentially damped sinusoid centered around the proton Larmor frequency. The decay envelope (here maintained unchanged for all NVs in the figure) is qualitatively captured by the presented model, with a decay time of 20 µs.

Sample B -PDMS
For the PDMS coating we examined a set of 85 NV centers, out of which 13 showed proton signal via an XY8--10 sequence. Only 2 of these 13 NVs reveal a noticeable signal after an XY8--3 correlation sequence.
The two correlation curves shown in Supplementary Figure 2 show the curve presented in Supplementary Figure 3 of the main text (NV1) and a second curve for a second NV (NV78). As explained above, experiments with this sample suffered from poor collection efficiency and required particularly long integration times, the reason why the second data set was measured only for evolution times < 20 µμs. Both measured signals are overlaid with the same calculated decay envelope, assuming a slow molecular dynamics of the PDMS molecules as described in the main text (see also Supplementary Note 4 below).

Sample C -Immersion oil
For the immersion oil a set of 83 NV centers was measured. Among them, nine NVs exhibited proton signal after an XY8--10 sequence. Five of those NVs showed an XY8--3 correlation signal. We measured two of those NVs over a longer evolution time window (> 80µμs, NV6 is the one shown in the main text and related SI notes). As shown in Supplementary Figure 3, we measure virtually identical responses characterized by a fast initial decay and a long--lived tail as highlighted in the main text.
Remarkably, the two data sets are so consistent that the same set of parameters (adsorbed layer thickness, translational diffusion constant, and rotational diffusion constant) can be used to describe both observations.

Origin of the signal
Recently it was pointed out that the XY8--N noise spectroscopy is prone to measure spurious effects of higher Larmor frequency harmonics 1 . Especially the 4 th harmonic of the 13 C nuclei intrinsic in the diamond lattice is spectrally located very close to the proton Larmor frequency, and can be easily misinterpreted as the signal of external protons. Our protocol is immune to this problem because it directly records the proton spin precession during , thus leading to resonance spectra at the specific nuclear Larmor frequency 2 . In particular, any residual effect from carbon spins during the XY8--N segments of the correlation protocol would lead to a peak at the carbon spin resonance frequency, which can be easily separated from proton--induced contributions.

Modeling the correlation signal envelope
The envelope of the correlation signal is described by (see Eq. (2) of the main text) where !,! is the conditional probability that the i--th nuclear spin of the j--th molecule remains in the detection volume over the evolution time , and ! !,! is the transverse relaxation time of said nuclear spin. Here we neglect any additional contribution to the envelope arising from longitudinal relaxation of the NV electron spin (see below Supplementary Note 6). This is a reasonable assumption as typical longitudinal relaxation times are in the order of hundreds of µs, while we focus here on effects on timescales of a few tens of µs, i.e. the time range of our measurements.
Before describing how we model the probability and relaxation time in detail, we introduce below the underlying assumptions and simplifications of the presented theory.
First, we assume that the value of the sum in Supplementary Eq. (1) is time independent, i.e. the total number of protons/molecules inside the detection volume stays approximately constant over time. Therefore the amplitude and decay components of the correlation signal depend solely on the spatial distribution of the molecular mobility with respect to the surface distance. Further, we neglect any hydrodynamic interactions among molecules other than those implicit when assigning a candidate translation or rotation diffusion coefficient.
The molecules are treated as rigid spheres of radius , and their molecular mobility varies on the length scale of the detection volume. This heterogeneity is described by position dependent diffusion constants ! ( ) and ! ( ) for the translational and rotational diffusion respectively, which takes into account the possible presence of adsorbate molecules on the diamond surface. We estimate that the transition between adsorbate and free molecules is sharp and that the latter obey classical diffusion equations, as well as the Stokes--Einstein--(Debye) equations for the diffusion where ! is the Boltzmann constant, is the medium temperature, and is the viscosity of the sample medium.
We further assume that the diffusion heterogeneity is restricted to the z--axis perpendicular to the where !"#$ , and !"#$%&' are the respective values of the diffusion constants in the bulk of the volume and at the surface, and is a transition parameter, which we keep fixed at an empirical value of = 0.25 nm. We note that this value is in the same order of magnitude as that found for AFM experiments measuring the transition of surface to bulk fluid dynamics 3 . We find that the exact value of has a negligible impact on the results of the simulation, as long as the transition takes place on a scale small compared to the radial size of the detection volume.
The coherences of the nuclear spin bath characterized by the transverse relaxation time ! !,! can be calculated using established means known from classical NMR 4 . In the following we will briefly outline the calculation steps of the relaxation time, before describing the model we use to estimate the conditional probability !,! .

Nuclear spin relaxation time
The nuclear spin dynamics mainly depends on position dependent correlation times of motion ! .
The only interaction considered among the nuclei is homonuclear dipolar coupling, so that one can distinguish different time regimes by the ratio of the rate of nuclear spin reorientation (∝ ! !! ) versus the dipolar coupling constant ! .
Among the interactions responsible for the correlation decay one can distinguish between intramolecular interactions inside a molecule, and intermolecular interactions among spins of different molecules. In the interior of a molecule, the variation of the dipolar coupling between spins arises almost exclusively from the rotation of the molecule (neglecting distance variations due to vibrations), while for the interactions between spins in different molecules their relative translation must be considered.

1H--1H interaction
The homonuclear dipolar coupling between two proton spins is described by the Hamiltonian where ! is the magnetic constant, ℏ is the reduced Planck constant, !" is the gyromagnetic ratio of a proton spin, ! are the spin vector operators of the i--th proton, and !! is the distance vector between two proton spins. Squaring Supplementary Eq. (5), averaging over the azimuthal and polar angles, and taking the trace over the spin matrices yields an expression for the rms--dipolar coupling where we use a representation of !! in spherical coordinates. In Supplementary Eq. (6) !! is the average distance to the nearest nuclear neighbor, which is dependent on the proton density and is distributed according to This results in an average proton--proton distance of The relaxation rate of the axial magnetization component for a 2--spin cluster due to their mutual dipolar interaction is calculated via 4 Since we assume a constant proton density throughout the detection volume, the average interaction strength and dipolar relaxation rate are therefore position independent. Using a typical value for the proton density in organic compounds of = 50 nm !! for all samples, yields a fluctuation rate of !"# ≈ 20 kHz, and average interaction strength of ! ! ≈ 45 kHz.

Translational diffusion of 1H spins
The translational diffusion of molecules generates a small fluctuating magnetic field. The fluctuation rate due to translational diffusion is given by

Rotational diffusion of 1H spins
For a spherical molecule of radius rotating in a liquid of viscosity the rotational fluctuation rate is described by Stokes' law Using the Stokes--Einstein relation (Supplementary Eq. ( 3)) we find

Combined dynamics and relaxation rate
In order to describe the nuclear spin relaxation dynamics of the sample molecules we use an autocorrelation function of the fluctuating magnetic field B(t) experienced by the nuclear spins where ! is the rms field strength experienced by the protons spins, and ! is the characteristic correlation time of the magnetic fluctuations, which includes all above mentioned mechanisms, i.e.
. The spectral density ( , ! ) is twice the Fourier transform of Supplementary Eq. (14) , The normalized spectral density ( , ! ) is The transverse relaxation is caused by the dephasing of the proton spins due to their mutual dipolar interactions. In the case of intramolecular (homonuclear) dipole--dipole relaxation, the transverse relaxation rate ! !! between two spin ½ nuclei is given by 4 where ! is the nuclear Larmor precession frequency, and ! ! is the average dipole--dipole coupling constant, given by Supplementary Eq. (6).

Probability ,
We model the probability !,! using a semiquantitative geometric model. Simply put, we draw a 'diffusion sphere' around a starting position with a radius which is the rms--diffusion distance in three dimensions for Brownian motion. We only consider translational diffusion here, and the translational diffusion constant ! only depends on the starting position (Supplementary Figure  4a,  Supplementary  Figure  5a). The conditional probability to find the molecule inside the detection volume ( !" ) after the time is then approximately given by the ratio of the overlap between the diffusion volume ( !"## ) and the detection volume, normalized to the detection volume (Supplementary Figure 5b): The detection volume is simplified as a hemisphere with radius !" = 4 nm , which roughly translates to an expected volume of (5 nm)³.
The overlap !"## ∩ !" can be calculated via the intersection of two spheres where is the absolute value of the starting position vector , and !"## is obtained from Supplementary Eq. (18). The intersection volume will have a negative value if either the diffusion sphere is completely inside the detection volume ( !,! = 1), or if the diffusion sphere becomes larger than the detection volume. In the latter case the probability !,! is determined by the ratio of the diffusion sphere volume to the detection volume.
Since the detection volume is only a hemisphere, we have to subtract any intersection or diffusion volume contributions extending below the surface before calculating the value for the probability.
The volume beneath the surface can be calculated as the volume of a sphere cap ( !"# ). Therefore the volume below the surface can be calculated as with !"# = ! ! ℎ ! (3 !"## − ℎ) , where ℎ = !"## − is the height of the cap, and is the z-component of the position vector . The volume fraction to be subtracted is approximately given by where = !"## ! − ! is the radius of the circular base of the cap volume, and = + − !" =

Monte Carlo (MC) simulation of the molecular dynamics in the detection volume and comparison to the experimental data
For the simulation we pick random starting positions inside the detection volume, then calculate the probability to remain within this volume for varying as described above of the signal, and we find best agreement for a value of !,!"#$ ~0.3 nm 2 /µs ( Supplementary Fig. 6).
We note that this is in the order of magnitude that one would expect from Supplementary Eq. (2), which yields a value of !,!"#$ ~ 1.0 nm 2 /µs, corresponding to a molecule of radius a = 0.5 nm, a medium of viscosity = 437 mPa.s (as listed in the compound datasheet), and room temperature conditions. The deviations could be readily explained by varying the molecule radius a, as the immersion oil contains various different molecules of different chain length etc., which are treated identically in the simulations.
Our model also manages to describe the decay times of the purely exponential signal envelopes for the solid polymers quite well. Here we assume no translational diffusion throughout the detection volume, i.e. the whole volume is treated as an 'adsorption layer'. We then vary the value of the rotational diffusion constant D R,surface until we obtain a reasonable match with the experimental data.
In the case of Merckoglas (sample A) we can only reproduce the signal shape by removing the rotational diffusion as well, i.e. molecules remain static over the course of the measurement. In this case the probability !,! becomes unity, and the signal decay is mainly governed by NMR relaxation due to inter--and intramolecular dipolar interactions, yielding a fast single exponential decay (Supplementary Figure 7a) as expected for a solid.
For PDMS the correlation signal also decays exponentially but the decay time (of order ~30 µs) is found to be longer than anticipated for an ensemble of static, dipolarly--coupled protons (at least assuming a standard proton density of = 50 nm --3 ). We can reproduce the observed behavior by assuming a rotational diffusion constant of D R,surface = 0.001 rad 2 /µs throughout the detection volume. The latter is supported by prior NMR experiments showing that PDMS can exhibit significant molecular dynamics at room temperature, thus leading to partial motional narrowing 7 . The self-diffusion depends on the PDMS chain length/molecular weight, i.e., on the mixing ratio of the components. Even if partly, this enhanced mobility can be the reason of the line narrowing captured in the non--negligible value of the rotational diffusion constant in the simulation. For completeness, we note that the observed data set can also be reproduced by a static ensemble of protons with a 20% lower density than normal (40 nm --3 as opposed to 50 nm --3 , lighter green trace in Supplementary Figure 7b). However, in light of the experimental evidence mentioned above, this alternative seems less realistic.

Model precision and accuracy
As described above the key parameters in our model are the bulk translational diffusion coefficient !,!"#$ (which governs the initial signal decay for a liquid sample), the rotational diffusion coefficient in the semi--mobile adsorption layer !,!"#$%&' (describing the decay time of the long lasting signal tail), and the thickness of the adsorption layer !"# (which sets the amplitude of this long signal component). Supplementary Figure 6 schematically visualizes the role of each parameter on the resulting decay envelope.
As shown in Supplementary Figure 8 for Sample C, we find that the modeled response of NV 6 is quite sensitive to the values assigned to these parameters. For example, Supplementary Figure 8a shows Over the measured proton evolution time we find that the correlation amplitude stays approximately constant, with no distinct signal decay. Therefore we can only estimate a lower bound for the rotational diffusion coefficient for which we can reproduce this behavior. We find that this signal shape is reproduced by rotational diffusion coefficients larger than !,!"#$%&' = 0.05 , and that a notable deviation from the data can be observed for coefficients below !,!"#$%&' = 0.02 . Therefore within the measurement region the model can distinguish changes below 0.05 with a precision of ∼ 60%.
Supplementary Figure  example, one of the underlying assumptions is that the sample molecules are spherical and that they diffuse according to a free Brownian motion. The effect of anisotropic diffusion due to molecular non--sphericity is difficult to estimate. We can, however, qualitatively examine the effects of non--Brownian motion by considering the case of anomalous diffusion. In this case, the diffusion radius is given by where Γ is a scaling factor, and characterizes the diffusion process. In the case of a typical (Brownian) diffusion = 1, while > 1 is often associated with super--diffusion processes (e.g., due to cellular transport processes).
The change to super--diffusion leads to a faster signal decay, while maintaining an overall shape qualitatively similar to free diffusion. This leads to a certain ambiguity in determining the diffusion regime, as the faster signal decay can be partly compensated by reducing the value of Γ and/or !,!"#$ .

Role of NV spin--lattice relaxation
The observation of a correlation signal depends on the NV ability to store the phase picked up during the first interrogation segment throughout the evolution interval . A question of interest is, therefore, whether NV spin--lattice relaxation !!" must be taken into account in our modeling.
While we did not specifically determine !!" for the NVs we used in our correlation measurements, observations from virtually identical NVs are presented in Supplementary  Figure 9. The plotted data points correspond to the difference between pump--probe and inversion--recovery protocols. Out of the 25 NVs we studied, we determine an average spin--lattice relaxation time of ~350 µs with a maximum dispersion between !!" values of about 100 µs. These relaxation times are substantially longer than the time constants describing the decay of the correlation signal in the solid samples (~20 µs and ~30 µs for Samples A and B, respectively), confirming the above conclusion that NV relaxation can be ignored in these two cases. Sample C, on the other hand, deserves some special consideration: The first segment of the correlation signal (governed by molecular diffusion and decaying on a time scale of !"## (!!!"#)~1 5 µs) is clearly insensitive to NV spin lattice relaxation, meaning that no corrections are required in our numerical estimates of the bulk diffusion constant segment of the correlation signal, which decays on a time scale !"#! (!"#$)~2 50 µs comparable to !!" (Supplementary Figure 10b). Therefore, our estimate of the rotational diffusional constant in the adsorbed layer !,!"# must be interpreted as a lower bound.

Comparison with bulk NMR data
Relevant information on the sample dynamics near the diamond surface can be gained by comparing the measured correlation signals with bulk 1 H NMR data. Supplementary Figure  10 shows the NMR proton signal for Sample A (Merckoglas) and Sample B (PDMS) upon application of a π/2--pulse at 600 MHz (corresponding to a 14.1 T magnetic field). In both cases we find that the NV--detected correlation signals (Supplementary Figs. 1 and 2) are shorter--lived than their inductively detected counterparts ( 1 H NMR FIDs in Supplementary Figs. 10a and 10c). The difference is starker in the case of Sample B, where the nuclear spin coherence time is seen to persist beyond 0.5 ms. We interpret these differences as a manifestation of the singular molecular dynamics near the diamond surface, where motion is likely more restricted than in the bulk.
Indirect support for this idea is provided by the results of Supplementary Figs. 10b and 10d where we expose the importance of molecular packing in the dynamics of both Sample A and Sample B: In the first case ( Supplementary Figure 10b), we let a sample of 'as--purchased' Merckoglas dry to form a solid crust, which we then analyze via inductive NMR. This system differs from the one used in Supplementary Fig. 10a only in its preparation protocol: Unlike the practice in our NV experimentswhere Mercoglas is dissolved in toluene to form a 1:1 solution prior to gluing the diamond crystal to the sample holder-no solvent was added. The result is an 1 H FID substantially longer than that in Supplementary Figure 10a, the increased mobility likely originating from the reduced number of interstitial molecules after solvent evaporation. Interestingly, we observe the converse effect in Sample B where we shorten the FID duration by bringing the ratio between PDMS and the curing agent to a value lower than that used for our NV experiments. The latter leads to a more rigid form of the resulting polymer and thus a shorter--lived proton FID.
For completeness, Supplementary Figure 11 presents 1 H NMR data corresponding to Sample C (immersion oil), where we find a proton spin transverse relaxation time of order 25 ms, limited by field inhomogeneity (an FID from acetone is included as a reference).