Hybrid quantum-classical simulations of magic angle spinning dynamic nuclear polarization in very large spin systems

The sensitivity of nuclear magnetic resonance (NMR) spectroscopy can be enhanced by transferring polarization from more magnetic electron spins, an approach known as dynamic nuclear polarization (DNP).^1–3 The last decade has seen a dramatic rise in magic angle spinning (MAS)-DNP, in large part due to the development and commercialization of high-frequency gyrotron microwave (MW) sources.^4,5 With the possible exception of the gyrotron, the greatest leaps in sensitivity have been made through the development of improved radical polarization sources. The most notable of which have been the invention of bisnitroxide polarizing agents that facilitate electron-to-nuclear polarization transfers by way of the three-spin cross-effect (CE) mechanism,^6–8 the development of large molecular weight bisnitroxides that have increased relaxation times,^9–11 and, most recently, the exploration of asymmetric biradicals.^12–14 

Theory has played a large role in the development of new polarization sources. For instance, early on, it was realized that biradicals needed to be designed with relative radical orientations that limited the number of inactive spin pairs.^15,16 A deeper understanding of the cross-effect mechanism during MAS revealed new ways to design DNP polarizing agents. Under MAS, the electron paramagnetic resonance (EPR) frequencies of the radicals are modulated due to the g-anisotropy. This rapid change in EPR frequency can lead to four types of anticrossing events,^17–19 namely, (1) microwave (MW) events wherein a radical’s EPR frequency crosses the frequency of the applied microwaves, (2) dipolar and exchange (D/J) events that occur when the EPR frequencies of two radicals cross one another, (3) cross-effect⁶ (CE) events that occur when the difference between two radicals’ EPR frequencies equals the Larmor frequency of the nuclear spins, and finally, (4) solid-effect^20,21 (SE) events that occur when the microwave frequency equals the sum or difference of the EPR and NMR frequencies. MW events lead to the saturation or inversion of an EPR resonance, D/J events lead to the exchange of magnetization between electron spins, CE events transfer magnetization to the nuclear spins of up to the difference in magnetization between the two electrons, and SE events lead to simultaneous flips of electron and nuclear spins.

Theoretical investigations revealed a particularly notorious interplay between the D/J and CE events.^17,18 Given that D/J events are never 100% efficient, particularly for inter-molecular events, they can equilibrate the polarization of the electron spins, leading to the draining of nuclear magnetization via the cross-effect. This effect has been termed electron-induced depolarization and often leads to the overestimation of DNP enhancement factors.^22,23 Inefficient D/J events thus lead to a lower than optimal difference in electron polarization and an overall lower level of hyperpolarization. This revelation has guided the development of improved polarizing agents that aimed to either increase the efficiency of D/J events^24,25 or avoid them altogether.^{12–14,26–28}

The role that internuclear interactions play in the overall DNP process is relatively underexplored. It is well understood that spin diffusion mediates the transfer of polarization to the observable nuclei;^29–32 however, experiments show that depleting the nuclear spin bath in DNP surface-enhanced NMR spectroscopy (SENS) experiments counterintuitively leads to an increase in DNP enhancements,³³ as does the perdeuteration of the polarizing agent,^34,35 however, at the cost of a slower DNP build-up rate. It is thought that the increased performance in DNP SENS when using deuterated solvents is made possible by radicals adsorbing onto surfaces,³⁶ coupled with a reduction in the overall number of spins needing to be polarized. Theoretical models to probe such mechanisms in silico are unfortunately difficult to design.

While MAS-DNP processes in small spin systems composed of two electrons and a handful of nuclei can be simulated exactly in Liouville space, the exponential scaling and need for very fine time steps of nanosecond order prevent the expansion of such methods to large spin systems.^37–41 The most successful solution to this dilemma has been the use of classical Fickian diffusion to model the transfer of polarization either away from a point source (a radical molecule) or toward the interior of a nanoparticle.^27,42–46 Such models have yielded a greater understanding of experimental results but do not provide a means to test new experiments or polarizing agents outside of the laboratory.

