Derivatives of spin dynamics simulations

Bayesian analysis of 1D 1H-NMR spectra

Extracting spin system parameters from 1D high resolution 1H-NMR spectra can be an intricate task requiring sophisticate methods. With a few exceptions methods to perform such a total line shape analysis commonly rely on local optimization techniques which for increasing complexity of the underlying spin system tend to reveal local solutions. In this work we propose a full Bayesian modeling approach based on a quantum mechanical model of the spin system. The Bayesian formalism provides a global optimization strategy which allows to efficiently include prior knowledge about the spin system or to incorporate additional constraints concerning the parameters of interest. The proposed algorithm has been tested on synthetic and real 1D 1H-NMR data for various spin systems with increasing complexity. The results show that the Bayesian algorithm provides accurate estimates even for complex spectra with many overlapping regions, and that it can cope with symmetry induced local minima. By providing an unbiased estimate of the model evidence the proposed algorithm furthermore offers a way to discriminate between different spin system candidates.

Optimal control gradient precision trade-offs: Application to fast generation of DeepControl libraries for MRI

We have recently demonstrated supervised deep learning methods for rapid generation of radiofrequency pulses in magnetic resonance imaging (https://doi.org/10.1002/mrm.27740, https://doi.org/10.1002/mrm.28667). Unlike the previous iterative optimization approaches, deep learning methods generate a pulse using a fixed number of floating-point operations - this is important in MRI, where patient-specific pulses preferably must be produced in real time. However, deep learning requires vast training libraries, which must be generated using the traditional methods, e.g., iterative quantum optimal control methods. Those methods are usually variations of gradient descent, and the calculation of the gradient of the performance metric with respect to the pulse waveform can be the most numerically intensive step. In this communication, we explore various ways in which the calculation of gradients in quantum optimal control theory may be accelerated. Four optimization avenues are explored: truncated commutator series expansions at zeroth and first order, a novel midpoint truncation scheme at first order, and the exact complex-step method. For the spin systems relevant to MRI, the first-order midpoint truncation is found to be sufficiently accurate, but also significantly faster than the machine precision gradient. This makes the generation of training databases for the machine learning methods considerably more realistic.

Tunable, Flexible, and Efficient Optimization of Control Pulses for Practical Qubits

Quantum computation places very stringent demands on gate fidelities, and experimental implementations require both the controls and the resultant dynamics to conform to hardware-specific constraints. Superconducting qubits present the additional requirement that pulses must have simple parameterizations, so they can be further calibrated in the experiment, to compensate for uncertainties in system parameters. Other quantum technologies, such as sensing, require extremely high fidelities. We present a novel, conceptually simple and easy-to-implement gradient-based optimal control technique named gradient optimization of analytic controls (GOAT), which satisfies all the above requirements, unlike previous approaches. To demonstrate GOAT's capabilities, with emphasis on flexibility and ease of subsequent calibration, we optimize fast coherence-limited pulses for two leading superconducting qubits architectures-flux-tunable transmons and fixed-frequency transmons with tunable couplers.

Rotation relaxation splitting for optimizing parallel RF excitation pulses with T1- and T2-relaxations in MRI

Exact solutions of the Bloch equations with T₁- and T₂-relaxation terms for piecewise constant magnetic fields are numerically challenging. We therefore investigate an approximation for the achieved magnetization in which rotations and relaxations are split into separate operations. We develop an estimate for its accuracy and explicit first and second order derivatives with respect to the complex excitation radio frequency voltages. In practice, the deviation between an exact solution of the Bloch equations and this rotation relaxation splitting approximation seems negligible. Its computation times are similar to exact solutions without relaxation terms. We apply the developed theory to numerically optimize radio frequency excitation waveforms with T₁- and T₂-relaxations in several examples.

First and second order derivatives for optimizing parallel RF excitation waveforms

For piecewise constant magnetic fields, the Bloch equations (without relaxation terms) can be solved explicitly. This way the magnetization created by an excitation pulse can be written as a concatenation of rotations applied to the initial magnetization. For fixed gradient trajectories, the problem of finding parallel RF waveforms, which minimize the difference between achieved and desired magnetization on a number of voxels, can thus be represented as a finite-dimensional minimization problem. We use quaternion calculus to formulate this optimization problem in the magnitude least squares variant and specify first and second order derivatives of the objective function. We obtain a small tip angle approximation as first order Taylor development from the first order derivatives and also develop algorithms for first and second order derivatives for this small tip angle approximation. All algorithms are accompanied by precise floating point operation counts to assess and compare the computational efforts. We have implemented these algorithms as callback functions of an interior-point solver. We have applied this numerical optimization method to example problems from the literature and report key observations.

Real-time 2D spatially selective MRI experiments: Comparative analysis of optimal control design methods

There is an increasing need for development of advanced radio-frequency (RF) pulse techniques in modern magnetic resonance imaging (MRI) systems driven by recent advancements in ultra-high magnetic field systems, new parallel transmit/receive coil designs, and accessible powerful computational facilities. 2D spatially selective RF pulses are an example of advanced pulses that have many applications of clinical relevance, e.g., reduced field of view imaging, and MR spectroscopy. The 2D spatially selective RF pulses are mostly generated and optimised with numerical methods that can handle vast controls and multiple constraints. With this study we aim at demonstrating that numerical, optimal control (OC) algorithms are efficient for the design of 2D spatially selective MRI experiments, when robustness towards e.g. field inhomogeneity is in focus. We have chosen three popular OC algorithms; two which are gradient-based, concurrent methods using first- and second-order derivatives, respectively; and a third that belongs to the sequential, monotonically convergent family. We used two experimental models: a water phantom, and an in vivo human head. Taking into consideration the challenging experimental setup, our analysis suggests the use of the sequential, monotonic approach and the second-order gradient-based approach as computational speed, experimental robustness, and image quality is key. All algorithms used in this work were implemented in the MATLAB environment and are freely available to the MRI community.

BEBE(tr) and BUBI: J-compensated concurrent shaped pulses for 1H-13C experiments

Shaped pulses designed for broadband excitation, inversion and refocusing are important tools in modern NMR spectroscopy to achieve robust pulse sequences especially in heteronuclear correlation experiments. A large variety of mostly computer-optimized pulse shapes exist for different desired bandwidths, available rf-field strengths, and tolerance to B1-inhomogeneity. They are usually derived for a single spin 1/2, neglecting evolution due to J-couplings. While pulses with constant resulting phase are selfcompensated for heteronuclear coupling evolution as long as they are applied exclusively on a single nucleus, the situation changes for concurrently applied pulse shapes. Using the example of a (1)H,(13)C two spin system, two J-compensated pulse pairs for the application in INEPT-type transfer elements were optimized: a point-to-point pulse sandwich called BEBE(tr), consisting of a broadband excitation and time-reversed excitation pulse, and a combined universal rotation and point-to-point pulse pair called BUBI, which acts as a refocusing pulse on (1)H and a corresponding inversion pulse on (13)C. After a derivation of quality factors and optimization protocols, a theoretical and experimental comparison with conventionally derived BEBOP, BIBOP, and BURBOP-180° pulses is given. While the overall transfer efficiency of a single pulse pair is only reduced by approximately 0.1%, resulting transfer to undesired coherences is reduced by several percent. In experiments this can lead to undesired phase distortions for pairs of uncompensated pulse shapes and even differences in signal intensities of 5-10% in HSQC and up to 68% in more complex COB-HSQC experiments.

Tailored real-time scaling of heteronuclear couplings

Heteronuclear couplings are a valuable source of molecular information, which is measured from the multiplet splittings of an NMR spectrum. Radiofrequency irradiation on one coupled nuclear spin allows to modify the effective coupling constant, scaling down the multiplet splittings in the spectrum observed at the resonance frequency of the other nuclear spin. Such decoupling sequences are often used to collapse a multiplet into a singlet and can therefore simplify NMR spectra significantly. Continuous-wave (cw) decoupling has an intrinsic non-linear offset dependence of the scaling of the effective J-coupling constant. Using optimal control pulse optimization, we show that virtually arbitrary off-resonance scaling of the J-coupling constant can be achieved. The new class of tailored decoupling pulses is named SHOT (Scaling of Heteronuclear couplings by Optimal Tracking). Complementing cw irradiation, SHOT pulses offer an alternative approach of encoding chemical shift information indirectly through off-resonance decoupling, which however makes it possible for the first time to achieve linear J scaling as a function of offset frequency. For a simple mixture of eight aromatic compounds, it is demonstrated experimentally that a 1D-SHOT {(1)H}-(13)C experiment yields comparable information to a 2D-HSQC and can give full assignment of all coupled spins.

Sensitivity of the NMR density matrix to pulse sequence parameters: a simplified analytic approach

We present a formalism for the analysis of sensitivity of nuclear magnetic resonance pulse sequences to variations of pulse sequence parameters, such as radiofrequency pulses, gradient pulses or evolution delays. The formalism enables the calculation of compact, analytic expressions for the derivatives of the density matrix and the observed signal with respect to the parameters varied. The analysis is based on two constructs computed in the course of modified density-matrix simulations: the error interrogation operators and error commutators. The approach presented is consequently named the Error Commutator Formalism (ECF). It is used to evaluate the sensitivity of the density matrix to parameter variation based on the simulations carried out for the ideal parameters, obviating the need for finite-difference calculations of signal errors. The ECF analysis therefore carries a computational cost comparable to a single density-matrix or product-operator simulation. Its application is illustrated using a number of examples from basic NMR spectroscopy. We show that the strength of the ECF is its ability to provide analytic insights into the propagation of errors through pulse sequences and the behaviour of signal errors under phase cycling. Furthermore, the approach is algorithmic and easily amenable to implementation in the form of a programming code. It is envisaged that it could be incorporated into standard NMR product-operator simulation packages.

Designing optimal universal pulses using second-order, large-scale, non-linear optimization

Recently, RF pulse design using first-order and quasi-second-order pulses has been actively investigated. We present a full second-order design method capable of incorporating relaxation, inhomogeneity in B(0) and B(1). Our model is formulated as a generic optimization problem making it easy to incorporate diverse pulse sequence features. To tame the computational cost, we present a method of calculating second derivatives in at most a constant multiple of the first derivative calculation time, this is further accelerated by using symbolic solutions of the Bloch equations. We illustrate the relative merits and performance of quasi-Newton and full second-order optimization with a series of examples, showing that even a pulse already optimized using other methods can be visibly improved. To be useful in CPMG experiments, a universal refocusing pulse should be independent of the delay time and insensitive of the relaxation time and RF inhomogeneity. We design such a pulse and show that, using it, we can obtain reliable R(2) measurements for offsets within ±γB(1). Finally, we compare our optimal refocusing pulse with other published refocusing pulses by doing CPMG experiments.

Implementing quantum gates by optimal control with doubly exponential convergence

We introduce a novel algorithm for the task of coherently controlling a quantum mechanical system to implement any chosen unitary dynamics. It performs faster than existing state of the art methods by 1 to 3 orders of magnitude (depending on which one we compare to), particularly for quantum information processing purposes. This substantially enhances the ability to both study the control capabilities of physical systems within their coherence times, and constrain solutions for control tasks to lie within experimentally feasible regions. Natural extensions of the algorithm are also discussed.

Molecular structure refinement by direct fitting of atomic coordinates to experimental ESR spectra

An attempt is made to bypass spectral analysis and fit internal coordinates of radicals directly to experimental liquid- and solid-state electron spin resonance (ESR) spectra. We take advantage of the recently introduced large-scale spin dynamics simulation algorithms and of the fact that the accuracy of quantum mechanical calculations of ESR parameters has improved to the point of quantitative correctness. Partial solutions are offered to the local minimum problem in spectral fitting and to the problem of spin interaction parameters (hyperfine couplings, chemical shifts, etc.) being very sensitive to vibrational excursions from the equilibrium geometry.

Second order gradient ascent pulse engineering

We report some improvements to the gradient ascent pulse engineering (GRAPE) algorithm for optimal control of spin ensembles and other quantum systems. These include more accurate gradients, convergence acceleration using the Broyden-Fletcher-Goldfarb-Shanno (BFGS) quasi-Newton algorithm as well as faster control derivative calculation algorithms. In all test systems, the wall clock time and the convergence rates show a considerable improvement over the approximate gradient ascent.

Fast and accurate algorithm for the simulation of NMR spectra of large spin systems

The computational cost for the simulation of NMR spectra grows exponentially with the number of nuclei. Today, the memory available to store the Hamiltonian limits the size of the system that can be studied. Modern computers enable to tackle systems containing up to 13 spins [1], which obviously does not allow to study most molecules of interest in research. This issue can be addressed by identifying groups of spins or fragments that are not or only weakly interacting together, i.e., that only share weakly coupled spin pairs. Such a fragmentation is only permitted in the weak coupling regime, i.e., when the coupling interaction is weak compared to the difference in chemical shift of the coupled spins. Here, we propose a procedure that removes weak coupling interactions in order to split the spin system efficiently and to correct a posteriori for the effect of the neglected couplings. This approach yields accurate spectra when the adequate interactions are removed, i.e., between spins only involved in weak coupling interactions, but fails otherwise. As a result, the computational time for the simulation of 1D spectra grows linearly with the size of the spin system.

Spinach--a software library for simulation of spin dynamics in large spin systems

We introduce a software library incorporating our recent research into efficient simulation algorithms for large spin systems. Liouville space simulations (including symmetry, relaxation and chemical kinetics) of most liquid-state NMR experiments on 40+ spin systems can now be performed without effort on a desktop workstation. Much progress has also been made with improving the efficiency of ESR, solid state NMR and Spin Chemistry simulations. Spinach is available for download at http://spindynamics.org.