Christina Graf1, Armin Rund2, Christoph Stefan Aigner1, Karl Kunisch2, and Rudolf Stollberger1
1Institute of Medical Engineering, Graz University of Technology, Graz, Austria, 2Institute of Mathematics and Scientific Computing, Karl-Franzens-University Graz, Graz, Austria
Synopsis
Accurate Bloch and Bloch-McConnell model simulation is essential for proper modelling and becomes crucial for iterative RF pulse design methods. This work demonstrates the advantages of a piecewise analytical solution using eigenvalues and eigenvectors and an asymmetric and symmetric operator splitting scheme. Those three solvers are compared to five other numerical solution methods regarding accuracy and runtime. The symmetric operator splitting prevails in both compared to the others. Moreover, it is quadratic convergent with respect to the time step and therefore seems to be a good choice for optimal control approaches for Bloch and for Bloch-McConnell models.
INTRODUCTION
Iterative RF pulse optimization1,2 requires multiple Bloch or extended Bloch model simulations with the need for fast and accurate solution methods. Especially in the case of chemical exchange saturation transfer (CEST) the typically long pulse durations are a challenge for iterative optimization methods3. In this work, different analytical and numerical methods for solving the Bloch and Bloch-McConnell equations are derived and compared in extensive numerical studies regarding accuracy and runtime. THEORY
Both the Bloch and Bloch-McConnell equations can mathematically be described as a system of differential equations $$\dot{M}(t) = A(t) \cdot M(t) + b, \: \: \: M(0)=M^{0} \: \: \text{for} \: \: t \in (0,T).$$
Therein, M(t) is the magnetization vector, M0 the initial magnetization and b the inhomogeneity including the longitudinal relaxation time. The coefficient matrix A(t) is non-constant, therefore a full analytic solution does not exist in general. However, if a piecewise constant discretization is chosen for RF pulse and slice-selective gradient (Gs) amplitude, an analytical solution per time step can then be derived using eigenvalues and eigenvectors (EV). On the other hand, by splitting the system matrices into relaxation/chemical exchange and rotation parts, an asymmetric (ASY)4 and a symmetric (SY) operator splitting scheme can be defined, where all subsystems are solved exact. The splitting error is known to be of first-order for asymmetric and of second-order for symmetric operator splitting.
METHODS
The three methods (EV, ASY and SY) are implemented in MATLAB (The MathWorks, Inc., Natick, USA) and are analyzed in extensive numerical studies, both for the Bloch equations and the Bloch-McConnell equations. For the former, the methods are compared to five other numerical solution methods including the spin domain solution with neglected relaxation terms SD5, the MATLAB built-in functions ode45 and expM, an implicit Euler scheme IE and the Crank-Nicholson method CN. The numerical comparison is performed for a complex 180o RF pulse and a trapezoidal Gs shape with a duration of T=1.9ms discretized with 1104 points. The spatial domain is set to (-50,50)mm using 101 points. Different relaxation times are investigated, T1=T2=$$$\infty$$$, T1=1331ms, T2=80ms, and T1=400ms, T2=5ms. The proposed Bloch-McConnell simulation methods (ASY and SY) are compared to the analytical EV solver for a series of Gaussian and block RF pulses with a total duration of T=2100ms discretized with 4421 points in the first place for an amide proton system in white matter6. Later, a convergence study with respect to time discretization is performed. The z-spectrum (-5,5)ppm is discretized using 101 points. The relaxation times of the water proton pool are set to T1,w=1048ms and T2,w=69ms and of the solute proton pool to T1,s=1048ms and T2,s=15ms. RESULTS AND DISCUSSION
Figure 1 and 2 show the accuracy of all Bloch solvers. In absence of relaxation effects, most solvers yield exact results in machine precision, larger errors are observed for CN and IE. A closer look reveals that IE results in a magnetization degradation over time (not shown). Furthermore, we analyze the behavior in the presence of relaxation effects. IE, CN, expM and SD yield larger errors, the latter since it does not account for relaxation effects. In contrast, ode45, ASY and SY perform with higher precision. Regarding runtime, ASY is slightly faster than SY, which in turn is 3 orders of magnitude faster than ode45 (for the chosen accuracy), and EV. In general, SY outperforms ASY in the accuracy for the chosen time-step. This is also true in presence of chemical exchange saturation transfer shown in Figure 3, 4 and 5. Figure 5 lists the numerically verified linear and quadratic convergence with respect to the time step for ASY and SY in high precision for the series of RF block pulses.CONCLUSION
We presented a piecewise analytical solver EV for the Bloch and Bloch-McConnell equations that allows exact simulation at the cost of an increased runtime. The symmetric operator splitting SY shows an excellent accuracy together with quadratic convergence in the time-step. Since its runtime is only slightly increased compared to ASY, it seems to be a good choice for optimal control approaches for Bloch as well as for Bloch-McConnell models.Acknowledgements
Supported by the Austrian Science Fund (FWF): SFB F32‐N18 and BioTechMed-Graz.References
1) Christoph Aigner, Christian Clason, Armin Rund, Rudolf Stollberger. Efficient high-resolution RF pulse design applied to simultaneous multi-slice excitation. J Magn Reson., 263:33-44, 2016.
2) Armin Rund, Christoph Aigner, Karl Kunisch, Rudolf Stollberger. Simultaneous multislice refocusing via time optimal control. Magn. Reson. Med. 80(4):1416-1428, 2018.
3) Gang Xiao, Renhua Wu, Phillip Zhe Sun. Fast simulation and optimization of pulse-train chemical exchange saturation transfer (CEST) imaging. Phys Med Biol, 60(12):4719-4730, 2015.
4) Jacques Bittoun, Jacques Taquin, Michel D. Sauzade. A computer algorithm for the simulation of any nuclear magnetic resonance (NMR) method. J Magn Reson, 2(2):113-120, 1984.
5) Matt A. Bernstein, Kevin F. King, Xiaohong Joe Zhou. Handbook of MRI Pulse Sequences. Elsevier Academic Press, 2004.
6) Moritz Zaiss and Peter Bachert. Exchange-dependent relaxation in the rotating frame for slow and intermediate exchange - modeling off-resonant spin-lock and chemical exchange saturation transfer. NMR in Biomedicine, 26:5, 2013.