To increase the scale of predictive ab initio DNP simulations, Köckenberger proposed the use of restricted state-space methods,^47,48 first focusing on static SE DNP. We later expanded this work to MAS-DNP by combining a restricted state-space approach to simulating spin diffusion, known as low-order correlations in Liouville space (LCL),^50–52 with a Landau–Zener treatment of the electron spin dynamics. This increased the scale of the simulations to ∼50 nuclei, at the expense of the propagator, given that it would overrun the memory of the computer.⁵³ Since the propagator could not be reused to calculate arbitrary time points, an optimization routine was used to find the steady-state DNP enhancements.⁵³ The use of more aggressive locally restricted (LR) basis sets,⁵⁴ similar to the IK-1 basis set implemented in Spinach,^47,55 enabled for the increase of the spin system size to ∼1000 spins in ab initio simulations.⁵⁶ These simulations revealed a steep decline in DNP enhancements within the first nanometer from the polarizing agent due to the so-called spin diffusion barrier.⁵⁷ This abrupt loss in hyperpolarization was found to be dependent on the geometry of the biradical molecule, suggesting the possibility of engineering spin diffusion highways to direct polarization away from the radical.

Despite these successes, LR-LCL-based multinuclear MAS-DNP simulations are highly computationally demanding, requiring on average a week of computation on a computer cluster. Given that the bulk of the compute time is taken up by the treatment of the spin diffusion, we wondered whether this task could be replaced by a rate matrix. Such an approach has been successful in describing proton spin diffusion under MAS⁵⁸ and has also been applied to large-scale MAS-DNP simulations.⁵⁹ Rate constants could then be parameterized to reproduce the ab initio result in a small fraction of the time. This idea is closely related to the kinetic Monte Carlo methods developed by Köckenberger for describing static DNP^61–62 and Corzilius’ treatment of spin diffusion in MAS-DNP.³¹ In these methods, however, estimates of the exchange rate constants are calculated from first-principles. The dramatic reduction in matrix sizes would once again enable the use of propagators and the calculation of DNP build-up times in addition to the steady-state enhancements. Here, we will describe the parameterization of the spin diffusion rate matrix and apply the new hybrid quantum-classical MAS-DNP model to simulate the MAS rate and concentration dependence of DNP enhancements and build-up times. Then, we apply this model to study the transfer of polarization into core–shell metal-organic-framework (MOF) nanoparticles.

An ensemble of spins can be described using the density operator, $σ̂$⁠. Its time evolution, due to the spin Hamiltonian $Ĥ(t)$⁠, is described using the Liouville–von-Neumann equation,

The commutator between the Hamiltonian and the density operator is known as the Liouvillian (⁠$L̂̂$⁠), and the Hamiltonian is generally expressed as a sum of independent Hamiltonian operators that describe the various interactions affecting the spins. In the case of MAS-DNP, the Hamiltonian takes the following form in the microwave rotating frame:

where

$ĤZ$ is the Zeeman Hamiltonian and describes the interaction between the spins and the magnetic field. In this expression, g_zz,i is the element of the “i” electron’s g tensor that is aligned with the magnetic field, μ_B is the Bohr magneton, ω_MW is the microwave frequency, σ_zz,j is the component of the “j” nucleus’ magnetic shielding tensor that is aligned with the magnetic field, and ω_n is the nuclear Larmor frequency. $Ĥe,e$⁠, $Ĥe,n$⁠, and $Ĥn,n$ describe the spin–spin coupling interactions between the electrons and nuclei either through space (d, A, and D) or through the exchange mechanism in the case of the electrons (J). Finally, $ĤMW$ describes the interaction between the electrons and the microwaves, which has a (usually constant) Rabi frequency of ω₁.

The evolution of $σ̂$ in Liouville space uses a propagator, $Û̂$⁠,

which is time-dependent due to the MAS rotation. As such, it is necessary to divide time into P short increments of duration $Δ$t over which $Ĥ(t)$ can be assumed to be time-independent,

where $T̂̂$ is the Dyson time-ordering operator. If $Δ$t is small, one can also apply the Suziki–Trotter approximation and propagate each of the Hamiltonians independently. To first order, this corresponds to⁶³ 

The exact solution of Eqs. (9) or (10) is typically beyond reach in all but the simplest of spin systems both due to the exponential scaling of $σ̂$ and the short time steps needed to treat the electron Hamiltonian as time-independent. The following sections will describe the approximations that are made to accelerate the treatment of the nuclei (Sec. II A) and electrons (Sec. II B), the optimization algorithm used to determine the steady-state nuclear polarizations (Sec. II C), and finally the hybrid quantum-classical model (Sec. II D).

The density operator is expressed in a basis set consisting of direct products of single-spin operators,^47,50

where b_r is the coefficient for the basis operator $B̂r$ and

In the above summation, q_r corresponds to the spin order of the particular basis operator, which is the number of nuclei involved in the product. $Îi,r$ is the operator associated with spin “i,” taking values of

where $Êi$ is the identity operator.

In this basis, it is possible to distinguish operators by their spin order and coherence order and by the identity of the spins involved in the product. Given that no pulses are applied to the nuclei, we can freely remove the non-zero-quantum operators from the basis set without sacrificing accuracy.⁶⁴ Given that we can only directly observe single-spin operators and that operators with high-q_r must sequentially decrease in q_r to contribute to an observable operator (q_r = 1), it also follows that operators with low q_r values have a larger impact on the spin dynamics.⁶⁵ Indeed, it has been demonstrated, through comparisons with exact simulations, that limiting q_r to a maximum value (q_max) of 4 is sufficient to reproduce an exact simulation in the case of MAS spin diffusion⁵⁰ and solid effect DNP.⁶⁷ These two state-space reductions describe the LCL basis set, which scales polynomially in size,^47,50 rather than exponentially. At the limit where q_max is equal to the number of nuclei in system (N), the LCL result corresponds to the exact result.

Edwards et al. demonstrated that it is possible to further eliminate product operators involving distant nuclei for the calculation of spin diffusion in liquids without much loss in accuracy,⁶⁶ an approach which we later adapted for the treatment of spin diffusion in rotating solids.⁵⁴ In this locally restricted (LR)-LCL basis set, basis operators are only kept if they involve spins that are among the N_max nearest neighbors of a given spin. A N_max value of 27 is usually sufficient to reproduce the LCL result. Note that if N_max equals N − 1, it corresponds to the LCL basis set.

High computational efficiency in LR-LCL calculations can be achieved with the use of a hash map to store the density operator and then propagate one dipolar interaction at a time using a Suzuki–Trotter propagation,⁵⁰ 

There, a list of the indices of a spin’s nearest neighbors, n_i, is used. In practice, the propagation is achieved with the use of a series of four-dimensional subspaces with the corresponding rotation matrices

for subspaces of the form [Î_iz$B̂r$⁠, Î_jz$B̂r$⁠, I_i+Î_j−$B̂r$⁠, Î_i−Î_j+$B̂r$] and

for subspaces of the form [Î_i±$B̂r$⁠, Î_j±$B̂r$⁠, Î_i±Î_jz$B̂r$⁠, Î_j±Î_i$B̂r$], as described previously,⁵⁰ where $B̂r$ is an operator excluding spins i and j and ω_D,i,j is the instantaneous dipolar frequency for the interaction between spins i and j,

R_DD,i,j is the dipolar coupling constant between spins i and j, and three sets of Euler angles in Eq. (17) are used to rotate the spin system in the rotor around the magic angle, from the rotor frame to the interaction frame and, finally, to the lab frame.

The chemical shift propagator in this basis is a diagonal matrix with elements of

where p_i is the coherence order of the i spin’s operator in the basis operator (i.e., 0 for $Êi$ and $Îi,z$⁠, 1 for $Îi,+2$⁠, and −1 for $Îi,−2$⁠). Equation (18) can also be used to model the effects of hyperfine interactions on the spin diffusion rates by applying an offset corresponding to the strength of the electron-nuclear dipolar coupling at that particular orientation.

An important thing to note is that the number of 4 × 4 operations that need to be performed [Eqs. (14)–(16)] is proportional to N, and as such, this type of calculation scales linearly with the size of the spin system, enabling ab initio calculations of spin diffusion in very large spin systems.⁵⁴ 

The explicit treatment of the electron spins has a significant impact on DNP calculations, given that each addition of an electron increases the basis set size by a factor of 4. Fortunately, electron non-zero quantum coherences are largely unpopulated and can be eliminated from the basis set.^47,67 Further simplification was made possible in 2012 by Thurber and Tycko who showed that polarization transfers involving electrons occurred in extremely short level anticrossing events brought forth by the MAS rotation.¹⁷ By directly calculating the polarization transfers occurring during these events, using the Landau–Zener formula, and making the assumption that the crossings are infinitely short, the electrons could be effectively ignored, save for relaxation, for the bulk of the calculation with minimal losses in accuracy.

Mentink-Vigier later refined this work by representing the Landau–Zener events in the matrix form, which enabled longer calculations that take advantage of the periodicity of the MAS rotation and propagators.^18,58 We will reiterate the different rotor events that can occur along with their corresponding propagators. Generally, such a calculation involves the pre-computation of the MAS dependence of the electron Larmor frequencies, ω_ei, to identify the time steps when a level anticrossing event occurs. The velocities of these events are stored to memory, along with the indices corresponding to the time increments when they occur for use in the later calculation.

When two electron spins’ EPR frequencies cross one another such that ω_ei = ω_ej, they will exchange magnetization in a quantity that is proportional to the strength of the dipolar (D) and exchange (J) interactions between them and inversely proportional to the rate of the following crossing:

where

if a D/J event occurs and is otherwise zero.

Note that in Eq. (19), the basis operators are ordered as Ê, followed by the electron Î_z operators, and then the nuclear operators. Equation (19) corresponds specifically to a D/J event between electrons 1 and 2.

When the EPR frequency of an electron comes into resonance with the frequency of the applied microwaves (ω_e,i = ω_MW), it is partially saturated or inverted. The corresponding propagator for such an event is as follows in the case of electron 1:

where

if a MW event occurs and is otherwise zero.

When the microwave frequency matches either the electron-nuclear zero-quantum or double-quantum frequencies (ω_ei ± ω_nj = ω_MW), a solid effect event occurs, and polarization is transferred from the electron to the nuclei. The propagator for such an event is as follows in the case that the irradiated solid-effect condition involves electron 1:

where

when a SE event occurs and is otherwise zero.

Finally, if two electrons have a difference in EPR frequencies that equals the nuclear Larmor frequency [±(ω_ei − ω_ej) = ω_n], a cross effect event occurs and a portion of the difference in electron polarization is transferred to the nuclear spins. As an example, if we consider a cross-effect event between electrons 1 and 2 and nucleus 1, the propagator would equal,

where

when a CE event occurs and is otherwise zero.

We can thus define a global Landau–Zener propagator, which describes all interactions involving the N_el electrons,

MAS-DNP simulations are generally performed using the following strategy. We first simulate the evolution of the spin system for the duration of a rotor period and store the propagator,

where

and b_r,eq is the equilibrium coefficient of a given operator: the relative Boltzmann polarization for I_z-type operators, otherwise zero, and T_n is the relaxation rate of a given operator in the product: T₁ for I_z-type operators and T₂ for raising and lowering operators. $Û̂r$ can then be exponentiated to calculate the spin polarizations at a later time point and predict DNP build-up curves, including the steady-state that yields the DNP enhancement factors in a fully relaxed spin system,

A problematic feature of this approach is that it requires the use of matrices with dimensions of n, where n is the number of basis operators. In particular, propagators are dense complex matrices and thus require the storage of 2n² double-precision floating-point numbers. As such, a basis set only needs to contain 11k vectors for the propagator to cost over 2 GB to store, which corresponds to only 11 spins for LCL simulations. It is, nevertheless, possible to simulate far larger basis sets by working in low-dimensional subspaces via Suzuki–Trotter propagation (see Sec. II A), but these approaches are limited in timespan, given that they cannot take advantage of the periodicity of the rotor [Eq. (28)]. For instance, a 10 s simulation takes ten times longer to perform than a 1 s simulation.

