一种电磁场非匹配棱边元的插值方法

技术领域

本发明涉及计算辅助工程技术领域，具体而言，尤其涉及一种电磁场非匹配棱边元的插值方法。

背景技术

跨区域建模方式是一种灵活的建模方式，常用于工程中多组件、多工况的多仿真场景，也常见于电磁场运动滑移和区域分解等问题中。在跨区域建模中，可以对各个组件的子域独立建模、网格划分甚至求解。当子域的几何或网格尺寸发生改变时，只需更换相应子域的几何或网格而不影响其它子域。在跨区域建模中，需要在不同子域的网格分界面处建立耦合约束，以保证未知量在边界处的连续性。以电磁场运动滑移问题为例，一般先将求解区域分为运动区域和静止区域，对它们分别进行独立的网格剖分，在两套网格分界面处建立滑移面，当网格发生运动时，滑移面两侧的网格节点一般不匹配，此时需要通过自由度耦合技术来施加滑移面处物理量的连续性条件，常见的自由度耦合方法有Lagrangian乘子法、罚方法及多点约束法。

Lagrangian乘子法是处理含约束类问题的一种通用方法，尤其适用于约束条件不能显式表示时。为了保证滑移面上的非匹配网格上未知量的连续性，通过Lagrangian乘子在泛函中引入这个约束关系表达式，从而将有附加条件的变分原理变成无附加条件的变分原理。然后再对这个修正后的变分原理取驻值，就可以得到关于主、从自由度以及Lagrangian乘子的有限元方程。需要注意的是，对于三维电磁棱边元，是引入Lagrangian矢量乘子来施加未知场量切向分量的连续性条件。Lagrangian法的缺点是会导致病态非定的离散系统。

罚函数法是另外一种有效处理约束的方法，它将耦合约束条件以乘积的形式引入泛函，得到修正的变分原理，并通过罚因子来控制约束刚度。然后对修正后的变分原理取驻值，就可以得到关于主、从自由度的有限元方程。相比于Lagrangian乘子法，罚函数方法不需要引入额外的耦合界面变量，求解矩阵的性态跟罚因子的选取有较大关系。

多点约束法是另外一种处理非匹配网格的有效方法，广泛应用于实际工程。基于主从面上网格的映射插值技术，它显式地将从面节点自由度表示为它所映射的主面单元上节点自由度的线性组合。文献中的节点插值法是多点约束法在电机等问题中的具体应用。不过对于三维棱边元来说，问题变得有些复杂，需要考虑棱边的映射插值，将从面棱边自由度表示为主面棱边自由度的线性组合。

Golovanov等人在中应是第一次提出连接棱边非匹配网格的插值方法，该方法对于未知从棱边在对应主面棱边上实现线性插值，可以认为是多点约束方法从节点元到棱边元上的推广。但是该方法对于一般性网格存在棱边分割和搜寻相关主棱边的困难，因此该方法限定于含规则网格的数值案例的计算分析。Kometani等人基于该方法结合规则网格提出了一种插值算法，成功应用于三维带斜槽的感应电机的分析；Okamoto等人也基于该方法使用规则网格分析了含旋转运动的电磁波问题。Muramatsu等人还提出了在分界面处细分主面和从面上的网格来生成匹配网格的有趣思路。但是，这种细分程序在拓展到任意三维非匹配网格中可能还存在一些技术难题，在所提供的数值算例中，使用的仍然是规则网格。

Ito等人在非匹配网格问题中对从棱边提出了一种棱边分割算法，他们将从棱边和主单元的交叉形式分为四类，并在插值公式中解析地计算了分割棱边的切线积分表达式。但是他们并未考虑工程中经常出现的交界面为曲面的情况。因此，该方法在应用于曲面时可能会由于缝隙或重叠区域的存在而失效。另一个困难是缺乏一个高效地搜索算法，当非匹配网格数量变得非常庞大时，将消耗较大的计算资源。

由于非匹配棱边元的插值算法尚存在的一些不足以及技术实现本身的难度，最终实现的非匹配棱边插值算法还应满足通用性，尤其是适应于交界面为曲面的情形；还应满足复杂情形下插值的精度要求，为了保证非匹配网格中磁矢势求解，应当跨区域构造合理的生成树；还应考虑效率问题，当网格数庞大时，应当保证棱边的映射插值是高效的。

