开放式磁共振系统中梯度磁场的快速计算方法

朱旭晨， 王秋良， 王厚生

(1. 中国科学院电工研究所， 北京 100190； 2. 中国科学院大学， 北京 100049)

开放式磁共振系统中梯度磁场的快速计算方法

朱旭晨1,2， 王秋良1， 王厚生1

(1. 中国科学院电工研究所， 北京 100190； 2. 中国科学院大学， 北京 100049)

梯度线圈是磁共振系统的核心部件之一，其产生的梯度磁场是成像空间定位的物理基础。在开放式磁共振系统中，大量铁磁介质会导致实际梯度磁场偏离设计值，因此需要对实际梯度磁场进行校验。本文针对开放式磁共振的结构特点，提出了一种基于扇环形状单元的积分方程法，可以计算二维/三维模型中的梯度磁场。针对两种常见的开放式磁共振结构，通过将本文方法的计算结果与其他三种方法(直接法、镜像法和有限元法)进行对比，证明了本文方法在保证计算精度的同时，可以大幅减少计算时间，在重复计算和迭代优化中具有显著优势。实验结果证明了该算法具有良好的可靠性。

开放式磁共振； 铁磁介质； 梯度线圈； 积分方程法

1 引言

作为一种无创诊断工具，磁共振成像(Magnetic Resonance Imaging，MRI)在医学和科研领域均起着重要作用[1-4]。开放式磁共振是常见的磁共振类型之一。和传统螺线管式磁共振相比，开放式磁共振能够提供一个开放的诊断空间，可以有效缓解患者的紧张情绪，并且给临床的介入治疗提供一个良好的平台。

在MRI系统中，梯度线圈在成像区域产生一个用于图像空间位置编码的梯度磁场，其精度决定了磁共振信号向图像转换时的空间定位精度。因此，扭曲的梯度磁场会导致最终成像出现伪影、扭曲和形变等问题，直接导致图像采集的失败。目前开放式磁共振系统中较为常用的梯度线圈设计方法为连续电流密度空间的方法[5-10]，在梯度线圈设计阶段通常不考虑铁磁介质影响。然而实际开放式磁共振系统中，常采用大量铁磁介质作为导磁回路和磁场调节部件，而铁磁介质受外界磁场磁化，从而在成像区域产生额外的磁场，进而对梯度磁场造成影响。

为了估算铁磁介质给梯度磁场造成的误差，目前较为常用的方法为镜像法[11,12]和有限元法[13,14]。镜像法的假设条件为铁磁介质是一个磁导率无穷大的无限大理想平面，将铁磁介质对空间磁场的影响等效为一个镜像电流，最终空间磁场由源电流和镜像电流叠加计算。镜像法的优点在于计算简单、快速，但缺点是铁磁介质几何结构及材料属性较为复杂时，误差较大。有限元法是目前最常用的电磁数值计算方法，可以计算各种复杂的电磁问题。由于开放式磁共振存在较大的空气域，因而可视为三维开域问题，大量计算资源和时间会被浪费在空气域的网格划分和计算。同时，由于网格在线圈结构变化后需要重新划分和计算，因而计算效率较低。

本文根据开放磁共振系统中铁磁介质具有轴对称的几何特点，提出了一种基于扇环形状单元的积分方程法，并推导了圆柱坐标系中铁磁介质扇环单元在成像区域上采样点的磁场系数矩阵计算公式。基于该方法，对两种常见开放式磁场共振结构下的z方向梯度线圈和x方向梯度线圈进行了仿真计算，并与多种传统方法在计算精度和计算效率两方面进行了对比。最后为了进一步验证本文方法的准确性，利用实际开放式磁体和实际梯度线圈搭建实验平台，并对两种梯度线圈的梯度磁场进行了测量。仿真结果和实验数据均表明本文方法在保证计算精度的前提下大大提高了计算效率，有助于梯度线圈的校验和后续优化工作。

2 理论模型

2.1 扇环单元数学模型