It is, nevertheless, possible to exploit rotor periodicity in large-scale calculations of steady-state DNP enhancements by considering the following relationship:⁵³ 

Namely, at the steady-state, the DNP enhancements are invariant at the start and end of a rotor period. The calculation of DNP enhancements can thus be rephrased into a search for a set of nuclear polarizations that satisfy Eq. (31). Whether the nuclear polarizations increase or decrease after a rotor period reflects whether the initial guess of nuclear polarizations was too low or high. We first utilized this feature, in combination with LCL, to calculate DNP enhancements in moderately large spin systems using a simple Monte Carlo optimization algorithm where the nuclear polarizations, P, are stepped in their direction of change by a random step size between 0 and λ, tabulated in vector $Λ$⁠,⁵³ 

In Eq. (32), the subscript m corresponds to the optimization step. The procedure for calculating steady-state DNP enhancements is then to propagate the density operator from an arbitrary initial guess for a period of t_r. The polarizations are then updated according to Eq. (32), and the calculation is continued. This process is repeated until the $Δ$P values are below a certain threshold. Assuming that λ is sufficiently small, this process is guaranteed to converge to the steady-state DNP enhancements.

We, however, found Eq. (32) to scale poorly in very large spin systems due to the time it takes for a change in polarization from a given nucleus to affect that from a distant one. Far faster convergence is achieved if, in addition to the random jumps from Eq. (32), all nuclear polarizations are stepped in unison in a manner that is proportional to the total change in the nuclear polarizations,⁵⁶ 

where the factor μ corresponds to the step size. The same process can also be applied to smaller groups of neighboring N_max spins, leading to

where

The step sizes λ, ν, and μ are set to obtain the fastest convergence; for the DNP enhancement of ¹H spins, we have found that the following values work well: λ = 0.5, ν = 0.25, and μ = 0.5, with the value of λ decreasing by a factor of 0.999 at each optimization step. Note that the step sizes are normalized to the nuclear Boltzmann polarization.

The combination of LR-LCL and the Landau–Zener model for MAS-DNP has enabled the calculation of DNP enhancements in periodic spin systems containing upward of 1000 spins through the use of the optimization approaches described in Sec. II C.⁵⁶ In this case, the propagation can be expressed as

which uses exclusively matrix-vector operations in small subspaces of up to four dimensions. By far, the most expensive part of the calculation is the LR-LCL part, which is needed to treat spin diffusion and the effects that hyperfine interactions have on spin diffusion rates (i.e., the so-called spin diffusion barrier). If the LR-LCL calculation could be replaced with a kinetic model, trained to reproduce the LR-LCL result, the space size could be dramatically reduced by a factor of about 2000 since the basis would only effectively contain the Î_z operators and the identity operator. In this case, Eq. (36) would be simplified to

Here, $K̂̂$ is a matrix containing the spin diffusion rates between the different nuclei in the system, k_i,j. The reduced matrix dimensions also enable the application of Eq. (30) and the calculation of DNP build-up curves in large atomistic models, hopefully, with the same level of accuracy as obtained using the ab initio model.

The application of the simulation approach summarized in Eq. (37) requires that the elements of $K̂̂$⁠, namely, the spin diffusion rate constants k_i,j, be predicted accurately. From perturbation theory, we expect that these rate constants will be proportional to $RDD2$^69–71 and inversely proportional to the MAS frequency, ν_r.^72,73 The spin diffusion rates should also have some sort of inverse relationship with the strength of the hyperfine interactions, A_zz,PAS, with the different electrons in the sample.⁷⁴ We, thus, tested the following expression for the calculation of the rate constants:

The expression has three variable parameters, A, B, and C, which can be tuned to improve the agreement between an LR-LCL-based ab initio calculation and the hybrid model. We first optimized the value of A by setting B to zero and comparing the result to an ab initio calculation that does not apply a hyperfine shift operator. In other words, the spin diffusion barrier is removed. Next, B and C are optimized together against simulations that included the spin diffusion barrier. We were able to obtain a good agreement between the hybrid and ab initio simulations with values of A = 0.000 21, B = 0.0027 s C⁻¹, and C = 1.2; see Fig. 1. The fitting was performed with the MAS rate in a range of 10–40 kHz on 16 mM periodic solutions of TEKPol¹¹ and bTbK¹⁶ in 1,1,2,2-tetrachloroethane (TCE).⁵⁶ 

Given that these parameters agree well with the ab initio results in TCE solutions of two different radicals, we expect that they are relatively reliable for TCE solutions, in general. In systems involving different NMR-active isotopes or differences in chemical shifts, it may be necessary to refine the fitted parameters or alter Eq. (38) to suit the particular situation. In certain circumstances, it may also be necessary to include high-order spin diffusion terms, for instance, at very fast MAS frequencies.^75,76 In any case, we expect that ab initio simulations of DNP enhancements could be performed on all or part of a system of interest and used to refine Eq. (38) for the particular case, before calculating the build-up rates, or a larger spin system with the hybrid model.

Having optimized the hybrid MAS-DNP model for the reproduction of ab initio results in TCE solutions, we were interested to see whether it could be used to predict DNP enhancements and build-up curves as the radical concentration and MAS rate are altered. In particular, our previous ab initio simulations were limited to periodic models containing only a single biradical molecule. As such, it was not possible to properly assess the concentration dependence, given that inter-radical dipolar events could not be included in the simulation. In the following sections,^1–3 we will compare the predicted DNP behavior of single- and multi-biradical solutions to that obtained from the experiment.

We previously used molecular dynamics simulations of biradical solutions to obtain representative ensembles of radical conformations and solvent proton distributions, both of which are expected to impact the calculated DNP enhancements.⁵⁶ As we are now interested in simulating DNP in solutions containing upward of 20 000 atoms, this approach proved to be impractical. As such, we instead opted to use the software Balloon⁷⁷ to search for the lowest-energy TEKPol conformers using a multiobjective genetic algorithm. 86 conformers were identified, and the five lowest-energy, non-duplicate, conformers were optimized using Density functional theory (DFT) in a polarizable continuum model (COSMO, 1,1,2,2-tetrachloroethane).⁷⁸ The energy difference between the first and fifth lowest-energy conformers was −2.1 kcal/mol. Using these computed energies, we calculated the expected Boltzmann populations of the different conformers in solution at 100 K (Fig. 2).

To generate a simulation model, we then defined a periodic unit cell with its dimensions set such that the radical concentration at room temperature would equal 2, 4, 8, 16, or 32 mM, with a given number of biradical molecules from 1 to 6. Radical molecules were inserted into the unit cell in randomized locations, with randomized orientations. The remainder of the cell was then populated with randomly located solvent protons with a density matching that found in frozen TCE.⁷⁹ No interatomic contacts below 2 Å were allowed in the model generation.

LR-LCL calculations scale linearly with the number of nuclei in the spin system since the number of four-dimensional matrix-vector products also scales linearly. In practice, optimization-based DNP calculations scale slightly worse due to the need for performing a greater number of optimization steps in larger spin systems (for example, 250 steps for 1250 nuclei vs 20 steps for 12 nuclei). A 1250-nuclei, two-electron, LR-LCL MAS-DNP calculation may take on the order of 2 days to complete on an AMD Epyc 7662 processor, but the same spin system can be calculated with the hybrid approach [Eq. (37)] in only 10 min. The hybrid calculations, however, scale as O(N³) due to the use of matrix–matrix products for the calculation of a rotor-periodic propagator. This places the practical limit of this simulation approach to roughly 10 000 nuclei, close to the memory limit for the storage of the propagator.

The main bottleneck of the calculations, however, shifts from the treatment of spin diffusion to the application of the numerous cross-effect propagators. As the number of electron spins in the system increases, so does the number of cross-effect events as well, and as such, the calculation also scales quadratically in the number of electron spins.