此外，处理电磁场非匹配网格问题还有许多其它方法，如体素建模方法，Galerkin投影方法，在电机旋转问题的建模中还有气隙单元法，根据定、转子在气隙处的未知自由度的解析关系，将气隙处理为一种特殊的宏单元，从而耦合定、转子间的未知自由度，此外还有有限元-边界元耦合方法，非重叠Mortar有限元法方法等。

发明内容

根据上述提出由于非匹配棱边元的插值算法尚存在的一些不足以及技术实现本身的难度，最终实现的非匹配棱边插值算法还应满足通用性等的技术问题，而提供一种电磁场非匹配棱边元的插值方法。本发明主要利用一种电磁场非匹配棱边元的插值方法，其特征在于，包括以下步骤：

步骤S1：对于电磁场仿真中的非匹配网格，在分界面处将从棱边分割为若干段，每一段分割棱边隶属所述耦合面上的一个主面单元；

步骤S2：对所述的每一段分割棱边上的场量线积分使用该主面单元内的棱边自由度的线性组合表示，从棱边自由度等于所有分割棱边的场量线积分的总和。

进一步地，所述步骤S1还包括以下步骤：

步骤S11：通过桶搜索算法将从面节点映射到主面单元上得到投影节点，根据投影节点的连线获得投影棱边；

步骤S12：计算每条投影棱边与主面网格的交叉点，依次连接投影棱边的交叉点获得若干段分割棱边。

更进一步地，所述的桶搜索算法，首先将主面单元按照坐标范围分配到若干数量的桶中，且桶内的主面单元数目大致相等；对所述从面节点，按照坐标范围找到它所从属的桶，再在该桶内逐个搜索出匹配的主面单元。

进一步地，所述步骤S12中的计算每条投影棱边与主面网格的交叉点，通过使交叉点与投影棱边的末端点的距离最小的方法来确定交叉点。

更进一步地，对分割棱边上的场量线积分的计算使用两点Gauss积分的计算方法。

较现有技术相比，本发明具有以下优点：

本发明可以进行电磁场跨区域建模，如电机、电磁发射技术中的滑移面耦合分析，实现不同子域的网格分界面处棱边场量的连续性条件。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图做以简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。

图1为本发明从面单元PQR的棱边在主面网格上的分割示意图。

图2为本发明从节点在主面网格上的映射示意图。

图3为本发明从棱边PQ在主面网格上的分割示意图。

图4为本发明对分割棱边的场量采用数值积分示意图。

图5为本发明方腔模型中的非匹配网格和磁场分布示意图；其中，(a)为方腔非匹配网格有限元模型；(b)为磁场强度矢量图。

图6为本发明非匹配网格界面上从棱边的分割示意图；其中，(a)为从棱边的分割；(b)为局部放大图。

图7为本发明旋转直线圈模型示意图；其中，(a)为初始时刻网格；(b)为初始时刻磁场分布。

图8为本发明初始时刻非匹配网格界面上的从棱边分割示意图。

图9为本发明旋转直线圈模型的磁场能量随时间变化曲线示意图。

图10为本发明三维电机有限元模型示意图。

图11为本发明转子的磁力转矩随时间的变化示意图。

具体实施方式

为了使本技术领域的人员更好地理解本发明方案，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分的实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都应当属于本发明保护的范围。

需要说明的是，本发明的说明书和权利要求书及上述附图中的术语“第一”、“第二”等是用于区别类似的对象，而不必用于描述特定的顺序或先后次序。应该理解这样使用的数据在适当情况下可以互换，以便这里描述的本发明的实施例能够以除了在这里图示或描述的那些以外的顺序实施。此外，术语“包括”和“具有”以及他们的任何变形，意图在于覆盖不排他的包含，例如，包含了一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元，而是可包括没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。

如图1-11所示，本发明提供了一种电磁场非匹配棱边元的插值方法，包括以下步骤：步骤S1：对于电磁场仿真中的非匹配网格，在分界面处将从棱边分割为若干段，每一段分割棱边隶属所述耦合面上的一个主面单元；步骤S2：对所述的每一段分割棱边上的场量线积分使用该主面单元内的棱边自由度的线性组合表示，从棱边自由度等于所有分割棱边的场量线积分的总和。

