基于解析法的无铁芯永磁同步直线电机磁场分析

Magnetic field analysis of ironless permanent magnet synchronous linear motor based on analytical method

摘要

本文旨在实现双边无铁芯永磁同步直线电机（DSILPMSLM）磁场参数的快速精确计算， 为DSILPMSLM的设计与优化提供高效的解析工具. 有限元分析方法结果精确， 但计算时间长且需大量计算资源， 无法建立磁场与电机结构参数之间的数学关系， 解析法可以有效避免有限元法求解时间过长的问题， 通过建立准确的数学模型， 实现磁场参数的快速计算. 本文首先基于等效磁化电流原理， 采用解析法对一台DSILPMSLM进行了磁场解析分析， 推导出了空载气隙磁通密度表达式. 在此基础上， 获得了非重叠绕组DSILPMSLM的空载相反电动势的精确解析式， 完成了对电机推力的建模求解. 其次， 将解析法计算结果与有限元仿真进行了对比分析. 最后， 通过搭建实验平台， 对研究样机的空载反电动势以及电机推力进行了实验测试. 对比显示， 解析法计算得到的气隙磁通密度、三相反电动势及电机平均推力与有限元法结果和实验数据高度吻合， 证明了该方法在解析DSILPMSLM磁场时的准确性. 本文采用的解析法在分析DSILPMSLM磁场时具有很高的精度， 显著缩短了计算时间， 对高性能永磁同步直线电机的设计与优化具有重要的指导意义.

Abstract

This study aims to achieve fast and accurate calculation of the magnetic field parameters in a double-sided ironless permanent magnet synchronous linear motor （DSILPMSLM） and to provide an efficient analytical tool for its design and optimization. Although the finite element method （FEM） provides accurate results， it requires substantial computational resources and time， and is unable to construct a mathematical relationship between the magnetic field and motor structure parameters. The analytical method can overcome the drawback of long computation time associated with FEM by establishing a precise mathematical model for rapid calculation. Firstly， based on the equivalent magnetization current principle， the magnetic field of a DSILPMSLM was analytically modeled. Expressions for the no-load air-gap flux density and no-load back electromotive force （EMF） of non-overlapping winding DSILPMSLM were derived， and a thrust model was developed and solved. Secondly， the results obtained from the analytical method were compared with those from FEM simulations and experimental measurements. Thirdly， an experimental platform was set up to test the no-load back EMF and thrust of the research prototype. The air-gap flux density， three-phase back EMF， and average motor thrust calculated by the analytical method show close agreement with both FEM results and experimental data， confirming the accuracy of the proposed method in predicting the magnetic field of DSILPMSLM. The analytical method employed in this study demonstrates high accuracy and significantly reduces computation time. It provides important guidance for the design and optimization of high-performance permanent magnet synchronous linear motors.

Graphical abstract

关键词

磁场分析 / 有限元法 / 数学模型 / 永磁电机

Key words

magnetic field analysis / finite element method / mathematical models / permanent magnet motors

永磁同步直线电机具有结构简单，高加速度、高推力密度、高定位精度等优点^［1］，目前广泛应用于工业生产、轨道交通、医疗系统等各种场合^［2］.双边无铁芯永磁同步直线电机（double-sided ironless permanent magnet synchronous linear motor， DSILPMSLM）因其无铁芯结构，具有零齿槽力、零法向吸力等独特优势，在高精度直线运动伺服系统中有巨大实用价值^［3］.

目前对永磁同步直线电机的研究集中在拓扑结构设计、电机控制算法^［4］、仿真建模分析、电磁特性分析等方面.永磁同步直线电机拓扑结构多种多样，磁场分布又较为复杂，为了得到精确的磁路计算结果，需要对电机建立准确的电磁模型.当前，广泛应用于电磁场的计算方法是等效磁路法^［5-8］，但是等效磁路法的缺点是仅对磁场参数进行单一的线性化求解， 未考虑漏磁效应和饱和效应， 且难以处理电机复杂的几何结构， 这将导致磁场的计算结果与磁场实际分布存在较大偏差. 有限元法^［9］支持各种复杂的边界条件， 将仿真模型剖分为细小的单元， 计算了每一个单元的近似解， 极大地提高了磁场的解析精度. 但是有限元法求解过程耗费大量时间， 计算资源需求过高. 可以采用解析法^［10］以避免上述方法的缺点， 解析法计算效率远高于有限元法，其通过数学公式揭示了磁场参数之间的物理规律， 从而方便对磁场进行优化分析. 此外， 解析法更适合电机设计中的初始迭代计算. 文献［11］用改进型的等效磁路法对有铁芯永磁同步直线电机的气隙磁场建模分析， 考虑了铁芯的饱和情况和端部效应， 提高了解析法的计算精度.

文献［12］用等效磁荷法分析了平面电机永磁阵列的三维气隙磁通密度分布， 简化了磁场分析模型. 文献［13］提出了一种三边型无铁芯永磁同步直线电机， 用等效磁化电流法对气隙磁场进行建模分析， 得到磁通密度、反电动势和推力表达式. 文献［14］用等效磁网格法计算了有铁芯永磁同步直线电机不同负载条件下的气隙磁通密度，研究了饱和电枢反应对气隙磁通密度分布的影响.文献［15-18］用解析法计算了有铁芯永磁同步直线电机的气隙磁场， 得到了气隙磁通密度的表达式.文献［19］用Maxwell应力张量法对无铁芯动磁式永磁同步直线电机建立了推力模型，搭建了推力实验测试系统，验证了该方法计算推力的准确性.文献［20］提出一种子域叠加模型分析电机电枢反应磁场，验证了该方法的可行性.文献［21］提出了一种快速三维分析模型来解析具有倾斜拓扑结构永磁体的磁场分量.文献［22］采用绕组函数法对永磁电机磁场进行分析，并得到负载和空载下的不同特性. 文献［23］提出了一种磁场半解析技术来模拟三维永磁体的力分量.文献［24］采用圆环电流法对气隙磁场进行解析计算， 在保证精度的同时提高了计算速度. 文献［25］提出了一种极坐标系下磁场解析方法来解析三维永磁体的气隙磁场分量.