One way to accelerate the calculations is to reduce the number of cross-effect events such that electrons will only transfer polarization to their immediate neighbors. This concept is based on the fact that the efficiency of DNP transfers has an r⁻⁶ dependence, see Eq. (20), while the number of nuclei at a given distance from the radical only has an r² dependence. While this approach has been used before in DNP simulations,⁵⁸ it was not clear how many electron-nuclear interactions are needed to retain the accuracy of the calculation. We performed MAS-DNP simulations in a 2 mM TEKPol/TCE solution (2 electrons, 9544 nuclei) as a function of the number of nuclei to which each electron could directly transfer polarization using the cross-effect. If each electron was only allowed to interact with its closest nuclear spin, then the predicted DNP enhancements were lower by 72%. This error is quickly reduced, however, with the inclusion of more nuclei and is below 1% when the electrons each interact with 35 nuclear spins. The total computational time was reduced by a factor of 10 (see Fig. 3).

It would, nevertheless, appear that a full half of the hyperpolarization that is fed to the 9544 nuclei in the spin system was initially generated on only four nuclear spins. These nuclei are situated on the biradical molecule itself, which highlights the importance of properly situating these innermost nuclei and providing them with efficient spin diffusion paths away from the electrons and into the bulk.⁸⁰ 

DNP enhancements, relative depolarization levels, and build-up times were measured as a function of the MAS frequency, from 10 to 39 kHz in 5 kHz increments, for TEKPol/TCE solutions with biradical concentrations of 0, 2, 4, 8, 16, and 32 mM. The results are tabulated in Table I. I_on and I_off are used to refer to the integrated signal volumes in the presence and absence of microwaves, respectively. I_on,10kHz and I_off,10kHz correspond to the signal intensities at the lowest MAS frequency, 10 kHz, and are used to follow the MAS-dependence of the enhancement and depolarization factors. T_on and T_off are the build-up time constants in the presence and absence of microwaves, respectively. Due to challenges associated with introducing controlled quantities of the solution into the small 1.3 mm rotors, only relative values are reported.

Generally, we observe that while the vertical signal intensity increases with increasing MAS rate, the integrated intensity is fairly constant in the absence of microwaves and radicals. Magnetization losses as the MAS rate is increased, both in the presence and absence of microwaves, are seen to increase slightly with the radical concentration. This signal decrease is similar at the lower radical concentrations but is worse in the presence of microwaves in the higher-concentration samples, leading to a moderate overall reduction in the DNP enhancement factors (I_on/I_off). Counter to prior reports where narrow-line radicals were used,^27,81 we did not observe an increase in DNP efficiency as the MAS rate was increased when using a bisnitroxide. Integrated intensities generally decreased by 25%–40% as the MAS rate was increased from 10 to 39 kHz. Vertical signal intensities, however, did increase, although this was caused by the overall reduction in linewidths (see Table I, Fig. S3).

In all cases, build-up time constants were seen to increase with the MAS frequency and decrease with the radical concentration. Similar time constants were measured in both the presence and absence of microwaves, suggesting that cross-effect-induced nuclear (de)polarization dominates the relaxation processes in both cases. In the absence of radicals, a strong increase in build-up times was again observed as the MAS rate was increased, in agreement with the source–sink model introduced by Chaudhari and co-workers (see Table I).²⁷ 

We then performed MAS-DNP calculations, both in the presence and absence of microwaves and in models containing single or multiple biradical molecules per unit cell, to determine how inter-radical interactions impact the accuracy of the simulations. Nuclear relaxation times were set according to the following expression:

with T_1,bulk set to the value measured at the corresponding MAS rate in the 0 mM solution. The electron relaxation times were set to the values measured experimentally on each of the solutions using pulsed X-band EPR spectroscopy (see the supplementary material). The field and orientational dependences of the electron relaxation times were neglected, which may lead to some errors in the calculated DNP behavior.^38,82 The microwave power was set to 0.85 MHz, in accordance with previous studies.^38,53,83 The results are tabulated in Tables II and III for the single- and multi-biradical models, respectively, and correspond to the average over five randomly generated solutions of a given radical concentration.