根据开放式磁共振轴对称的结构特点，引入了扇环单元对铁磁单元进行划分。如图1所示的圆柱坐标系中，扇环单元的空间位置由其内外径r1、r2，上下平面坐标z1、z2和圆心角φ1、φ2决定。

假设单元内存在均匀分布的三维磁化强度M，其在单元外空间中任一点P(r,φ,z)处产生三维空间磁场H，则可建立H与M的关系式为：

(1)

式中，A为磁化强度与磁场强度之间的系数；下标r、φ、z分别表示径向分量、周向分量和轴向分量。

图1 柱坐标下扇环单元和空间场点示意图Fig.1 Sector-ring unit and field point in cylindrical coordinates

为了得到系数A，用等效磁流法[15]将磁化强度等效为磁化电流，再用毕奥-萨法尔定律进行计算：

(2)

式中，J为磁化电流密度。以计算系数Arr，Arφ，Arz为例：

利用相同方法，可求得其他系数，联立可得系数矩阵A：

系数矩阵A由铁磁单元的几何尺寸、空间位置以及目标场点的空间位置决定[16]。在进行梯度线圈优化设计和结果校验时，铁磁介质和目标场点是固定不变的，因此在此过程中，该系数矩阵也是固定不变的。因而在相同结构磁体中，涉及梯度磁场计算时，铁磁介质的网格划分和系数矩阵只需进行一次计算，便可进行后续所有的优化和校验。

2.2 积分方程法

空间磁场可以视为由源区(电流源、永磁体等)和磁化区域综合产生[17]：

H=HS+HM

(6)

式中，H为空间中任一点磁场强度；HS为源区产生的磁场；HM为磁化区域产生的磁场。假设铁磁介质共划分L个单元，则式(6)可扩充为：

(7)

式中，HabM表示第b个铁磁介质单元在第a个铁磁单元内产生的磁场强度。将式(1)代入式(6),并整理可得：

(8)

式中，A为上文所求的系数矩阵；χ为该单元铁磁介质磁化率，由铁磁介质磁化曲线和工作点决定。

通过求解式(8)，可求得所有铁磁单元内部的磁场强度，进而可以求得铁磁单元内部的磁化强度。最后结合源区产生的磁场，便可求得空间任意一点的磁场强度。

当假设铁磁介质工作区域接近线性时，磁化率可视为常数，式(8)为线性方程组；当铁磁介质工作区域位于非线性段时，需利用磁化曲线通过迭代计算单元内磁场强度。

3 仿真计算

本节中，采用第2节的方法对梯度线圈产生的梯度磁场进行仿真计算，并与直接法(不考虑铁磁介质)、镜像法和有限元法结果进行比较，从而验证本文方法的计算精度和计算效率。

将现有开放磁共振磁体的结构归纳为两种常见结构(如图2所示)：结构1的特点在于调节成像区域磁场分布的铁磁介质部件位于整体结构的边缘(通常称为匀场环)，通过调节该部件的长度和宽度来修正磁场，是目前最常用的结构；结构2的特点在于调节成像区域磁场分布的铁磁介质位于整体结构的中央，通过调节中央无铁磁介质区域的大小来修正磁场。这两种结构不仅涵盖了常用开放式磁共振的结构特点，其对梯度磁场的影响也有较大差别。因此对两种结构下的梯度磁场进行仿真计算，可以有效地验证计算方法的普适性和计算精度。

图2 两种开放式磁共振结构示意图Fig.2 Two structures of biplanar superconducting MRI

由于在磁共振成像中，梯度磁场仅轴向(z方向)分量参与空间定位，因此下文仅计算梯度线圈的轴向磁场分量，但铁磁介质中的磁化强度是三维空间向量。

3.1 计算精度

3.1.1z方向梯度线圈

z方向梯度线圈的结构如图3所示，其由多个不同匝数和电流方向的同轴电流环组成。

图3 z方向梯度线圈结构图Fig.3 Structure of gradient coil z

现采用四种方法计算z方向梯度线圈在两种结构中的梯度磁场，并进行误差分析，所得误差结果见表1。

表1 z方向梯度线圈磁场计算结果比较Tab.1 Results of gradient coil z

直接法的结果说明在两种铁磁介质结构中，直接法的计算结果均存在较大的误差。镜像法的结果在两种结构中存在较大差别。分析认为这是由于铁磁介质和线圈的结构特点造成的。图3所示的z方向梯度线圈的外侧电流较少，因此该线圈对结构1中的匀场环不敏感，因而结构1整体与镜像法假设的无限大平面接近，所以计算误差较小；而结构2中铁磁介质的中心空孔会极大地破坏镜像法的假设，因此误差较大。本文方法在两种结构中的结果均较为理想，说明本文方法不会由于铁磁介质的结构差异导致计算结果出现明显偏差。相对而言结构1的结果略差于结构2的结果，其原因可能为结构1边缘凸起处的磁化强度变化较大，从而导致磁化强度的计算出现一定的误差。

3.1.2x方向梯度线圈

x方向梯度线圈的结构如图4所示，其由不规则的“指纹”型电流回路组成。

图4 x方向梯度线圈结构图Fig.4 Structure of gradient coil x

现采用四种方法计算该结构x方向梯度线圈在两种结构中的梯度磁场，并进行误差分析，误差结果见表2。

表2 x方向梯度线圈磁场计算结果比较Tab.2 Results of gradient coil x

由表2可知，直接法的结果在两种铁磁介质结构中同样不理想。x方向梯度线圈镜像法的结果与z方向梯度线圈镜像法的结果相反，结构1的结果要差于结构2的结果。其原因为图4所示的x方向梯度线圈，其电流主要集中在边缘区域，因此对于结构1中匀场环的部分比较敏感，而对结构2中铁磁介质中心空孔不太敏感，从而导致结构1的结果要差于结构2的结果。分析本文方法的结果可以发现其在两种结构中依旧能保证精度。

结合z方向和x方向两种梯度线圈的分析结果可知：①直接法的结果通常都不理想，因此用常规设计方法得到的梯度线圈无法直接在实际磁体中使用，需要进行校验和优化；②镜像法的结果往往受到线圈和铁磁介质结构的影响，计算精度无法得到保证；③本文方法不受线圈和铁磁介质结构影响，计算精度可以得到有效控制。

3.2 计算效率

在同一台计算机上用Matlab运行本文所用方法和Ansoft有限元软件计算同一组案例。计算机配置如下：4核CPU i5-3740 3.2GHz；内存8GB。

采用有限元法时，对于z方向梯度线圈，由于其具有轴对称结构，因此采用二维建模，计算速度很快，仅需2min；对于x方向梯度线圈，三维模型的计算时间较长，约需12min，由于线圈形状特殊，因此需要联合Solidworks等软件进行建模，流程比较复杂。

采用本文方法，两种梯度线圈均采用三维模型，总用时小于1min。其中系数矩阵计算时间为50.7s，实际用于磁场计算小于5s。

通过以上数据对比可以发现，本文方法即使考虑系数矩阵计算，在计算效率上相对于有限元法依然有明显优势；如果再考虑系数矩阵无需重复计算的情况，则计算速度方面的优势将更加显著。

4 实验结果

为了验证本文理论方法和仿真结果的正确性，以开放式超导磁共振实际磁体为实验平台，通过对安装在实验平台中的梯度线圈产生的梯度磁场进行测量，并与仿真结果进行对比，从而对本文方法进行验证。

实验平台和测量设备如图5所示。图5(a)为0.7T开放式超导磁共振的主磁体；图5(b)为测量工装，测量工装安装在上下极面之间，测量工装上距离中心不同距离的测量点为成像区域上的磁场采样点；图5(c)和图5(d)分别为z方向梯度线圈和x方向梯度线圈。

图5 实验平台和测量设备Fig.5 Experimental platform and measuring equipment

成像区域梯度磁场的测量流程为：①将主磁体励磁至工作电流，此时中心点磁场为0.7 T；利用测量工装测试成像区域测量点处的磁场强度(测量仪器为Metrolab PT2025核磁共振仪)，该磁场强度为主磁场强度；②将梯度线圈外接一个恒定电流源，并施加1A电流，再次测量成像区域测量点处的磁场强度，该磁场强度为主磁场强度和梯度磁场强度之和；③将步骤②中的磁场强度减去步骤①中的磁场强度，得到在1A工作电流下，梯度线圈产生的实际梯度磁场强度。通过上述流程，分别对z方向梯度线圈和x方向梯度线圈产生的梯度磁场进行测量，并将实验测量结果与四种方法的计算结果进行对比，结果见表3。

表3 仿真结果与实验数据比较Tab.3 Comparison between experimental result and simulation result

由于开放式磁共振实际磁体采用的是结构1的铁磁介质结构，因此实验结果与结构1的结果相近。直接法和镜像法结果均存在较大误差，有限元法和本文方法相对精度较高，但依旧与实验结果有一定的偏差。造成这种偏差的原因有两方面，一方面是铁磁介质仿真时采用的材料属性和实际材料属性有一定偏差，同时实际铁磁介质材料本身存在差异；另一方面梯度线圈在实际制作过程中无法达到理想模型的精度，并且还需要额外添加进出线。

5 结论

本文基于开放式磁共振的轴对称结构特点，提出了一种基于扇环形状单元的积分方程法计算开放式磁共振的梯度磁场。利用该方法计算开放式磁共振的梯度磁场相比现有方法具有以下优点：

(1)计算精度高且不受磁体结构影响。当模型中铁磁介质具有轴对称特点或可以采用扇环单元进行近似划分时，均可采用本文方法计算。在两种开放式磁共振结构中，本文方法的计算精度已证明具有实际应用价值。

(2)计算效率高。本文方法只对铁磁介质划分网格，相对于有限元法中存在大量的空气域网格，可以节省大量计算资源，有效提高计算效率。

(3)有利于重复校验和后续优化算法。在铁磁介质结构不变时，本文方法只需计算一次系数矩阵。在重复验算和优化的迭代计算中，具有明显的优势。

仿真结果和实验测量结果都很好地验证了本文方法在计算精度和计算效率上的优势，因此本文方法有利于开放式磁共振梯度线圈的实际验算，并对后续的优化计算具有一定的参考意义。

[1] Qiuliang Wang. Practical design of magnetostatic structure using numerical simulation[M]. Singapore: John Wiley & Sons, 2013.

[2] Q Wang, Y Dai, B Zhao, et al. A superconducting magnet system for whole-body metabolism imaging[J]. IEEE Transactions on Applied Superconductivity, 2012, 22 (3): 4400905.

[3] Z Wang, T Song, Z R Dong, et al. A method on decreasing magnetic moment of Halbach cylinder magnets[J]. IEEE Transactions on Applied Superconductivity, 2008, 18 (2): 961-964.

[4] 籍勇亮, 何为, 吴高林, 等(Ji Yongliang, He Wei, Wu Gaolin, et al.). 单边核磁共振仪中射频线圈的优化(Optimal of RF coil in single-side NMR system) [J]. 电工电能新技术( Advanced Technology of Electrical Engineering and Energy ), 2015, 34 (9):75-80.

[5] Q M Tieng, V Vegh, I M Brereton. Globally optimal superconducting magnets - Part I: Minimum stored energy (MSE) current density map [J]. Journal of Magnetic Resonance, 2009, 196 (1): 1-6.

[6] Geli Hu, Junsheng Chen, Zhiping Ni, et al. A novel target field approach to design of biplanar gradient coils for permanent MRI system [A]. 2013 IEEE International Conference on Applied Superconductivity and Electromagnetic Devices (ASEMD) [C]. 2013. 446-449.

[7] Larry K Forbes, Stuart Crozier. Novel target-field method for designing shielded biplanar shim and gradient coils [J]. IEEE Transactions on Magnetics, 2004, 40 (4): 1929-1938.

[8] F Romeo, D I Hoult. Magnet field profiling: Analysis and correcting coil design [J]. Magnetic Resonance in Medicine, 1984, 1 (1): 44-65.

[9] R Turner．A target field approach to optimal coil design [J]．IEEE Transactions on Applied Physics，1986，19 (8)： 147-151．

[10] R Turner．Gradient coil design：A review of methods [J]．IEEE Transactions on Magnetic Resonance Imaging，1993, 11 (7)：903-920．

[11] L K Forbes，S Crozier．A novel target-field method for finite-length magnetic resonance shim coils - Part I: Zonal shims[J]．Journal of Physics D: Applied Physics，2001，34 (24)：3447-3455．

[12] 俎栋林, 洪礼明, 曹学明, 等 (Zu Donglin, Hong Liming, Cao Xueming, et al.). 对永磁 MRI 脉冲梯度产生涡流的背景磁场的解析(Analysis on background magnetic field to generate eddy current by pulsed gradient of permanent-magnet MRI)[J]. 中国科学: 技术科学(Science China Technological Sciences), 2010,40 (10): 1221-1226.

[13] 刘若茜, 黄常纲(Liu Ruoqian, Huang Changgang). 永磁磁共振成像系统中磁体极板尺寸与梯度线圈的关系(The effects of gradient coil on pole size of magnet for permanent MRI system)[J]. 电工电能新技术(Advanced Technology of Electrical Engineering and Energy), 2000, 19 (1): 14-18.

[14] 郑建立, 张志伟, 聂生东,等(Zheng Jianli，Zhang Zhiwei，Nie Shengdong, et al.). 基于 Ansoft 的台式磁共振成像仪双平面梯度线圈的有限元分析(Finite element analysis of the bi-planar gradient coils of the desktop-magnetic resonance imaging based on Ansoft)[J]. 中国医学物理学杂志(Chinese Journal of Medical Physics), 2012, 29 (1): 3133-3138.

[15] 宁毅(Ning Yi). 新型永久磁体及双平面纵向屏蔽梯度线圈设计(Designs of novel permanent magnet and bi-planar shielded longitudinal gradient coils)[D]. 武汉：华中科技大学(Wuhan: Huazhong University of Science and Technology), 2005.

[16] E P Furlani. Permanent magnet and electromechanical devices: Materials, analysis, and applications [M]. New York: Academic Press, 2001. 209-264.

[17] J Caldwell, A Zisserman, R Saunders. A GFUN approach to include the effects of iron on coil systems with cylindrical symmetry[J]. Journal of Physics D: Applied Physics, 1984, 17 (9): 1759-1772.

A fast calculation method of gradient field in biplanar MRI

ZHU Xu-chen1,2, WANG Qiu-liang1, WANG Hou-sheng1

(1. Institute of Electrical Engineering， Chinese Academy of Sciences， Beijing 100190， China; 2. University of Chinese Academy of Sciences， Beijing 100049， China)

Gradient coil is one of the key components in the magnetic resonance imaging (MRI) system. The effect of gradient magnetic field is to provide the spatial orientation of imaging. In biplanar MRI system, a large number of ferromagnetic materials would destroy the gradient magnetic field, so it is necessary to calibrate the gradient magnetic field. Aiming at the structural characteristics of biplanar MRI, a new method to calculate the gradient magnetic field in three dimensional model based on the integral equation method and sector-ring element is presented in this paper. In two different MRI structures, simulation results of the proposed method is compared with the method without ferromagnetic materials, mirror method and finite element method. It is proved that the proposed method can effectively reduce the computing time and has a significant advantage in the repeated calculation and iterative optimization while ensuring the accuracy. Experimental results also show that the proposed method has good reliability.

biplanar MRI; ferromagnetic material; gradient coil; integral-equation method

2016-01-20

国家自然科学基金项目(50925726； 50577063)

朱旭晨(1989-)， 男， 上海籍， 博士研究生， 研究方向为磁共振成像系统电磁设计与优化; 王秋良(1965-)， 男， 湖北籍， 研究员， 博士， 研究方向为超导磁体与强磁场应用。

TM153

A

1003-3076(2016)10-0043-06