具体的，如图1所示，三角形PQR为从面单元，棱边PQ在从面上的分割棱边由AB、BC、CD、DE这四个线段组成，交叉点为B、C、D。

为了确定主面上发生交叉的棱边以及交叉点的位置，最简单的算法是对每条从面棱边遍历搜索所有的主面网格。但是当网格量很多时，这种方法的计算代价较大；另一个困难是当耦合面为曲面时，主、从网格将存在间隙或重叠，应如何判断交叉问题。

现以图1中棱边PQ的分割为例进行算法说明。步骤S11：通过桶搜索算法将从面节点映射到主面单元上得到投影节点，根据投影节点的连线获得投影棱边；即图2中PQ在主面上的映射点A和E。

步骤S12：计算每条投影棱边与主面网格的交叉点，依次连接投影棱边的交叉点获得若干段分割棱边。且所述步骤S12中的计算每条投影棱边与主面网格的交叉点，通过使交叉点与投影棱边的末端点的距离最小的方法来确定交叉点。

如图3所示，从棱边AE的A端出发，第一个交叉点B应在三角形123的棱边上。我们通过最小化距离和|AB|+|BE|的方法计算出B在三角形123棱边上的位置，图3所示的B点在主棱边23上，并获得第一段分割棱边AB。根据主面单元之间的连接关系以及B点所在的共享棱边编号，可以查找出与三角形123相邻的面单元234。接着从B点出发，在三角形234内确定新的交叉点C，使距离和|BC|+|CE|最小，从而确定如图3所示的C点位置，同时获得第二段分割棱边BC。重复此步骤，可以获得余下的交叉点D及分割棱边CD。当交叉点进入末端附着点E所在单元456时，直接提取最后一段分割棱边DE，结束从棱边的分割程序。最终，我们获得从棱边PQ在主面网格上的全部分割棱边，它们是图3所示的一系列相连线段AB、BC、CD、DE，分别隶属主面单元123、234、345、456。

更进一步地，所述的桶搜索算法，首先将主面单元按照坐标范围分配到若干数量的桶中，且桶内的主面单元数目大致相等；对所述从面节点，按照坐标范围找到它所从属的桶，再在该桶内逐个搜索出匹配的主面单元。

如图4，每段分割棱边的积分采用两点Gauss积分，其中“×”是Gauss点的位置标识。对分割棱边上的场量线积分的计算使用两点Gauss积分的计算方法。

实施例1

以1×1×1的正方体区域沿着某个面法向通以大小为1.0×10

表1方腔的磁场能量的对比

实施例2

如图7所示，为含解析解的直线圈电流激励磁场的案例，直线圈和外部空气的相对磁导率均为1，直线圈通以恒定电流密度1A·m

实施例3

如图10所示为一个4极6槽表贴永磁同步电机的三维有限元模型。对该三维电机建模时，将电机分为定子相关区域Ω

定子中的线圈通有三相电流U-V-W，匝数是100，采用串联方式，电流有效值是11.31A，频率为100Hz；转子中放置四块永磁铁，矫顽力大小为898404A/m。定子和转子设为相对磁导率为1000的铁性材料。对该电机稳定工作时的磁场进行仿真分析。设定分析的总时间为0.005s，时间步长为0.0001s，转子的转速为3000r/min。

图11给出了基于滑移网格耦合的电机旋转过程的转矩履历，并与ANSYS-Maxwell的计算结果进行对比。可以看出两者的结果基本一致，说明了三维滑移网格耦合算法的有效性。

上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。

在本发明的上述实施例中，对各个实施例的描述都各有侧重，某个实施例中没有详述的部分，可以参见其他实施例的相关描述。

在本申请所提供的几个实施例中，应该理解到，所揭露的技术内容，可通过其它的方式实现。其中，以上所描述的装置实施例仅仅是示意性的，例如所述单元的划分，可以为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些接口，单元或模块的间接耦合或通信连接，可以是电性或其它的形式。

所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。

另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。上述集成的单元既可以采用硬件的形式实现，也可以采用软件功能单元的形式实现。

所述集成的单元如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的全部或部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可为个人计算机、服务器或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：U盘、只读存储器(ROM，Read-Only Memory)、随机存取存储器(RAM，Random Access Memory)、移动硬盘、磁碟或者光盘等各种可以存储程序代码的介质。

最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。