本文在上述磁场解析的基础上，引入等效磁化电流模型，采用解析法分析双边无铁芯永磁同步直线电机的磁场.具体工作如下：首先，介绍了DSILPMSLM的拓扑结构并分析了其工作原理；其次，引入矢量磁位，建立了永磁体区域和气隙区域的等效磁化电流解析模型，用分离变量法推导了气隙磁通密度（下文简称磁密）的表达式；再次，以一台8极6槽非重叠绕组的DSILPMSLM为研究对象，结合磁密表达式获得了其空载相反电动势的精确解析式，基于磁场解析结果， 对该款电机推力进行建模求解；最后， 将有限元法和解析法进行对比分析， 通过空载反电动势及推力实验证明了该方法在分析DSILPMSLM时的准确性.

1 DSILPMSLM的基本结构

本文所研究的DSILPMSLM基本结构如图1所示， 该电机主要由动子和定子两个模块组成. 动子由三相绕组、环氧树脂及铝壳组成. 三相绕组的单个线圈绕制后经浸漆处理并采用模具压缩使之成型， 线圈经过环氧树脂封装固定， 通过铝壳引出绕组引线. 三相绕组的空间电角度相差120°. 定子由导磁性较好的背铁和表贴式长方形永磁体构成， 永磁体以均匀间隙安装在背铁上， 具有交错的磁极和对称的拓扑结构， 三相绕组与上下永磁体之间留有相同的机械气隙. 在整体装配中， 动子绕组模块被安装在上下两层永磁体阵列之间， 两者保持对称的气隙间距， 构成“动子-气隙-永磁体”的层叠耦合结构. 当对称的三相正弦交流电流通入三相绕组后， 电流产生的磁场和永磁体产生的励磁磁场相互作用形成气隙磁场， 动子线圈在气隙磁场的作用下在电机内部做直线运动.

2 磁场解析分析

2.1 等效磁化电流解析模型

等效磁化电流法是电磁场分析中用于简化磁性材料磁场计算的一种物理方法， 其核心思想是将永磁体产生的励磁磁场等效为束缚面电流或者体电流， 通过引入矢量磁位， 然后根据已知的边界条件求解矢量磁位满足的拉普拉斯方程或者泊松方程. 本文所研究的DSILPMSLM具有对称结构， 分析电机一侧的磁场物理量便可以等效成整体. 本文选取电机x轴以上永磁体的部分进行研究分析， DSILPMSLM的等效磁化电流空间分布函数M（x）如图2所示.

图2中， {\tau }_{\mathrm{m}}为永磁体的宽度， \tau为极距， M_{\mathrm{0}}为磁化电流的磁化强度.

根据图2， 等效磁化电流分布函数M（x）由对称性可知：

M(x)=\left\{\begin{array}{l}\mathrm{0},k\tau -\frac{\tau -{\tau }_{\mathrm{m}}}{\mathrm{2}}\le x\le k\tau +\frac{\tau -{\tau }_{\mathrm{m}}}{\mathrm{2}}\\ {(-\mathrm{1})}^k{M}_{\mathrm{0}},k\tau +\frac{\tau -{\tau }_{\mathrm{m}}}{\mathrm{2}}\le x\le k\tau +\frac{\tau +{\tau }_{\mathrm{m}}}{\mathrm{2}}\end{array}\right.（1）

式中：k=\mathrm{0},\pm \mathrm{1},\pm \mathrm{2},\dots.

等效电流密度分布函数J_{\mathrm{m}}(x)为：

J_{\mathrm{m}}(x)=\frac{\partial M(x)}{\partial x}（2）

M_{\mathrm{0}}=\frac{B_{\mathrm{r}}}{{\mu }_{\mathrm{0}}}（3）

由图2可知M（x）为周期性方波函数， 可以得到M（x）的傅里叶级数展开式：

M(x)={\displaystyle \sum _{k=\mathrm{1}}^{\infty }}{M}_{k}sin\left[\frac{\left(\mathrm{2}k-\mathrm{1}\right)\mathrm{\pi }}{\tau }x\right]（4）

根据傅里叶级数中三角函数展开式正弦项的计算公式（5）可以求得M（x）对应的傅里叶系数M_k：

b_n=\frac{\mathrm{2}}{T_{\mathrm{0}}}{\int }_{\frac{-T_{\mathrm{0}}}{\mathrm{2}}}^{\frac{T_{\mathrm{0}}}{\mathrm{2}}}X(t)\mathrm{s}\mathrm{i}\mathrm{n}{\omega }_{\mathrm{0}}t\mathrm{d}t（5）

M_k=\frac{\mathrm{2}}{\tau }{\int }_{\frac{\tau -{\tau }_{\mathrm{m}}}{\mathrm{2}}}^{\frac{\tau +{\tau }_{\mathrm{m}}}{\mathrm{2}}}{M}_{\mathrm{0}}\mathrm{s}\mathrm{i}\mathrm{n}\frac{(\mathrm{2}k-\mathrm{1})\mathrm{\pi }x}{\tau }\mathrm{d}x（6）

经过积分可以得到：

M_k={(-\mathrm{1})}^{k+\mathrm{1}}{M}_{\mathrm{0}}\frac{\mathrm{4}}{(\mathrm{2}k-\mathrm{1})\mathrm{\pi }}\mathrm{s}\mathrm{i}\mathrm{n}\frac{(\mathrm{2}k-\mathrm{1})\mathrm{\pi }{\tau }_{\mathrm{m}}}{\mathrm{2}\tau }（7）

2.2 气隙磁密解析计算

为了方便对DSILPMSLM的空载气隙磁密进行解析计算， 将整个区域划分为气隙区域和永磁体区域， DSILPMSLM的解析模型如图3所示，其中h_{\mathrm{m}}为永磁体高度， \delta为气隙总高度，g为有效气隙长度.

实际情况下，DSILPMSLM的内部磁场对外呈现三维状态，但由于三维磁场描述较为复杂且方程的未知数过多，给求解带来了较大的困难.为了方便对DSILPMSLM气隙磁场的求解，将电机的三维问题转化为二维问题，提出以下假设：永磁体均匀充磁， 退磁曲线和回复线重合且都为直线；电机沿x方向无限长， 磁场强度沿x方向无衰减；背铁无磁阻， 磁导率无限大；电机端部磁场不发生畸变.

在双边无铁芯永磁同步直线电机中， 关于磁场的研究主要包括两个部分， 即永磁体产生的励磁磁场和气隙磁场.

在气隙区域（区域1）：

B_{\mathrm{1}}={\mu }_{\mathrm{0}}{H}_{\mathrm{1}}（8）

在永磁体区域（区域2和区域3）：

B_{\mathrm{2}}={\mu }_{\mathrm{0}}(M+{\mu }_{\mathrm{r}}{H}_{\mathrm{2}})（9）

式中：B_{\mathrm{1}}和H_{\mathrm{1}}分别是气隙区域的磁感应强度和磁场强度；B_{\mathrm{2}}和H_{\mathrm{2}}分别是永磁体区域的磁感应强度和磁场强度. {\mu }_{\mathrm{0}}是真空磁导率即空气磁导率；{\mu }_{\mathrm{r}}是永磁体相对磁导率；M是永磁体剩余磁化强度.

引入矢量磁位A_{z\mathrm{1}}(x,y)， A_{z\mathrm{2}}(x,y)， A_{z\mathrm{3}}(x,y)， 由于在气隙区域1中不存在磁化电流， 而在永磁体区域2和区域3中存在磁化电流， 因此使用拉普拉斯方程求解气隙区域1中的矢量磁位A_{z\mathrm{1}}(x,y)， 使用泊松方程求解区域2和区域3的矢量磁位A_{z\mathrm{2}}(x,y)和A_{z\mathrm{3}}(x,y)：

\frac{{\partial }^{\mathrm{2}}{A}_{z\mathrm{1}}}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}{A}_{z\mathrm{1}}}{\partial y^{\mathrm{2}}}=\mathrm{0}（10）

