CN113325350B - 一种基于离散网格的高性能梯度线圈设计方法 - Google Patents

一种基于离散网格的高性能梯度线圈设计方法 Download PDF

Info

Publication number
CN113325350B
CN113325350B CN202010129685.5A CN202010129685A CN113325350B CN 113325350 B CN113325350 B CN 113325350B CN 202010129685 A CN202010129685 A CN 202010129685A CN 113325350 B CN113325350 B CN 113325350B
Authority
CN
China
Prior art keywords
matrix
gradient coil
vector
current density
design method
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010129685.5A
Other languages
English (en)
Other versions
CN113325350A (zh
Inventor
平学伟
殷兴辉
刘海韵
陈嘉琪
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hohai University HHU
Original Assignee
Hohai University HHU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hohai University HHU filed Critical Hohai University HHU
Priority to CN202010129685.5A priority Critical patent/CN113325350B/zh
Publication of CN113325350A publication Critical patent/CN113325350A/zh
Application granted granted Critical
Publication of CN113325350B publication Critical patent/CN113325350B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/28—Details of apparatus provided for in groups G01R33/44 - G01R33/64
    • G01R33/38—Systems for generation, homogenisation or stabilisation of the main or gradient magnetic field
    • G01R33/385—Systems for generation, homogenisation or stabilisation of the main or gradient magnetic field using gradient magnetic field coils

Abstract

本发明公开了一种基于离散网格的高性能梯度线圈设计方法。该方法将成像区域内的磁场误差平方和与功耗的加权和作为待优化的目标函数，首先将梯度线圈所在的面进行网格剖分，然后在每个网格内将电流密度用基函数展开，并将展开式带入目标函数，可以得到一个稠密的线性系统。本发明采用迭代法求解该线性系统，并提出了一种加速迭代过程中的矩阵矢量乘的方法，从而能够极大的提升优化方法的效率。

Description

一种基于离散网格的高性能梯度线圈设计方法
技术领域
本发明涉及一种高性能梯度线圈设计方法，更具体地涉及梯度线圈优化方法中的高效线性系统求解技术。
背景技术
梯度线圈是磁共振设备中的核心部件，决定了磁共振成像的很多重要指标。从上世纪80年代至今，如何设计高性能梯度线圈一直是众多学者研究的一个重要问题。目前最常用的梯度线圈设计方法是流函数方法。在流函数法中，基于网格离散的方法获得了很多学者的关注与研究。该类方法具有如下一些优点。首先，该类方法可以适用于复杂的线圈形状；其次，对于不同形状的梯度线圈，只要采样同种类型的剖分网格，就可以采用同一种基函数，有利于编写梯度线圈设计通用代码。其次，该类方法与其他方法相比具有更多的自由度，从而能够保证所设计的线圈具有最佳的电磁性能。但是该类方法在实际应用中往往生成一个具有数千甚至数万未知量的稠密线性系统。如果采用直接解法求解该线性系统，则一次求解需要的时间会很长。在进行梯度线圈设计时，往往需要反复修改输入参数来确定最合适的线圈参数。因此，在进行梯度线圈设计时需要很多次求解线性系统的操作。线性系统求解方法的性能会严重影响优化方法的效率。因此，有必要研究快速求解该线性系统的方法。
发明内容
本发明目的在于采用离散网格法设计梯度线圈，并给出一种求解离散网格法产生的线性系统的快速迭代求解方法与加速矩阵矢量乘的方法，以提升线圈设计的效率。
为实现上述发明目的，本发明所述的一种基于离散网格的高性能梯度线圈设计方法，包括如下步骤：
(1)将梯度线圈布线区域进行网格剖分并对网格进行编号；
(2)在每个网格内定义电流密度基函数，并将每个电流密度基函数分别用局部编号与全局编号表示，在每个网格内，电流密度用基函数展开式表示为：
Figure BDA0002395453680000021
式中m为网格编号，Nm为第m个网格内的基函数个数，
Figure BDA0002395453680000022
为第m个网格内任意一点的坐标矢量，
Figure BDA0002395453680000023
为
Figure BDA0002395453680000024
处的电流密度，
Figure BDA0002395453680000025
为待求系数,
Figure BDA0002395453680000026
为电流密度基函数。
(3)构建如下的目标函数：
Figure BDA0002395453680000027
式中α为权重系数，S为梯度线圈所在的面，P为采样点个数，
Figure BDA0002395453680000028
为编号为p的采样点处的坐标矢量，
Figure BDA0002395453680000029
为在第p个坐标点处梯度线圈产生的z方向磁场与理想的z方向磁场；
(4)将电流密度展开式带入目标函数并令
Figure BDA00023954536800000210
得到以下的线性系统：
Aa＝b
其中a是由待求系数构成的矢量，A为一个矩阵，其表达式为：
Figure BDA00023954536800000211
上式中，
Figure BDA00023954536800000212
为一个列矢量，其第i个元素表达式为：
Figure BDA00023954536800000213
上式中，μ0为真空磁导率，
Figure BDA00023954536800000214
为全局编号为i的电流密度基函数；
矩阵B第i行j列元素表达式为：
Figure BDA00023954536800000215
矢量b的表达式为：
Figure BDA0002395453680000031
(5)采用迭代法求解(4)中的线性系统，并采用如下方法计算迭代方法中的矩阵矢量乘运算：
Figure BDA0002395453680000032
p＝1,2…P
Q0＝BQ，
Figure BDA0002395453680000033
上式中，Q为需要与矩阵A进行矩阵矢量乘运算的矢量，R为矩阵矢量乘运算得到的结果；
(6)根据步骤(5)得到的解求得梯度线圈的结构。
进一步，采用三角形网格进行网格剖分。
进一步，在每个三角形内，采用如下的电流密度基函数；
Figure BDA0002395453680000034
上式中，m为三角形网格编号，
Figure BDA0002395453680000035
为第m个三角形内的
Figure BDA0002395453680000036
点对应的第i个面积坐标，
Figure BDA0002395453680000037
为三角形的单位法向矢量，▽×表示取旋度。
进一步，采用Krylov子空间迭代法求解步骤(4)中的线性系统。
进一步，采用CSR稀疏存储格式存储B矩阵。
进一步，采用对角预条件加速迭代法的收敛，并采用如下方法计算对角预条件矩阵G的元素：
Figure BDA0002395453680000038
进一步，采用如下方法施加预条件矩阵：
(1)计算矩阵F：F＝G-1/2；
(2)更新矢量c：
Figure BDA0002395453680000039
(3)更新矩阵B：B＝FBF；
(4)更新右侧矢量b：b＝Fb；
(5)用新的c、B、b进行迭代过程中的矩阵矢量乘运算；
(6)迭代结束后，对解a进行如下运算：a＝Fa。
与现有技术相比，本发明设计方法的优点是：采用本发明中的控制方程设计的梯度线圈光滑性好，易于加工；与直接求解方法相比，本发明中的线性系统求解方法求解速度快了上千倍，能节省大量的设计时间。
附图说明
图1为柱面梯度线圈布线区域内的网格剖分示意图。
图2为采用离散网格法设计的柱面梯度线圈示意图。
具体实施方式
本发明是关于一种梯度线圈设计方法。发明内容是提出了一种基于离散网格的梯度线圈设计方案。在该方案中，首先将布线区域进行网格剖分，并将每个网格内的电流密度表示为基函数与待求系数的组合，然后定义待优化的目标函数，令目标函数对待求系数的导数为0，得到一个稠密的线性系统。采用Krylov子空间迭代法，并采用快速的矩阵矢量乘求解技术，可以快速求解该线性系统。下面我们针对磁共振系统中的梯度线圈设计来说明该方法的具体实施方案。
在磁共振系统中，梯度线圈分布在平面、柱面或椭圆柱面骨架上。本发明实施例公开的一种基于离散网格的高性能梯度线圈设计方法如下：
(1)将梯度线圈布线区域进行网格剖分并对每个网格进行编号；
(2)在每个网格内定义电流密度基函数，并将每个电流密度基函数分别用局部编号与全局编号表示，在每个网格内，电流密度用基函数展开式表示为：
Figure BDA0002395453680000041
式中m为网格编号，
Figure BDA0002395453680000042
为第m个网格内任意一点处的坐标矢量，
Figure BDA0002395453680000043
为
Figure BDA0002395453680000044
处的电流密度，
Figure BDA0002395453680000045
为待求系数，
Figure BDA0002395453680000046
为第m个三角形内的
Figure BDA0002395453680000047
点对应的第i个面积坐标(关于面积坐标的定义可以参阅金建铭的《电磁场有限元方法》)，
Figure BDA0002395453680000048
为三角形的单位法向矢量，
Figure BDA0002395453680000049
表示取旋度。
(3)构建如下的目标函数：
Figure BDA00023954536800000410
式中α为权重系数，S为梯度线圈所在的面，P为采样点个数，
Figure BDA0002395453680000051
为编号为p的采样点处的坐标矢量；
(4)将电流密度展开式带入目标函数并令
Figure BDA0002395453680000052
得到以下的线性系统：
Aa＝b
其中a是由所有待求系数构成的矢量，表达式为：
a＝[ai],i＝1,2,...N
ai为全局编号为i的待求系数，N为待求系数的个数，这里等于网格内部节点的数目。A为一个矩阵，其表达式为：
Figure BDA0002395453680000053
上式中，
Figure BDA0002395453680000054
为一个列矢量，其第i个元素表达式为：
Figure BDA0002395453680000055
上式中，μ0为真空磁导率，
Figure BDA0002395453680000056
为布线区域内任意一点处的坐标矢量，
Figure BDA0002395453680000057
为全局编号为i的电流密度基函数，S为线圈所在的面；
矩阵B第i行j列元素表达式为：
Figure BDA0002395453680000058
矢量b的表达式为：
Figure BDA0002395453680000059
(5)采用Krylov子空间迭代法求解(4)中的线性系统。Krylov子空间迭代法指的是一类算法，常用的共轭梯度法(CG)、双共轭梯度法(BCG)、一般最小余量法(GMRES)都属于这类算法。这类算法的主要运算量是矩阵矢量乘运算，即矩阵A与某一矢量Q相乘的操作。只考虑乘法运算，该操作的运算量为N2。本发明根据矩阵的特点，采用如下方法计算迭代法中的矩阵矢量乘运算：
Figure BDA0002395453680000061
p＝1,2…P
Q0＝BQ，
Figure BDA0002395453680000062
上式中，R为矩阵矢量乘运算得到的结果，βp、Qp、Q0为中间变量。B为一个维度为N×N的矩阵，其大部分元素为0。因此，采用CSR稀疏存储格式存储B矩阵可以极大的节省存储空间与计算量。所谓CSR存储格式，在一些文献里有介绍(如Yousef Saad的书籍IterativeMethods for Sparse Linear Systems)，这里不再赘述。
在应用Krylov子空间迭代法时，一般采用预条件技术来加速迭代法的收敛。本实施例中采用对角预条件，采用如下方法计算对角预条件矩阵G的元素：
Figure BDA0002395453680000063
并采用如下方法施加预条件矩阵：
(1)计算矩阵F：F＝G-1/2；
(2)更新矢量c：
Figure BDA0002395453680000064
(3)更新矩阵B：B＝FBF；
(4)更新右侧矢量b：b＝Fb；
(5)用新的c、B、b进行迭代过程中的矩阵矢量乘运算；
(6)迭代结束后，对解a进行如下运算：a＝Fa。
由于G为对角阵，B为稀疏矩阵，上述操作所需运算量非常少，且只需在迭代前对矩阵进行一次处理，并在迭代结束后对解进行一次处理即可。
在求解上述的线性系统后，即可得到线圈骨架上的电流密度分布、根据电流密度分布可求出流函数分布，进而求得梯度线圈的结构。根据流函数求线圈结构的方法在很多关于梯度线圈设计的文献中有叙述，这里不再赘述。
下面根据本具体实施例的步骤设计一个结构对称的柱面梯度线圈。柱面半径为r＝0.3m，轴向长度为L＝1.2m。梯度场强约为G＝60.5uT/m/A，在直径为45cm的球形区域内，线性度为6％。采用离散网格法设计本梯度线圈，首先将布线区域采用三角形网格剖分。考虑线圈结构的对称性，只需对第1、2象限内的柱面进行网格剖分，如图1所示。网格内不在边界上的点为4395个。也就是说，采用离散网格法对前述目标函数进行离散，生成的线性系统的系数矩阵为一个维数为4395×4395的稠密矩阵。采用直接解法求解该线性系统，每次求解耗时约为1570秒。采用GMRES迭代法，如果不采用本发明中的矩阵矢量乘方法，则每次求解线性系统的时间约为7.9秒。采用本发明中的矩阵矢量乘方法，则每次求解时间约为0.9秒。可以看出，本发明中的线性系统求解方法大大节省了线性系统求解的时间，从而节省了线圈设计时间。
求解线性系统后，即可得到每个网格节点处的流函数值。根据该值可以求出梯度线圈的结构，如图2所示。
以上所述仅为本发明的实施例子，并不用于限制本发明，本发明中的方法对于各种形状的梯度线圈均同样适用。凡在本发明的原则之内，所作的等同替换，均应包含在本发明的保护范围之内。本发明未作详细阐述的内容属于本专业领域技术人员公知的已有技术。

Claims (7)

1.一种基于离散网格的高性能梯度线圈设计方法，其特征在于，包含以下几个步骤：
(1)将梯度线圈布线区域进行网格剖分并对网格进行编号；
(2)在每个网格内定义电流密度基函数，并将每个电流密度基函数分别用局部编号与全局编号表示，在每个网格内，电流密度用基函数展开式表示为：
Figure FDA0002395453670000011
式中m为网格编号，Nm为第m个网格内的基函数个数，
Figure FDA0002395453670000012
为第m个网格内任意一点的坐标矢量，
Figure FDA0002395453670000013
为
Figure FDA0002395453670000014
处的电流密度，
Figure FDA0002395453670000015
为待求系数,
Figure FDA0002395453670000016
为电流密度基函数；
(3)构建如下的目标函数：
Figure FDA0002395453670000017
式中α为权重系数，S为梯度线圈所在的面，P为采样点个数，
Figure FDA0002395453670000018
为编号为p的采样点处的坐标矢量，
Figure FDA0002395453670000019
为在第p个坐标点处梯度线圈产生的z方向磁场与理想的z方向磁场；
(4)将电流密度展开式带入目标函数并令
Figure FDA00023954536700000110
得到以下的线性系统：
Aa＝b
其中a是由待求系数构成的矢量，A为一个矩阵，其表达式为：
Figure FDA00023954536700000111
上式中，
Figure FDA00023954536700000112
为一个列矢量，其第i个元素表达式为：
Figure FDA00023954536700000113
上式中，μ0为真空磁导率，
Figure FDA0002395453670000021
为全局编号为i的电流密度基函数；
矩阵B第i行j列元素表达式为：
Figure FDA0002395453670000022
矢量b的表达式为：
Figure FDA0002395453670000023
(5)采用迭代法求解(4)中的线性系统，并采用如下方法计算迭代法中的矩阵矢量乘运算：
Figure FDA0002395453670000024
Q0＝BQ，
Figure FDA0002395453670000025
上式中，Q为需要与矩阵A进行矩阵矢量乘运算的矢量，R为矩阵矢量乘运算得到的结果；
(6)根据步骤(5)得到的解求得梯度线圈的结构。
2.根据权利要求1所述的梯度线圈设计方法，其特征在于，采用三角形网格进行网格剖分。
3.根据权利要求2所述的梯度线圈设计方法，其特征在于，在每个三角形内，采用如下的电流密度基函数；
Figure FDA0002395453670000026
上式中，m为三角形网格编号，
Figure FDA0002395453670000027
为第m个三角形内的
Figure FDA0002395453670000028
点对应的第i个面积坐标，
Figure FDA0002395453670000029
为三角形的单位法向矢量，
Figure FDA00023954536700000210
表示取旋度。
4.根据权利要求1所述的梯度线圈设计方法，其特征在于，采用Krylov子空间迭代法求解步骤(4)中的线性系统。
5.根据权利要求1所述的梯度线圈设计方法，其特征在于，采用CSR稀疏存储格式存储B矩阵。
6.根据权利要求1所述的梯度线圈设计方法，其特征在于，采用对角预条件加速迭代法的收敛，并采用如下方法计算对角预条件矩阵G的元素：
Figure FDA0002395453670000031
7.根据权利要求6所述的梯度线圈设计方法，其特征在于，采用如下方法施加预条件矩阵：
(1)计算矩阵F：F＝G-1/2；
(2)更新矢量c：
Figure FDA0002395453670000032
(3)更新矩阵B：B＝FBF；
(4)更新右侧矢量b：b＝Fb；
(5)用新的c、B、b进行迭代过程中的矩阵矢量乘运算；
(6)迭代结束后，对解a进行如下运算：a＝Fa。
CN202010129685.5A 2020-02-28 2020-02-28 一种基于离散网格的高性能梯度线圈设计方法 Active CN113325350B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010129685.5A CN113325350B (zh) 2020-02-28 2020-02-28 一种基于离散网格的高性能梯度线圈设计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010129685.5A CN113325350B (zh) 2020-02-28 2020-02-28 一种基于离散网格的高性能梯度线圈设计方法

Publications (2)

Publication Number Publication Date
CN113325350A CN113325350A (zh) 2021-08-31
CN113325350B true CN113325350B (zh) 2022-04-08

Family

ID=77412700

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010129685.5A Active CN113325350B (zh) 2020-02-28 2020-02-28 一种基于离散网格的高性能梯度线圈设计方法

Country Status (1)

Country Link
CN (1) CN113325350B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113805485B (zh) * 2021-10-13 2023-04-14 吉林大学 一种暖启动c/gmres方法、系统、设备及介质

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103018689A (zh) * 2012-08-27 2013-04-03 嘉恒医疗科技有限公司 基于流函数的磁共振射频线圈设计方法

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4456881A (en) * 1982-01-18 1984-06-26 Technicare Corporation Gradient-coil apparatus for a magnetic resonance system
DE4203582C2 (de) * 1992-02-07 1994-03-03 Siemens Ag Transversale Gradientenspule
CN101996273B (zh) * 2010-11-29 2012-07-04 浙江大学 Mri系统梯度线圈的有限差分设计方法
CN104007406B (zh) * 2014-06-13 2016-08-24 河海大学 基于空间变换的椭圆柱面横向梯度线圈设计方法
CN104007407B (zh) * 2014-06-13 2016-08-24 河海大学 结合泛函与模拟退火算法的梯度线圈设计方法
WO2016183284A1 (en) * 2015-05-12 2016-11-17 Hyperfine Research, Inc. Radio frequency coil methods and apparatus
CN105205227B (zh) * 2015-08-31 2018-05-08 河海大学 一种基于最小vrms误差的梯度线圈设计方法
CN105445683B (zh) * 2015-11-16 2018-05-11 河海大学 一种柱面横向自屏蔽梯度线圈设计方法
CN105718677A (zh) * 2016-01-22 2016-06-29 中国科学院电工研究所 自屏蔽超导核磁共振成像系统梯度线圈设计方法
CN107831461A (zh) * 2017-10-24 2018-03-23 中国计量大学 基于0‑1整数规划的纵向梯度线圈设计方法
CN108802646B (zh) * 2017-10-24 2021-11-02 中国计量大学 基于0-1整数规划的永磁型纵向梯度线圈
CN109885957A (zh) * 2019-03-01 2019-06-14 中国科学院长春光学精密机械与物理研究所 一种用于生物医学成像系统的梯度线圈设计方法及梯度线圈
CN110568390A (zh) * 2019-10-14 2019-12-13 惠仁望都医疗设备科技有限公司 一种双平面匀场线圈的设计方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103018689A (zh) * 2012-08-27 2013-04-03 嘉恒医疗科技有限公司 基于流函数的磁共振射频线圈设计方法

Also Published As

Publication number Publication date
CN113325350A (zh) 2021-08-31

Similar Documents

Publication Publication Date Title
CN110110413B (zh) 一种基于材料场缩减级数展开的结构拓扑优化方法
CN113325350B (zh) 一种基于离散网格的高性能梯度线圈设计方法
CN111563299A (zh) 一种旋翼噪声确定方法及系统
Wang et al. Application of tree-cotree splitting to the time-domain finite-element analysis of electromagnetic problems
CN109657288B (zh) 一种三维显隐时域电磁学数值方法
Zhang et al. A morphing wing with cellular structure of non-uniform density
CN111881629B (zh) 一种气动热-结构热传导耦合非线性降阶模型方法
WO2022126348A1 (zh) 一种基于各向异性体调和场的附面层网格生成方法
CN116013443B (zh) 一种传热特性预测方法、装置、设备及可读存储介质
CN111027250B (zh) 一种基于网格变形技术的异形曲面加筋壳建模方法
Toal et al. Geometric filtration using POD for aerodynamic design optimization
Zheng et al. An enhanced topology optimization approach based on the combined MMC and NURBS-curve boundaries
CN109948253B (zh) 薄板无网格Galerkin结构模态分析的GPU加速方法
Liu et al. A finite element domain decomposition combined with algebraic multigrid method for large-scale electromagnetic field computation
Namura et al. Multipoint design optimization of vortex generators on transonic swept wings
CN116362079A (zh) 一种基于新型插值模型的多材料结构拓扑优化方法
Chen et al. Efficient multiscale methods for the semiclassical Schrödinger equation with time-dependent potentials
CN108280309A (zh) 一种应用于pic静电模型的电位有限元求解算法
Lin et al. Three-Dimensional-Slice-Super-Resolution-Net: A Fast Few Shooting Learning Model for 3D Super-Resolution Using Slice-Up and Slice-Reconstruction
Krukovskiy et al. Implementing a fully conservative Lagrangian–Euler scheme for two-dimensional problems of magnetic gas dynamics
CN112182909A (zh) 一种用于工业cae方向的流动求解器建立方法
Dziekonski et al. GPU-accelerated finite element method
Yang et al. Domain decomposition combined radial basis function collocation method to solve transient eddy current magnetic problems with moving conductors
Kompiš et al. Computational simulation methods for fiber reinforced composites
CN115130337B (zh) 一种磁悬浮执行器解耦器的设计实现方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant