一种在电磁场棱边单元法中施加规范的方法.pdf

本发明公开了一种在电磁场棱边单元法中施加规范的方法，包括：建立待分析对象三维模型，并对所述三维模型进行网格剖分；计算剖分后的棱边单元的磁矢势双旋度项、库伦规范项以及电流密度项所对应的单元矩阵，并进行集成形成总体有限元方程；对总体有限元方程施加对应的边界条件；对施加对应的边界条件的总体有限元方程进行求解，获得唯一的磁矢势基本解，并对获得的磁矢势基本解进行后处理，得到对应的电磁场量的解。本发明在保证精度的前提下，有效的解决了棱边单元法中的磁矢势唯一性问题；避免了传统方法复杂生成树的构建问题；同时效率高于棱边单元中所采用的拉格朗日乘子法施加库伦规范的方法。

1.  一种在电磁场棱边单元法中施加规范的方法，其特征在于：
包括如下步骤：
⑴、建立待分析对象三维模型，并对所述三维模型进行网格剖分；
⑵、计算剖分后的棱边单元的磁矢势双旋度项所对应的单元矩阵、库伦规范项所对应的单元矩阵以及电流密度项所对应的单元矩阵，并对它们进行集成，形成总体有限元方程；
⑶、对步骤⑵所形成的总体有限元方程，施加对应的边界条件；即在总体有限元方程的强迫边界上迫使剖分后的棱边单元的磁矢势自由度为0；
⑷、对步骤⑶所形成的总体有限元方程进行求解，获得唯一的磁矢势基本解，并对获得的磁矢势基本解进行后处理，得到对应的电磁场量的解。

2.  根据权利要求1所述的一种在电磁场棱边单元法中施加规范的方法，其特征在于：
所述对磁矢势双旋度项的单元矩阵进行计算，包括：
不导电空间中的静磁场方程：
▿×[v]▿×A=Jsin Ω---(1)]]>
A×n=0on ∂Ω---(2)]]>
其中，[ν]:磁阻矩阵(磁导率矩阵的逆矩阵)；A:磁矢势；J_s:电流密度矢量源项；旋度(curl)计算符号；Ω:不导电的三维空间领域；n:分析领域的外法向；强迫边界；为磁矢势的双旋度项；
对于满足边界条件(2)的任意磁矢势A，建立对应方程(1)的泛函：
F(A)=12∫Ω[v]|▿×A|2dV-∫ΩJs·AdV]]>
最小化F(A)，即可得到(1)、(2)所对应的微分方程解；
同时单元内的磁矢势通过形函数表示为：
A＝W·A_e                  (4)
其中，W是单元磁矢势矢量形函数矩阵；A_e是单元棱边上的磁矢势自由度；
根据公式(3)以及(4)，
棱边单元的有限元方程表示为：
KAAAe=Jesin Ω---(5)]]>
其中，
KAA=∫Ve(▿×WT)T[v](▿×WT)dV---(6)]]>
即为所要计算的双旋度项单元矩阵，为单元内离散后的电流密度项。

3.  根据权利要求1所述的一种在电磁场棱边单元法中施加规范的方法，其特征在于：
所述库伦规范项的单元矩阵计算：
鉴于方程(1)、(2)中的磁矢势A具有不唯一性，需要施加额外的约束条件来保证其唯一性，采用库伦规范条件：
▿·A=0in Ω---(7)]]>
并将库伦规范条件(7)作为一个罚项引入到方程(1)、(2)所对应的泛函(3)中，最小化泛函，结合公式(4)，得到相应的有限元方程为：
(KAA+f·GAA)Ae=Jesin Ω---(8)]]>
其中，K^AA表示双旋度项的单元矩阵，f·G^AA为所引入的罚项即库伦规范项所对应的单元矩阵，f为罚因子，需要按照不同情况做适当的选取；
G^AA为磁矢库伦规范所对应的罚项的单元矩阵：其表达式为：
GAA=∫Ve(▿·WT)T[v](▿·WT)dV---(9).]]>

4.  根据权利要求1所述的一种在电磁场棱边单元法中施加规范的方法，其特征在于：
所述电流密度项对应的单元矩阵的计算：
计算公式(1)中的单元电流密度源项，即公式(5)中的
Jes=∫VeWTJSdV---(10)]]>
最后，对上述的单元矩阵K^AA、G^AA、进行集成，形成总体有限元方程。

一种在电磁场棱边单元法中施加规范的方法
技术领域
本发明涉及一种在电磁场棱边单元法中施加规范的方法，属于计算机辅助工程(CAE)软件中新的数值方法类，具体涉及电磁场棱边单元法求解的核心方法，可广泛应用于电气行业，如电机、电磁阀等等的电磁场分析。
背景技术
电磁场有限元方法是一种在电磁场数值求解中广泛应用的数值方法，然而在保证精度的情况下解决磁矢势唯一性问题一直是个难点。通过对已有文献和技术资料的详细查询，总结出目前在电磁场有限元中常常采用的4种技术手段：1采用棱边单元法并结合树-余树方法，2采用棱边单元法并结合CG法(或预处理CG法)，3采用棱边单元法并结合拉格朗日乘子法施加库伦规范的方法，以及4采用节点元法结合罚函数法施加库伦规范的方法；
第一种方法是棱边单元法结合树-余树方法，这种方法是在待分析的网格空间建立起Spanning Tree，如图1所示(粗实线为生成的树),对位于Spanning Tree上所对应的棱边上的磁矢势的自由度进行约束从而确保解的唯一性；采用棱边单元法并结合树-余树方法，虽然避免了采用节点元法的精度问题，但是采用树-余树方法，面临着生成树的方法不唯一，而不合理的生成树会导致不合理的解的问题，另外复杂问题的最优生成树的构建存在困难。
第二种方法是采用棱边单元法并结合CG法(或预处理CG法)，这种方法是对棱边单元法离散所得到的有限元方程，采用CG法(或预处理CG法)进行求解，在不唯一的解空间中获得合理的解；但是采用棱边单元法并结合CG法(或预处理CG法)需要对电流密度场进行重新求解，增大了计算量同时降低了计算精度，并且CG迭代方法常常存在着收敛性问题。
第三种是采用拉格朗日乘子法引入库伦规范条件，这种方法是将库伦规范作为约束条件，并通过拉格朗日乘子将其引入电磁场泛函表达式，从而建立含磁矢势和拉格朗日乘子自由度的电磁场棱边有限元方程，获得合理的解。采用拉格朗日乘子法引入库伦规范条件保证了解的唯一性，但是由于增加了拉格朗日乘子自由度，使得方程包含了拉格朗日乘子与磁矢势的耦合矩阵，使得矩阵 规模变大，求解效率降低。
第四种方法，是节点元法并结合罚函数法施加库伦规范的方法，这种方法，是对电磁场方程采用节点元法离散，为了保证磁矢势的唯一性，将库伦规范条件作为罚项引入电磁场泛函中，获得唯一的电磁场解，方程中没有增加多余的自由度，仅通过罚项对应的矩阵达到对磁矢势的约束作用；采用节点元法结合罚函数法施加库伦规范的方法虽然有效的保证了磁矢势的唯一性，但是采用节点元法时，在处理磁导率出现较大变化的材料分界面上解的精度较差。
发明内容
鉴于已有技术存在的缺陷，本发明的目的是要提供一种棱边单元法结合罚函数法施加库伦规范的方法，在保证精度的前提下，有效的解决了棱边单元法中的磁矢势唯一性问题。
为了实现上述目的，本发明的技术方案：
一种在电磁场棱边单元法中施加规范的方法，其特征在于：
包括如下步骤：
⑴、建立待分析对象三维模型，并对所述三维模型进行网格剖分；
⑵、计算剖分后的棱边单元的磁矢势双旋度项所对应的单元矩阵、库伦规范项所对应的单元矩阵以及电流密度项所对应的单元矩阵，并对它们进行集成，形成总体有限元方程；
⑶、对步骤⑵所形成的总体有限元方程，施加对应的边界条件；即在总体有限元方程的强迫边界上迫使剖分后的棱边单元的磁矢势自由度为0；
⑷、对步骤⑶所形成的总体有限元方程进行求解，获得唯一的磁矢势基本解，并对获得的磁矢势基本解进行后处理，得到对应的电磁场量的解。
所述对磁矢势双旋度项的单元矩阵进行计算，包括：
不导电空间中的静磁场方程：
▿×[v]▿×A=Jsin Ω---(1)]]>
A×n=0on ∂Ω---(2)]]>
其中，[ν]:磁阻矩阵(磁导率矩阵的逆矩阵)；A:磁矢势；J_s:电流密度矢量源项；旋度(curl)计算符号；Ω:不导电的三维空间领域；n:分析领域的外法向；强迫边界；为磁矢势的双旋度项；
对于满足边界条件(2)的任意磁矢势A，建立对应方程(1)的泛函：
F(A)=12∫Ω[v]|▿×A|2dV-∫ΩJs·AdV---(3)]]>
最小化F(A)，即可得到(1)、(2)所对应的微分方程解；
同时单元内的磁矢势通过形函数表示为：
A＝W·A_e                (4)
其中，W是单元磁矢势矢量形函数矩阵；A_e是单元棱边上的磁矢势自由度；
根据公式(3)以及(4)，
棱边单元的有限元方程表示为：
KAAAe=Jesin Ω---(5)]]>
其中，
KAA=∫Ve(▿×WT)T[v](▿×WT)dV---(6)]]>
即为所要计算的双旋度项单元矩阵，为单元内离散后的电流密度项。
所述库伦规范项的单元矩阵计算：
鉴于方程(1)、(2)中的磁矢势A具有不唯一性，需要施加额外的约束条件来保证其唯一性，采用库伦规范条件：
▿·A=0in Ω---(7)]]>
并将库伦规范条件(7)作为一个罚项引入到方程(1)、(2)所对应的泛函(3)中，最小化泛函，结合公式(4)，得到相应的有限元方程为：
(KAA+f·GAA)Ae=Jesin Ω---(8)]]>
其中，K^AA表示双旋度项的单元矩阵，f·G^AA为所引入的罚项即库伦规范项所对应的单元矩阵，f为罚因子，需要按照不同情况做适当的选取；
G^AA为磁矢库伦规范所对应的罚项的单元矩阵：其表达式为：
GAA=∫Ve(▿·WT)T[v](▿·WT)dV---(9)]]>
所述电流密度项对应的单元矩阵的计算：
计算公式(1)中的单元电流密度源项，即公式(5)中的
Jes=∫VeWTJSdV---(10)]]>
最后，对上述的单元矩阵K^AA、G^AA、进行集成，形成总体有限元方程。
对步骤⑵所形成的总体有限元方程，施加对应的边界条件；即在总体有限元方程的强迫边界上迫使剖分后的棱边单元的磁矢势自由度为0；
对步骤⑶所形成的总体有限元方程进行求解，获得唯一的磁矢势基本解，并对获得的磁矢势基本解进行后处理，得到对应的电磁场量的解。
与现有技术相比，本发明的有益效果：
本发明提出了一种新型的在电磁场棱边单元法中施加规范的方法即采用棱边单元法结合罚函数法施加库伦规范的方法，在保证精度的前提下，有效的解决了棱边单元法中的磁矢势唯一性问题。该方法精度高于节点元法，避免了传统方法如第一种方法中的复杂生成树的构建问题；同时效率高于棱边单元中所采用的拉格朗日乘子法施加库伦规范的方法，避免了第二种方法中的收敛性问题以及避免了第三种方法中的由于自由度增多，导致的求解效率低等问题。
附图说明
图1为现有技术棱边单元法结合树-余树方法网格空间Spanning Tree示意图；
图2为本发明方法步骤流程图；
图3本发明实施例待求解区域的三维几何模型示意图；
图4为本发明实施例空气包围铁区的有限元模型示意图；
图5为图4计算得到的磁感应强度矢量分布图；
图6为图4计算得到的磁感应强度矢量分布图--俯视图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图，对本发明进行进一步详细说明。
如图2所示，本发明所述解决现有技术存在的缺陷的技术方案具体是指在电磁场棱边单元法中施加规范的方法，包括如下步骤：
⑴、建立待分析对象三维模型，并对所述三维模型进行网格剖分；
⑵、计算剖分后的棱边单元的磁矢势双旋度项所对应的单元矩阵、库伦规范项所对应的单元矩阵以及电流密度项所对应的单元矩阵，并对它们进行集成，形成总体有限元方程；
进一步的，所述对磁矢势双旋度项的单元矩阵进行计算，包括：
不导电空间中的静磁场方程：
▿×[v]▿×A=Jsin Ω---(1)]]>
A×n=0on ∂Ω---(2)]]>
其中，[ν]:磁阻矩阵(磁导率矩阵的逆矩阵)；A:磁矢势；J_s:电流密度矢量源项；旋度(curl)计算符号；Ω:不导电的三维空间领域；n:分析领域 的外法向；强迫边界；为磁矢势的双旋度项；
对于满足边界条件(2)的任意磁矢势A，建立对应方程(1)的泛函：
F(A)=12∫Ω[v]|▿×A|2dV-∫ΩJs·AdV---(3)]]>
最小化F(A)，即可得到(1)、(2)所对应的微分方程解；
同时磁矢势的棱边单元形函数表示为：
A＝W·A_e                (4)
其中，W是磁矢势矢量形函数矩阵；A_e是单元棱边上的磁矢势自由度；
根据公式(3)以及(4)，
棱边单元的有限元方程表示为：
KAAAe=Jesin Ω---(5)]]>
其中，
KAA=∫Ve(▿×WT)T[v](▿×WT)dV---(6)]]>
即为所要计算的双旋度项单元矩阵，为离散后的电流密度项。
所述库伦规范项的单元矩阵计算：
鉴于方程(1)、(2)中的磁矢势A具有不唯一性，需要施加额外的约束条件来保证其唯一性，采用库伦规范条件：
▿·A=0in Ω---(7)]]>
并将库伦规范条件(7)作为一个罚项引入到方程(1)、(2)所对应的泛函(3)中，最小化泛函，结合公式(4)，得到相应的有限元方程为：
(KAA+f·GAA)Ae=Jesin Ω---(8)]]>
其中，K^AA表示双旋度项的单元矩阵，f·G^AA为所引入的罚项即库伦规范项，f为罚因子，需要按照不同情况做适当的选取；若罚因子过大，则导致得出的解精度下降，若罚因子过小，则导致方程组的病态。
G^AA为磁矢库伦规范所对应的罚项的单元矩阵：其表达式为：
GAA=∫Ve(▿·WT)T[v](▿·WT)dV---(9)]]>
所述电流密度项对应的单元矩阵的计算：
计算公式(1)中的单元电流密度源项，即公式(5)中的
Jes=∫VeWTJSdV---(10)]]>
最后，对上述的单元矩阵K^AA、G^AA、进行集成，形成总体有限元方程。
⑶、对步骤⑵所形成的总体有限元方程，施加对应的边界条件；即在总体 有限元方程的强迫边界上迫使剖分后的棱边单元的磁矢势自由度为0；
⑷、对步骤⑶所形成的总体有限元方程进行求解，获得唯一的磁矢势基本解，并对获得的磁矢势基本解进行后处理，得到对应的电磁场量的解。
下面以文献中的空气包围铁芯的算例(具体算例信息详见：M.Kaltenbacher,Member,IEEE,and S.Reitzinger Appropriate Finite-Element Formulations for3-D Electromagnetic-Field Problems IEEE TRANSACTIONS ON MAGNETICS,VOL.38,NO.2,MARCH2002)来说明上述的操作过程：
建立待求解区域的三维几何模型，并对其进行网格剖分，如图3所示，黑色区域为铁区，其余区域为空气区域；
创建对应的物理模型，并对铁区和空气区域赋材料属性(相对磁导率)，以及单元属性；
施加物理模型的边界条件，使通过铁区和空气区域垂直方向的磁通量为常数，其他方向为0，如图4所示；
采用棱边单元法对模型的求解区域进行棱边单元法离散，计算双旋度项、罚项以及电流项对应的单元矩阵，并进行集成，得到对应的总体有限元方程；
对集成得到的总体有限元方程上施加中物理边界条件所对应的约束条件；
求解该方程，获得磁矢势基本解，进行后处理，得到对应电磁场量，如图5、6所示，所得结果合理正确，显示出该方法的有效性。
以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，根据本发明的技术方案及其发明构思加以等同替换或改变，都应涵盖在本发明的保护范围之内。