\frac{{\partial }^{\mathrm{2}}{A}_{z\mathrm{2}}}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}{A}_{z\mathrm{2}}}{\partial y^{\mathrm{2}}}=-{\mu }_{\mathrm{0}}\frac{\partial M_{\mathrm{2}}(x)}{\partial x}（11）

\frac{{\partial }^{\mathrm{2}}{A}_{z\mathrm{3}}}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}{A}_{z\mathrm{3}}}{\partial y^{\mathrm{2}}}=-{\mu }_{\mathrm{0}}\frac{\partial M_{\mathrm{3}}(x)}{\partial x}（12）

采用分离变量法求得A_{z\mathrm{1}}(x,y)和A_{z\mathrm{2}}(x,y)的通解：

\begin{array}{l}{A}_{z\mathrm{1}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}[(A_{k\mathrm{1}}\mathrm{c}\mathrm{o}\mathrm{s}{m}_{k}x+B_{k\mathrm{1}}\mathrm{s}\mathrm{i}\mathrm{n}{m}_{k}x)\times \\ {}_{\begin{array}{cc}& \end{array}}^{}{}_{\begin{array}{cc}& \end{array}}^{}{}_{\begin{array}{cc}& \end{array}}^{}{}_{\begin{array}{cc}& \end{array}}^{}{}_{\begin{array}{cc}& \end{array}}^{}{}^{}(C_{k\mathrm{1}}\mathrm{s}\mathrm{h} m_{k}y+D_{k\mathrm{1}}\mathrm{c}\mathrm{h} m_{k}y)]\end{array}（13）

A_{z\mathrm{2}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}\left[\begin{array}{l}(A_{k\mathrm{2}}\mathrm{c}\mathrm{o}\mathrm{s}{m}_{k}x+B_{k\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}{m}_{k}x)\times \\ (C_{k\mathrm{2}}\mathrm{s}\mathrm{h}{m}_{k}y+D_{k\mathrm{2}}\mathrm{c}\mathrm{h}{m}_{k}y)+\\ ({\mu }_{\mathrm{0}}{M}_k/m_k)\mathrm{c}\mathrm{o}\mathrm{s}{m}_{k}x]\end{array}\right]（14）

式中：m_k为分离常数.

m_k=\frac{\left(\mathrm{2}k-\mathrm{1}\right)\mathrm{\pi }}{\tau }（15）

矢量磁位函数A_{z\mathrm{1}}(x,y)和A_{z\mathrm{2}}(x,y)关于y轴对称， 电机自身求解定义域关于x轴对称， 则有：

A_z(x,-y)=A_z(x,y)（16）

A_z(-x,y)=A_z(x,y)（17）

代入通解可以化简求得：

A_{z\mathrm{1}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}\mathrm{c}\mathrm{o}\mathrm{s}{m}_{k}x\times D_{k\mathrm{3}}\times \mathrm{c}\mathrm{h} m_{k}y（18）

\begin{array}{l}{A}_{z\mathrm{2}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}\mathrm{c}\mathrm{o}\mathrm{s} m_{k}x(C_{k\mathrm{4}}\mathrm{s}\mathrm{h} m_{k}y+\\ {}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{D}_{k\mathrm{4}}\mathrm{c}\mathrm{h} m_{k}y+\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k})\end{array}（19）

式中：D_{k\mathrm{3}}、C_{k\mathrm{4}}、D_{k\mathrm{4}}为待定系数.

根据DSILPMSLM满足的三个边界条件，代入上述方程求得待定系数，通过对矢量磁位求偏导数可以得到不同区域的磁感应强度表达式，边界条件为：

\left\{\begin{array}{l}y=\frac{g}{\mathrm{2}},B_{y\mathrm{1}}(x,y)=B_{y\mathrm{2}}(x,y){,}^{}{}^{}{B}_y=-\frac{\partial A_z}{\partial x}\\ y=\frac{g}{\mathrm{2}},H_{x\mathrm{1}}(x,y)=H_{x\mathrm{2}}(x,y){,}^{}{H}_x=\frac{\mathrm{1}}{{\mu }_{\mathrm{0}}}\frac{\partial A_z}{\partial y}\\ y=\frac{\delta }{\mathrm{2}},B_{x\mathrm{2}}(x,y)=\mathrm{0}{,}^{}{}^{}{}^{}{}^{}{}^{}{}^{}{}^{}{}^{}{}^{}{}^{}{}^{}{B}_{x\mathrm{2}}=\frac{\partial A_{z\mathrm{2}}}{\partial y}\end{array}\right.（20）

式中：B_{y\mathrm{1}}、B_{y\mathrm{2}}分别为气隙区域和永磁体区域y方向上的磁密； H_{x\mathrm{1}}、H_{x\mathrm{2}}分别为气隙区域和永磁体区域x方向上的磁场强度； B_{x\mathrm{2}}为永磁体区域x方向上的磁密.

将边界条件代入式（18）和式（19）可以得到：

\left\{\begin{array}{l}{D}_{k\mathrm{3}}\mathrm{s}\mathrm{h}{m}_k\left(\frac{\delta }{\mathrm{2}}-h_m\right)=C_{k\mathrm{4}}\mathrm{c}\mathrm{h}{m}_k\left(\frac{\delta }{\mathrm{2}}-h_{\mathrm{m}}\right)+D_{k\mathrm{4}}\mathrm{s}\mathrm{h}{m}_k\left(\frac{\delta }{\mathrm{2}}-h_{\mathrm{m}}\right)\\ C_{k\mathrm{4}}\mathrm{c}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}+D_{k\mathrm{4}}\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}=\mathrm{0}\\ D_{k\mathrm{3}}\mathrm{c}\mathrm{h}{m}_k\left(\frac{\delta }{\mathrm{2}}-h_{\mathrm{m}}\right)=C_{k\mathrm{4}}\mathrm{s}\mathrm{h}{m}_k\left(\frac{\delta }{\mathrm{2}}-h_{\mathrm{m}}\right)+D_{k\mathrm{4}}\mathrm{c}\mathrm{h}{m}_k\left(\frac{\delta }{\mathrm{2}}-h_{\mathrm{m}}\right)+\\ {\begin{array}{ccc}\begin{array}{ccc}& & \end{array}& & \end{array}}^{}{}^{}{}^{}{}^{}{}^{}\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k}\end{array}\right.（21）

\left\{\begin{array}{l}{C}_{k\mathrm{4}}=-\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k}\mathrm{s}\mathrm{h}{m}_k\left(h_{\mathrm{m}}-\frac{\delta }{\mathrm{2}}\right)\\ D_{k\mathrm{4}}=\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k}\frac{\mathrm{c}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}{\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}\mathrm{s}\mathrm{h}{m}_k\left(h_{\mathrm{m}}-\frac{\delta }{\mathrm{2}}\right)\\ D_{k\mathrm{3}}=\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k}\frac{\mathrm{s}\mathrm{h}{m}_k{h}_{\mathrm{m}}}{\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}\end{array}\right.（22）

将待定系数D_{k\mathrm{3}}、C_{k\mathrm{4}}、D_{k\mathrm{4}}反代入式（18）和式（19）可以得到：

A_{z\mathrm{1}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k}\frac{\mathrm{s}\mathrm{h}{m}_k{h}_{\mathrm{m}}}{\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}\mathrm{c}\mathrm{o}\mathrm{s}{m}_{k}x\times \mathrm{c}\mathrm{h}{m}_{k}y（23）

A_{z\mathrm{2}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}\mathrm{c}\mathrm{o}\mathrm{s}{m}_{k}x\frac{{\mu }_{\mathrm{0}}{M}_k}{m_k}\left[{}^{\begin{array}{l}\\ \\ \\ \end{array}}-\mathrm{s}\mathrm{h}{m}_k\left(h_{\mathrm{m}}-\frac{\delta }{\mathrm{2}}\right)\mathrm{s}\mathrm{h}{m}_{k}y+\right.

\left.\mathrm{s}\mathrm{h}{m}_k\left(h_{\mathrm{m}}-\frac{\delta }{\mathrm{2}}\right)\frac{\mathrm{c}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}{\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}\mathrm{c}\mathrm{h}{m}_{k}y+\mathrm{1}\right]（24）

求矢量磁位函数在x方向上的偏导数便得到气隙区域和永磁体区域在y方向上的磁密表达式：

B_{y\mathrm{1}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}{\mu }_{\mathrm{0}}{M}_k\frac{\mathrm{s}\mathrm{h}{m}_k{h}_{\mathrm{m}}}{\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}\mathrm{s}\mathrm{i}\mathrm{n}{m}_{k}x\times \mathrm{c}\mathrm{h}{m}_{k}y（25）

\begin{array}{l}{B}_{y\mathrm{2}}(x,y)={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}{\mu }_{\mathrm{0}}{M}_k\left[{}^{\begin{array}{l}\\ \\ \\ \end{array}}-\mathrm{s}\mathrm{h}{m}_k\left(h_{\mathrm{m}}-\frac{\delta }{\mathrm{2}}\right)\mathrm{s}\mathrm{h}{m}_{k}y\right.+\\ \left.\mathrm{s}\mathrm{h}{m}_k\left(h_{\mathrm{m}}-\frac{\delta }{\mathrm{2}}\right)\frac{\mathrm{c}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}{\mathrm{s}\mathrm{h}{m}_k\frac{\delta }{\mathrm{2}}}\mathrm{c}\mathrm{h}{m}_{k}y+\mathrm{1}\right]\mathrm{s}\mathrm{i}\mathrm{n}{m}_{k}x\end{array}（26）

2.3 反电动势计算

本文所研究的DSILPMSLM采用非重叠绕组， 在上述推导出的气隙磁密分布的基础上， 推导出空载相反电动势的精确解析式. 图4为非重叠绕组的平视图.

一匝线圈的磁链表达式为：

\psi =-{\displaystyle \underset{vt+\gamma }{\overset{vt+D-\gamma }{\int }}}{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}{B}_{k\mathrm{1}}\frac{\mathrm{2}\mathrm{s}\mathrm{h}(m_k{h}_{\mathrm{w}}/\mathrm{2})}{h_{\mathrm{w}}}{L}_{\mathrm{e}\mathrm{f}}\mathrm{s}\mathrm{i}\mathrm{n}(m_{k}x)\mathrm{d}x（27）

B_{k\mathrm{1}}=\frac{\mathrm{4}{B}_{\mathrm{r}}\mathrm{s}\mathrm{i}\mathrm{n}(m_k\tau /\mathrm{2})\mathrm{s}\mathrm{i}\mathrm{n}(m_k{\tau }_{\mathrm{m}}/\mathrm{2})\mathrm{s}\mathrm{h}(m_k{h}_{\mathrm{m}})}{\tau {m_k}^{\mathrm{2}}\mathrm{s}\mathrm{h}(m_k\delta /\mathrm{2})}（28）

式中： d为单个线圈边宽度； \gamma为一匝环形线圈的宽度； D为单个线圈总宽度； L_{\mathrm{e}\mathrm{f}}为线圈直线部分长度； h_{\mathrm{w}}为线圈高度.

对单个线圈的磁链积分可得：

\psi =-{\displaystyle \sum _{k=\mathrm{1}}^{\infty }}\frac{\mathrm{4}{B}_{k\mathrm{1}}{L}_{\mathrm{e}\mathrm{f}}N}{{m_k}^{\mathrm{2}}{h}_{\mathrm{w}}d}sh\left(m_k\frac{h_{\mathrm{w}}}{\mathrm{2}}\right)sin\left(m_k\frac{d}{\mathrm{2}}\right)\times

{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}{}_{}^{}\mathrm{s}\mathrm{i}\mathrm{n}\left(m_k\frac{D-d}{\mathrm{2}}\right)\mathrm{s}\mathrm{i}\mathrm{n}\left[m_k\left(vt+\frac{D}{\mathrm{2}}\right)\right]（29）

则空载相反电动势为：

E=-k_{\mathrm{q}}\frac{Q}{\mathrm{3}a}\frac{\mathrm{d}\psi }{\mathrm{d}t}=\frac{\mathrm{4}{k}_{\mathrm{q}}Q{L}_{\mathrm{e}\mathrm{f}}Nv}{\mathrm{3}a{h}_{\mathrm{w}}d}\times

{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}\frac{B_{k\mathrm{1}}}{m_k}\mathrm{s}\mathrm{h}\left(m_k\frac{h_{\mathrm{w}}}{\mathrm{2}}\right)\mathrm{s}\mathrm{h}\left(m_k\frac{d}{\mathrm{2}}\right)\mathrm{s}\mathrm{h}\left(m_k\frac{D-d}{\mathrm{2}}\right)\times

\mathrm{c}\mathrm{o}\mathrm{s}\left[m_k\left(vt+\frac{D}{\mathrm{2}}\right)\right]（30）

式中：k_{\mathrm{q}}为分布因数； Q为槽数； a为并联支路数；v为动子直线运动速度；t为时间.

2.4 推力建模分析

本文电机的槽极配合是8极6槽，单元电机为4极3槽，且绕组结构为非重叠绕组，因此只需对非叠式4极3槽结构的DSILPMSLM进行推力建模.图5为推力求解模型.

当电机绕组内部通入三相正弦交流电时， A相绕组的电流密度分布函数为：

（31）

式中A相电流密度J_{\mathrm{a}\mathrm{0}}的表达式为：

J_{\mathrm{a}\mathrm{0}}=\frac{\sqrt[]{\mathrm{2}}NI\mathrm{s}\mathrm{i}\mathrm{n}(\theta )}{h_w{t}_w}（32）

式中：N为匝数； I为电流有效值； \theta为电流初始相位角； t_{\mathrm{w}}为线圈有效宽度.

A相绕组产生的推力的微分方程可以表示为：

\mathrm{d}{F}_{\mathrm{A}}=J_{\mathrm{a}}{B}_{y}l\mathrm{d}x\mathrm{d}y（33）

式中：l为电机沿Z方向的有效长度.

对式（33）积分可得A相绕组推力， 通过对称性可知BC两相推力， 因此三相绕组的推力分别为：

\left\{\begin{array}{l}{F}_{\mathrm{A}}={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}-\mathrm{4}K\mathrm{s}\mathrm{i}\mathrm{n}\frac{m_k{t}_{\mathrm{c}}}{\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}\frac{m_k{t}_{\mathrm{w}}}{\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}\theta \mathrm{c}\mathrm{o}\mathrm{s}{m}_k(x_{\mathrm{0}}+t_{\mathrm{c}})\\ F_{\mathrm{B}}={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}-\mathrm{4}K\mathrm{s}\mathrm{i}\mathrm{n}\frac{m_k{t}_{\mathrm{c}}}{\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}\frac{m_k{t}_{\mathrm{w}}}{\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}\left(\theta -\frac{\mathrm{2}\mathrm{\pi }}{\mathrm{3}}\right)\times \\ \mathrm{c}\mathrm{o}\mathrm{s}{m}_k(x_{\mathrm{0}}+t_{\mathrm{c}}+t_{\mathrm{c}\mathrm{d}})\\ F_{\mathrm{C}}={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}-\mathrm{4}K\mathrm{s}\mathrm{i}\mathrm{n}\frac{m_k{t}_{\mathrm{c}}}{\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}\frac{m_k{t}_{\mathrm{w}}}{\mathrm{2}}\mathrm{s}\mathrm{i}\mathrm{n}\left(\theta +\frac{\mathrm{2}\mathrm{\pi }}{\mathrm{3}}\right)\times \\ \mathrm{c}\mathrm{o}\mathrm{s}{m}_k(x_{\mathrm{0}}+t_{\mathrm{c}}+\mathrm{2}{t}_{\mathrm{c}\mathrm{d}})\end{array}\right.（34）

K=\frac{\mathrm{2}\sqrt[]{\mathrm{2}}lNI{\mu }_{\mathrm{0}}{M}_k\mathrm{s}\mathrm{h}\left(m_k{h}_{\mathrm{m}}\right)}{{m_k}^{\mathrm{2}}\mathrm{s}\mathrm{h}(m_k\delta /\mathrm{2})h_{\mathrm{w}}{t}_{\mathrm{w}}}\mathrm{s}\mathrm{h}(m_k{h}_{\mathrm{w}}/\mathrm{2})（35）

将ABC三相推力相加即可得到电机推力.

3 有限元仿真分析

本文研究的8极6槽无铁芯永磁同步直线电机的主要结构参数见表1.

有限元仿真分析非重叠绕组无铁芯永磁同步直线电机材料分配见表2.

根据上述结构参数和材料分配，进行有限元仿真计算， 通过Ansys Maxwell仿真软件建立了该款电机的三维有限元模型（充分考虑磁饱和状态），如图6所示.

电机同步速度为1.524 m/s，通入绕组中的三相电流的表达式为：

\left\{\begin{array}{l}{I}_{\mathrm{A}}=\mathrm{3.6}\times \sqrt[]{\mathrm{2}}\times \mathrm{s}\mathrm{i}\mathrm{n}(\mathrm{T}\mathrm{i}\mathrm{m}\mathrm{e}\times \mathrm{100}\mathrm{\pi })\\ I_{\mathrm{B}}=\mathrm{3.6}\times \sqrt[]{\mathrm{2}}\times \mathrm{s}\mathrm{i}\mathrm{n}\left(\mathrm{T}\mathrm{i}\mathrm{m}\mathrm{e}\times \mathrm{100}\mathrm{\pi }-\mathrm{2}\mathrm{\pi }/\mathrm{3}\right)\\ I_{\mathrm{C}}=\mathrm{3.6}\times \sqrt[]{\mathrm{2}}\times \mathrm{s}\mathrm{i}\mathrm{n}\left(\mathrm{T}\mathrm{i}\mathrm{m}\mathrm{e}\times \mathrm{100}\mathrm{\pi }+\mathrm{2}\mathrm{\pi }/\mathrm{3}\right)\end{array}\right.（36）

电机磁密云图分布如图7所示，可以看出电机各部分未呈现红色和橙色部分，说明电机未达到磁饱和状态.

反电动势和电机推力是无铁芯永磁同步直线电机中两个重要的性能参数，反电动势波形（正弦性）直接反映气隙磁场的合理性，推力是电机能力最直接的体现.因此有限元仿真主要分析电机反电动势和平均推力.通过仿真计算，可以得到电机三相反电动势和推力随时间变化的波形.

三相反电动势波形如图8所示，可以看出，波形具有明显的对称性和周期性，因此只需单独对A相进行快速傅里叶变换（fast Fourier transform）， 分析谐波影响.

对A相反电动势进行傅里叶分析，得到的反电动势谐波含量如图9所示，不难看出，基波幅值在18.8 V左右，理论上DSILPMSLM的对称磁场不产生偶次谐波，产生的原因可能是永磁体充磁不均匀或绕组不对称引起的三相电感不平衡，所以仅考虑奇数次谐波的影响，三次谐波电压幅值在0.2 V左右， 占比不到2%，五次谐波、九次谐波含量更低，总谐波失真（total harmonic distortion，THD）为2.62%，因此基本上可以忽略高次谐波对此款电机性能的影响.

推力与时间的关系如图10所示，电机的平均推力稳定在95 N左右，采用取点求平均值的方法得到电机平均推力为94.83 N，整个时间段内，波动范围非常小，最大推力为97.34 N，最小推力为93.29 N， 相差为4.05 N，推力波动仅约为4%.

4 解析法与有限元法的对比分析

为了验证解析法建立的数学模型能否准确描述电机内部磁场情况， 分别对气隙区域磁密、永磁体区域磁密、空载反电动势以及电机推力开展解析计算和有限元仿真并进行对比分析.

对于气隙处的磁密而言， 本文分析气隙中心（y=0 mm）处的磁密分量B_y随轴向位置变化的关系. 0~15 mm范围内气隙磁密随轴向位置变化的解析值如表3所示. 如图11所示， 气隙中心处的法向磁密呈现良好的正弦分布， 可以看出解析法计算结果和有限元法所得波形基本一致.

接下来分析永磁体区域中心磁密法向分量B_y随轴向位置变化的波形， 因为电机关于x轴对称， 所以只需分析一侧永磁体的磁密分布情况. 图12为y=5.9 mm处（永磁体区域2）磁密分量B_y随轴向位置变化的波形， 可以看出永磁体区域波形接近于方波， 磁密的幅值为1.22 T， 波形的周期大约为一对极距， 解析法计算结果与有限元法所得波形基本一致.

结合图11和图12， 不难看出， 越靠近永磁体区域， 磁密幅值越大， 高次谐波的含量也越大， 气隙的正弦性越差， 波形越接近于方波.

速度为1.524 m/s时， A相空载反电动势如图13所示， 感应电动势的幅值大概为18 V， 可以看出解析法计算结果和有限元法所得波形基本一致， 二者之间存在误差的原因可能是解析法假设磁场沿轴向均匀分布而忽略了端部的磁场衰减.

对于推力而言， 首先利用式（34）分别计算ABC三相推力， 然后用不同电流计算出的平均推力来验证解析法的准确性. 以A相为例， 对一个周期内的时间采样进行计算，表4列出0~10 ms内解析法的计算结果.

ABC三相推力解析法的计算结果和有限元法所得波形如图14所示， 可以看出， 两种方法所得到的波形趋势基本一致.

对采样电流（0~10 A）下解析法和有限元法得到的平均推力进行对比分析，依据式（37）对有限元法所得推力求取平均值F_av，N为采样时间点的个数，F_N为对应推力.

F_{\mathrm{a}\mathrm{v}}=\frac{F_{\mathrm{1}}+F_{\mathrm{2}}+\cdot \cdot \cdot +F_N}{N}（37）

不同电流下解析法和有限元法所得平均推力如图15所示.

从上述分析可知，对于不同物理量，解析法和有限元法所得的结果基本一致，但两种方法之间还是稍有偏差，存在误差的主要原因是解析法计算时需要对磁场做路径积分，端部磁场的衰减会将积分误差放大，有限元法的计算结果与网格剖分的大小、边界条件的确定存在密切关系，同时有限元法忽略了电磁热耦合，这可能导致有限元法求解结果存在误差.

5 实验验证

空载反电动势的大小、正弦性畸变率与空载气隙磁密分布密切相关， 同时空载相反电动势波形正弦性畸变率直接反映DSILPMSLM的推力波动水平， 电机平均推力是电机能力的直接体现. 因此选择空载反电动势和平均推力为研究对象， 将解析法、有限元法和实验测试得到的波形进行对比， 来验证解析法的准确性.

实验测试平台如图16所示， 将初级三相绕组 U/V/W和示波器的三个信号输入口相连， 通过匀速推动电机（实验测得速度为0.34 m/s）， 可以得到电机反电动势的变化情况. 伺服驱动器通过调试软件ServoStudio给电机施加电流激励， 并用测力计记录电机推力变化情况.

由于三相对称， 选取A相为研究对象， 有限元法、解析法和实验所测得A相空载反电动势波形如图17所示.

从图17可以看出， 本文提出的解析法和有限元法、实验所测得波形基本吻合， 且在仿真速度和实验速度下解析法均能保持良好的一致性. 样机实际反电动势小于仿真值， 出现误差的原因可能是线圈通过绕线机绕制后出现松散的现象， 导致实际绕组因数小于理论值.此外，动子安装时可能出现的上下气隙不对称也会导致误差增加.

不同电流下三种方法得到的平均推力如图18所示. 可以看出实际所测推力平均值依然小于仿真值， 产生原因可能是加工导致实际绕组因数偏小和气隙不对称.

6 结 论

1）本文采用的解析法能准确反映电机磁场分布特性， 其计算结果与有限元仿真及实验测量的高度一致性表明该建模方法具有普适性， 适用于同类型无铁芯永磁同步直线电机的磁场计算.

2）相较于传统有限元法， 本文所提解析法在保证精度的同时显著提升了计算效率， 与已有解析研究相比， 通过引入等效磁化电流模型更完整地考虑了永磁体空间谐波的影响， 该方法对高性能永磁同步直线电机的设计与优化具有重要的指导意义.

参考文献

[1] 卢琴芬， 沈燚明， 叶云岳 ．永磁直线电动机结构及研究发展综述［J］．中国电机工程学报， 2019（9）： 2575-2588．

[2] LU Q F， SHEN Y M， YE Y Y ．Development of permanent magnet linear synchronous motors structure and research［J］. Proceedings of the Chinese Society for Electrical Engineering， 2019（9）：2575-2588．（in Chinese）

[3] YANG J T， CHEN J H， CAI Z Y， et al ．Design and optimization of a novel permanent magnet linear motor with the same number of poles and slots［J］．IEEE Transactions on Industry Applications， 2025， 61（1）：16-24.

[4] ZHANG Y Q， YANG Z L， YU M H， et al ．Analysis and design of double-sided air core linear servo motor with trapezoidal permanent magnets［J］．IEEE Transactions on Magnetics， 2011， 47（10）：3 236-3239.

[5] 徐鹏， 郭铖， 苏鑫， 等 ．永磁同步电机模糊积分终端滑模控制［J］．电工电能新技术， 2025， 44（4）：52-59．

[6] XU P， GUO C， SU X， et al ．Sliding mode control of fuzzy integral terminal of PMSM［J］. Advanced Technology of Electrical Engineering and Energy， 2025， 44（4）：52-59．（in Chinese）

[7] 王伟华， 张子娇， 周海波， 等. 一种同步直线电机推力波动特性的检测方法［J］. 中国电机工程学报， 2016， 36（2）：540-546．

[8] WANG W H， ZHANG Z J， ZHOU H B， et al. A method for detecting the thrust fluctuation characteristics of synchronous linear motors［J］. Proceedings of the Chinese Society for Electrical Engineering， 2016， 36（2）：540-546．（in Chinese）

[9] MENG S， YANG X Y， ZHAO S W， et al ．Development of the equivalent magnetic circuit model for a surface-interior permanent magnet synchronous motor［C］//2015 6th International Conference on Power Electronics Systems and Applications （PESA）， December 15-17， 2015， Hong Kong， China. IEEE， 2016： 1-4.

[10] SI J K， HE S， CAO W P， et al ．Electromagnetic characteristics analysis of surface-mounted and interior hybrid PMSM based on equivalent magnetic circuit method［C］//2014 17th International Conference on Electrical Machines and Systems （ICEMS）．October 22-25， 2014， Hangzhou， China．IEEE， 2015：1125-1131.

[11] LEE K D， LEE J， LEE H W. Inductance calculation of flux concentrating permanent magnet motor through nonlinear magnetic equivalent circuit［J］. IEEE Transactions on Magnetics， 2015， 51（11）： 8204304.

[12] CUI F R， SUN Z L， XU W， et al ．Optimization analysis of long primary permanent magnet linear synchronous motor［J］. IEEE Transactions on Applied Superconductivity， 2021， 31（8）：1-4.

[13] 王淑红， 熊光煜 ．永磁直线同步电机解析分析及有限元验证［J］．电机与控制学报， 2011， 15（7）： 8-13．

[14] WANG S H， XIONG G Y ．Analyses of permanent magnet linear synchronous motor and validated by FEM［J］．Electric Machines and Control， 2011， 15（7）： 8-13．（in Chinese）

[15] VAEZ-ZADEH S， ISFAHANI A H ．Enhanced modeling of linear permanent-magnet synchronous motors［J］. IEEE Transactions on Magnetics， 2007， 43（1）：33-39.

[16] RUI H， JIAN F ．Magnetic field analysis of permanent magnet array for planar motor based on Equivalent Magnetic Charge method［C］//Proceedings of the 10th World Congress on Intelligent Control and Automation.

[17] ZHANG Z J， LUO M Z， KOU B Q， et al. Design and analysis of a trilateral permanent magnet linear synchronous motor with slotless ring windings for transport systems［C］//2016 IEEE Vehicle Power and Propulsion Conference （VPPC）， October 17-20， 2016， Hangzhou， China． IEEE， 2016：1-6.

[18] SHEIKH-GHALAVAND B， VAEZ-ZADEH S， ISFAHANI A H ．An improved magnetic equivalent circuit model for iron-core linear permanent-magnet synchronous motors［J］. IEEE Transactions on Magnetics， 2010， 46（1）：112-120.

[19] WANG F R， LIAO Y Y， CHEN J H， et al ．Analytical calculation of air gap magnetic field in permanent magnet linear motors［C］//2017 12th IEEE Conference on Industrial Electronics and Applications （ICIEA）， June 18-20， 2017， Siem Reap， Cambodia． IEEE， 2018：1706-1711.

[20] 王淑红， 熊光煜 ．新型筒型永磁动圈式直线电动机气隙磁场解析分析［J］．电工技术学报， 2007， 22（5）： 40-44．

[21] WANG S H， XIONG G Y. Analysis of air-gap field in a new tubular linear PM motor with moving coil［J］. Transactions of China Electrotechnical Society， 2007， 22（5）：40-44．（in Chinese）

[22] 赵镜红， 张晓锋， 张俊洪， 等 ．圆筒永磁直线同步电机磁场和推力分析［J］．电机与控制学报， 2010， 14（1）： 12-17．

[23] ZHAO J H， ZHANG X F， ZHANG J H， et al. Field and thrust analysis of tubular permanent magnet linear synchronous motor［J］. Electric Machines and Control， 2010， 14（1）：12-17．（in Chinese）

[24] 赵镜红， 张俊洪， 方芳， 等. 径向充磁圆筒永磁直线同步电机磁场和推力解析计算［J］. 电工技术学报， 2011， 26（7）：154-160．

[25] ZHAO J H， ZHANG J H， FANG F， et al. Analytical calculation of magnetic field and thrust force in radial magnetized tubular permanent magnet linear synchronous motor［J］. Transactions of China Electrotechnical Society， 2011， 26（7）：154-160．（in Chinese）

[26] 宋玉晶， 张鸣， 朱煜 ．Halbach永磁阵列磁场解析求解及推力建模［J］．电工技术学报， 2014， 29（11）：61-67．

[27] SONG Y J， ZHANG M， ZHU Y ．Modeling of Halbach array magnetic field and motor thrust［J］. Transactions of China Electrotechnical Society， 2014， 29（11）： 61-67．（in Chinese）

[28] 李斌， 董川， 金世鑫 ．基于子域叠加模型的表贴式永磁电机电枢反应磁场分析［J］．电工电能新技术， 2017， 36（1）：10-15．

[29] LI B， DONG C， JIN S X ．Subdomain superposition model of armature reaction magnetic field of surface mounted permanent-magnent machine［J］．Advanced Technology of Electrical Engineering and Energy， 2017， 36（1）：10-15．（in Chinese）

[30] GYSEN B L J， MEESSEN K J， PAULIDES J J H， et al ．3-D analytical and numerical modeling of tubular actuators with skewed permanent magnets［J］. IEEE Transactions on Magnetics， 2011， 47（9）：2200-2212.

[31] 孙少男， 李岩， 张炳义， 等 ．集中绕组模块组合式定子结构对低速大功率永磁电机性能影响［J］．电工电能新技术， 2024， 43（6）：11-24．

[32] SUN S N， LI Y， ZHANG B Y， et al.Influence of modular combined stator structure with centralized winding on the performance of low speed high power permanent magnet motor［J］．Advanced Technology of Electrical Engineering and Energy， 2024， 43（6）：11-24．（in Chinese）

[33] MEESSEN K J， PAULIDES J J H， LOMONOVA E A. Force calculations in 3-D cylindrical structures using Fourier analysis and the Maxwell stress tensor［J］. IEEE Transactions on Magnetics， 2013， 49（1）：536-545.

[34] 王婷婷， 周睿， 郑常宝， 等 ．基于转矩解析模型和遗传算法的球形电机电流优化计算［J］．电工电能新技术， 2022， 41（2）：42-52．

[35] WANG T T， ZHOU R， ZHENG C B， et al. Optimal calculation of spherical motor current based on torque analytical model and genetic algorithm［J］. Advanced Technology of Electrical Engineering and Energy， 2022， 41（2）：42-52．（in Chinese）

[36] SHIN K H， LEE S H， CHO H W， et al ．Analysis on the pitching moment in permanent-magnet linear synchronous motor for linear motion stage systems［J］. IEEE Transactions on Magnetics， 2016， 52（7）：1-4.

基金资助

湖南省科技计划资助项目(2026JK2019)

Hunan Provincial Science and Technology Program(2026JK2019)