Fractional-Order Sliding Mode Control for Maglev Rotary Table Based on Disturbance Compensation
-
摘要:
针对存在非线性、耦合性和不确定性的磁悬浮转台的高精度运动控制问题,提出一种基于非线性干扰观测器的分数阶滑模控制方法以提高跟踪精度. 首先,基于系统电磁力模型和动态解耦方法,构建六自由度磁悬浮转台系统动力学模型;其次,设计非线性干扰观测器,对包含系统误差、六自由度间耦合项和外界干扰的集总扰动进行估计,证明了估计误差有界且可调节到任意小;然后,在离散域提出了一种分数阶滑模面,采用分数幂函数替代传统符号函数来抑制抖振,引入分数阶微积分来减小跟踪误差;最后,设计有限时间收敛的分数阶滑模控制策略,并利用李雅普诺夫稳定性理论证明闭环系统稳定性. 实验结果表明:与整数阶滑模控制方法相比,采用所提方法,2个水平自由度和绕竖直方向旋转自由度对三角波的跟踪误差均方根分别减小了12.8%、16.8%和23.7%,最大跟踪误差分别减小9.26%、13.00%和33.20%;跟踪圆形轨迹时,2个水平自由度的跟踪误差均方值分别减小6.39%和12.40%,最大跟踪误差分别减小9.90%和12.10%.
Abstract:In view of the high-precision motion control problem of the maglev rotary table with nonlinearity, coupling, and uncertainty, a fractional-order sliding mode control method based on a nonlinear disturbance observer was proposed to improve the tracking accuracy. Firstly, based on the electromagnetic force model of the system and the dynamic decoupling method, the dynamical model of the six-degree-of-freedom maglev rotary table system was constructed. Secondly, a nonlinear disturbance observer was designed to estimate the lumped disturbance including system error, coupling term between six degrees of freedom, and external interference. It was proved that the estimation error was bounded and could be made arbitrarily small. Then, a fractional-order sliding surface was proposed in the discrete domain, where the fractional power function was used instead of the traditional symbolic function to suppress jitter, and the fractional calculus was introduced to reduce the tracking error. Finally, a fractional-order sliding mode control strategy with finite time convergence was designed, and the stability of the closed-loop system was proved by Lyapunov stability theory. The experimental results reveal that compared to the integer-order sliding mode control method, the proposed method reduces the root mean square of tracking error for triangular waves by 12.8%, 16.8%, and 23.7% for the two horizontal degrees of freedom and the rotational degree about the vertical axis, respectively, while the maximum tracking errors are reduced by 9.26%, 13.00%, and 33.20% respectively. When tracking a circular trajectory, the mean square values of tracking errors for two horizontal degrees of freedom are decreased by 6.39% and 12.40%, and the maximum tracking errors are reduced by 9.90% and 12.10%, respectively.
-
平面电机式磁悬浮转台具有多自由度平移、旋转运动能力[1],且克服了传统机械传动平台存在的接触摩擦、回程误差等限制因素[2],在微细加工、复杂微结构成型等精密制造应用场景中具有极大优势. 然而,由于采用Halbach永磁阵列,磁悬浮转台气隙磁场的空间周期性使得输出磁力大小和方向随运动次级的位置变化而变化,且在不同运动自由度之间存在非线性的磁力耦合. 此外,磁悬浮系统固有的开环不稳定性使得工作台易受扰动影响,对运动控制算法的稳定性和抗扰性要求较高.
为实现磁悬浮工作台对指定轨迹的快速精准跟踪,研究人员基于PID(proportion integration differentiation)控制器对其进行控制器设计[3-4]. 但PID控制难以兼顾动态性能和运动精度,如增大微分增益会导致在提升响应速度的同时放大测量噪声,对进一步提高跟踪精度具有局限性. 此外,建模和测量的不精确以及外界扰动导致系统存在不确定性,易于引起控制系统品质恶化,从而导致运动精度变低. 滑模控制在系统存在参数摄动和不确定性扰动时表现优异[5],魏静波等[6]结合终端滑模干扰观测器和全局快速终端滑模控制优化了磁悬浮球系统的控制效果. 此外,采用分数幂函数代替原有的符号函数,利用非线性部分辅助线性部分削弱切换项带来的抖振问题[7-8]. 近年来,具有全局记忆性的分数阶微积分被引入滑模控制策略[9],相较于传统的整数阶滑模控制更具优越性. Sun等[10]针对直线电机提出一种离散时间分数阶终端滑模控制策略,提升了其运动跟踪性能;Kuang等[11]利用交叉耦合二阶分数阶滑模控制方法提高了直线电机的同步控制和跟踪精度,对磁悬浮转台系统运动控制方法设计具有借鉴意义. 此外,相比于直接针对系统的非线性模型进行控制器设计[12],采用动态解耦方法将六自由度运动模型转换为6个二阶线性模型可简化控制器设计过程[13-14]. 为充分考虑系统的耦合效应,采用干扰观测器进行扰动补偿[15].
为抑制系统非线性、强耦合性和不确定性对磁悬浮转台运动精度的影响,首先,采用非线性干扰观测器对包含模型误差、测量误差、六自由度间耦合项和其他外界干扰在内的系统集总扰动进行补偿;然后,设计了一种离散时间分数阶滑模控制策略,分数幂函数削弱抖振,分数阶微积分对跟踪误差的全局记忆性来提升系统的动态性能和轨迹跟踪性能;最后,通过实验验证了所提方法的优越性和实用性.
1. 磁悬浮转台动力学分析
图1中的磁悬浮转台具有中心对称结构,主体部分由8个矩形线圈绕组(定子,Coil 1~Coil 8)和1个搭载圆周型Halbach永磁阵列的圆台(动子)构成,适合执行平面运动和绕竖直方向的大范围旋转运动. 以线圈阵列上表面的几何中心为坐标原点O,通过设计控制器调节线圈通电电流,可使线圈长边区域的电流与Halbach永磁阵列产生的空间磁场相互作用,产生可控的磁力和磁力矩,从而驱动磁悬浮转台动子在x、y、z、α、β、γ方向上进行六自由度运动.
1.1 电磁力模型
圆周型Halbach永磁阵列由12个磁极组成,其中每个磁极包含4块永磁体. 在动子坐标系{m}中,以Halbach永磁阵列下表面的几何中心为原点mO,将该Halbach永磁阵列沿其切向方向展开,则其下方空间的磁通密度分布如图2所示,图中:mr、mθ、mz分别为坐标系{m}的3个坐标变量,mθ∈[0,2π );wm、hm分别为永磁体宽和高;ϕm为相邻永磁体间夹角; r为永磁阵列中心半径;τ为永磁阵列周期角.
忽略高次谐波和边端效应,Halbach永磁阵列下方空间的磁通密度mB可表示为[16]
mB=[mBrmBθmBz]=[0−4Brπsin(πwmτ)(1−e−wm/λc)sin(πmθ2ϕm)emz/λc4Brπsin(πwmτ)(1−e−wm/λc)cos(πmθ2ϕm)emz/λc], (1) 式中:mBr、mBθ、mBz分别为mB在mr、mθ、mz方向上的分量,Br为永磁体剩磁;λc=τ/2π ,为特征波长.
根据洛伦兹积分公式和牛顿第三定律,动子受到的驱动力为通电线圈所受电磁力的反作用力,即
{Fi=−R∫Vcoil iJi×mBdν,Ti=−R∫Vcoil iri×Ji×mBdν, (2) 式中: Fi,Ti∈R3,分别为第i相线圈提供的磁力和磁力矩,i = 1,2,…,8,R为不同坐标系间的转换矩阵;Vcoil i为第i相线圈所占的空间,为线圈的有效作用区域;dν=mrdmrdmθdmz;Ji为线圈体积元内的电流密度,方向为通入的电流方向,大小为|Ji|=NIirchc, N为线圈匝数,Ii为驱动电流,rc为线圈侧边宽度,hc为线圈厚度;ri为线圈体积元到点mO的力臂. 此外,线圈的长边宽度和短边宽度分别为lc和wc.
可采用磁节点和高斯积分方法对式(2)进行计算[1]. 将Fi、Ti与其驱动电流Ii的关系式简记为
{Fi=[Fi,xFi,yFi,z]T=[Γi,xΓi,yΓi,z]TIi,Ti=[Ti,αTi,βTi,γ]T=[Γi,αΓi,βΓi,γ]TIi, (3) 式中:Fi,x,Fi,y,Fi,z∈R,分别为第i相线圈提供的磁力在平动自由度x、y、z的分量,Γi,x,Γi,y,Γi,z∈R,为相应的磁力系数;Ti,α,Ti,β,Ti,γ∈R,分别为第i相线圈提供的磁力矩在转动自由度α、β、γ的分量,Γi,α,Γi,β,Γi,γ∈R,为对应的磁力矩系数,与动子的六自由度位移有关.
综合8相线圈提供的磁力和磁力矩,可得磁悬浮转台动子所受到的合驱动力F和力矩T为
w=[FT]=[Γ1,xΓ2,x⋯Γ8,xΓ1,yΓ2,y⋯Γ8,yΓ1,zΓ2,z⋯Γ8,zΓ1,αΓ2,α⋯Γ8,αΓ1,βΓ2,β⋯Γ8,βΓ1,γΓ2,γ⋯Γ8,γ][I1I2⋮I8]=ΓI, (4) 式中:w∈R6,为由动子所受的合磁力和磁力矩构成的向量;Γ∈R6×8,为系统的电流-动力传递矩阵;I∈R8,为通入各相线圈的电流构成的向量.
1.2 系统动力学
在磁悬浮转台系统设计中,通过采用合理的电流分配方法,可以实现不同运动自由度上的力和力矩的动态解耦[13],如式(5).
I=Γ†w=ΓT(ΓΓT)−1w, (5) 式中:Γ†为Γ的伪逆.
记磁悬浮转台相对于初始状态的六自由度位移x=[xmymzmαmβmγm]T,对动子进行动力学分析可得
¨x=J−1(w−[00mg000]T)+d=J−1u+d, (6) 式中:J=diag (m, m, m, Jα, Jβ, Jγ),为动子的惯性矩阵,m为动子质量, Jα, Jβ, Jγ∈R,分别为转动自由度α、β、γ的转动惯量;u=w−[00mg000]T;d∈R6,为系统集总扰动,包括模型误差、测量误差、六自由度间耦合项和其他外界干扰.
因此,可将u作为控制量设计闭环控制算法,根据式(5)得到电流向量I,从而调节线圈电流大小以实现磁悬浮转台运动控制.
假设1 扰动d缓慢变化,即||˙d||2⩽, \varepsilon \approx 0 ,其中||·||2表示二范数, \varepsilon 为{\boldsymbol{d}}的二范数的上界.
由于系统实际工作环境较为稳定且能量有限,且工作台稳态运动时各自由度间的耦合项也随时间缓慢变化,故假设1合理.
设 {{\boldsymbol{x}}_{\text{d}}} \in {\mathbb{R}^6} ,为磁悬浮转台的期望轨迹, {\boldsymbol{e}} = {\boldsymbol{x}} - {{\boldsymbol{x}}_{\text{d}}} = [ \begin{array}{*{20}{c}} {{e_1}}&{{e_2}}& \cdots &{{e_n}} \end{array} ]{{\text{ }}^{\text{T}}} ,为跟踪误差,则跟踪系统的动力学方程为
{\boldsymbol{\ddot e}} = {{\boldsymbol{J}}^{ - 1}}{\boldsymbol{u}} + {\boldsymbol{d}} - {{\boldsymbol{\ddot x}}_{\text{d}}}. (7) 采用前向欧拉法将式(7)离散化,可得跟踪误差的离散域状态方程为
{\Delta ^2}{\boldsymbol{e}}(k) = {{\boldsymbol{J}}^{ - 1}}{\boldsymbol{u}}(k) + {\boldsymbol{d}}(k) - {\Delta ^2}{{\boldsymbol{x}}_{\text{d}}}(k), (8) 式中: \Delta 为前向差分算子; \Delta {{\boldsymbol{x}}_{\text{d}}}(k) = \dfrac{{{{\boldsymbol{x}}_{\text{d}}}(k + 1) - {{\boldsymbol{x}}_{\text{d}}}(k)}}{h} , {\Delta }^{{n}_{ + } + 1}{{\boldsymbol{x}}}_{\text{d}}(k)=\Delta ({\Delta }^{{n}_{ + }}{{\boldsymbol{x}}}_{\text{d}}(k)),{n}_{ + }\in {\mathbb{N}}^{ + } , {{\boldsymbol{x}}_{\text{d}}}(k) = {{\boldsymbol{x}}_{\text{d}}}(kh) ,h \in {\mathbb{R}^ + },为采样间隔,且满足 0 < h \ll 1{\text{ s}} , k \in \mathbb{N} ,为采样步数.
2. 分数阶滑模控制器设计
2.1 非线性干扰观测器
采用非线性干扰观测器估计系统集总扰动{\boldsymbol{d}},以实现干扰补偿. 设计非线性干扰观测器为[15]
\left\{ \begin{array}{l} {\boldsymbol{\hat d}} = {\boldsymbol{r}} + {\boldsymbol{p}}({\boldsymbol{e}}), \\ {\boldsymbol{\dot r}} = - {\text{ }}{{\boldsymbol{L}}_0}({\boldsymbol{r}} + {\boldsymbol{p}}({\boldsymbol{e}}) + {{\boldsymbol{J}}^{ - 1}}{\boldsymbol{u}} - {{{\boldsymbol{\ddot x}}}_{\text{d}}}), \\ \end{array} \right. (9) 式中: {\boldsymbol{\hat d}} \in {\mathbb{R}^6} ,为扰动估计值; {\boldsymbol{r}},\dot{{\boldsymbol{r}}}\in {\mathbb{R}}^{6} ,为干扰观测器状态; {{\boldsymbol{L}}_0} \in {\mathbb{R}^{6 \times 6}} ,为观测器增益,是一个正定的常数矩阵;{\boldsymbol{p}}( {\text{•}} )为非线性函数,且\dot {\boldsymbol{p}}({\boldsymbol{e}}) = {{\boldsymbol{L}}_0}{\boldsymbol{\ddot e}}.
由式(9)可知,扰动估计值的变化率为
\begin{split} & {{\dot {\hat {\boldsymbol d}}}} = {{\dot {\boldsymbol r}}} + \dot p({\boldsymbol{e}}) =- {{\boldsymbol{L}}_0}({\boldsymbol{r}} + p({\boldsymbol{e}}) + {{\boldsymbol{J}}^{ - 1}}{\boldsymbol{u}} - {{{{\ddot {\boldsymbol x}}}}_{\text{d}}}) + \\ &\quad {{\boldsymbol{L}}_0}{\boldsymbol{\ddot e}}= {{\boldsymbol{L}}_0}[({{\ddot {\boldsymbol e}}} - {{\boldsymbol{J}}^{ - 1}}{\boldsymbol{u}} + {{{{\ddot {\boldsymbol x}}}}_{\text{d}}}) - \\ &\quad ({\boldsymbol{r}} + p({\boldsymbol{e}}))] = {{\boldsymbol{L}}_0}({\boldsymbol{d}} - {{\hat {\boldsymbol d}}}) = {{\boldsymbol{L}}_0}{{\tilde {\boldsymbol d}}}, \end{split} (10) 式中: {\boldsymbol{\tilde d}} = {\boldsymbol{d}} - {\boldsymbol{\hat d}} ,为扰动估计误差,其变化率如式(11).
{{\dot {\tilde {\boldsymbol d}}}} = {\boldsymbol{\dot d}} - {\boldsymbol{\dot {\hat d}}} = {\boldsymbol{\dot d}} - {{\boldsymbol{L}}_0}{\boldsymbol{\tilde d}}. (11) 定理1 在假设1的条件下,利用式(9)估计{\boldsymbol{d}},可以使估计误差 {\boldsymbol{\tilde d}} 有界,且其上界可调节到任意小.
证明 构造Lyapunov函数为
U = \frac{1}{2}{{\boldsymbol{\tilde d}}^{\text{T}}}{\boldsymbol{\tilde d}}. (12) 设 {{\boldsymbol{L}}_0} 的最小特征值为 {\lambda _{\min }}({{\boldsymbol{L}}_0}) , \sigma 为任意小正常数且满足 \varsigma = {\lambda _{\min }}({{\boldsymbol{L}}_0}) - \sigma > 0 ,使得杨氏不等式 {{\boldsymbol{\tilde d}}^{\text{T}}}{\boldsymbol{\dot d}} \leqslant \sigma {{\boldsymbol{\tilde d}}^{\text{T}}}{\boldsymbol{\tilde d}} + \dfrac{1}{{4\sigma }}{{\boldsymbol{\dot d}}^{\text{T}}}{\boldsymbol{\dot d}} 成立. 则
\begin{split} & \dot U = {{{\boldsymbol{\tilde d}}}^{\text{T}}}{\boldsymbol{\dot {\tilde d}}} = {{{\boldsymbol{\tilde d}}}^{\text{T}}}({\boldsymbol{\dot d}} - {{\boldsymbol{L}}_0}{\boldsymbol{\tilde d}})\; \leqslant - {{{\boldsymbol{\tilde d}}}^{\text{T}}}{{\boldsymbol{L}}_0}{\boldsymbol{\tilde d}} + \sigma {{{\boldsymbol{\tilde d}}}^{\text{T}}}{\boldsymbol{\tilde d}} + \\ &\quad \frac{1}{{4\sigma }}{{{\boldsymbol{\dot d}}}^{\text{T}}}{\boldsymbol{\dot d}} \leqslant - [{\lambda _{\min }}({{\boldsymbol{L}}_0}) - \sigma ]{{{\boldsymbol{\tilde d}}}^{\text{T}}}{\boldsymbol{\tilde d}} + \frac{{{\varepsilon ^2}}}{{4\sigma }} =\\ &\quad - 2\varsigma U + \frac{{{\varepsilon ^2}}}{{4\sigma }}, \end{split} (13) 解得
U(t) \leqslant \frac{{{\varepsilon ^2}}}{{8\sigma \varsigma }} + \left[U(0) - \frac{{{\varepsilon ^2}}}{{8\sigma \varsigma }}\right]{{\text{e}}^{-2\varsigma t}}. (14) 由式(12)可得
||{\boldsymbol{\tilde d}}|{|_2} \leqslant \sqrt {\frac{{{\varepsilon ^2}}}{{4\sigma \varsigma }} + \left[2U(0) - \frac{{{\varepsilon ^2}}}{{4\sigma \varsigma }}\right]{{\text{e}}^{ -2\varsigma t}}} . (15) 因此, ||{\boldsymbol{\tilde d}}|{|_2} 最终将收敛到以 \sqrt {{\varepsilon ^2}/4\sigma \varsigma } 为半径的球域内,即所设计的非线性干扰观测器对系统集总扰动{\boldsymbol{d}}的估计误差 {\boldsymbol{\tilde d}} 是有界的. 通过适当地选择 {{\boldsymbol{L}}_0} 和 \sigma ,可以使其上界任意小.
显然,采用式(10)中的干扰观测器,对于离散域中 {\boldsymbol{\tilde d}}(k) 的每个分量 {\tilde d_n}(k) ,存在常数 {d^*} \in {\mathbb{R}^ + } 使得 {\tilde d_n}(k) 的绝对值 |{\tilde{d}}_{n}(k)|\leqslant {d}^{*},k\in \mathbb{N} ,n = 1,2,…,6,即
-{d}^{*}\leqslant {\tilde{d}}_{n}(k)\leqslant {d}^{*}. (16) 2.2 离散时间分数阶滑模控制器
定义1 非光滑分数幂函数定义为
{ [\kern-0.15em[ e(k) ]\kern-0.15em] ^q} \triangleq {\text{ }}|e(k){|^q}{{\mathrm{sgn}}} (e(k)), (17) 式中: q 为正奇数之比, 0 < q < 1 .
对于 {\boldsymbol{e}}(k) = [ {{e_1}(k)}\;\,{{e_2}(k)}\;\, \cdots \;\,{{e_n}(k)} ]^{\text{T}} \in {\mathbb{R}^n} , { [\kern-0.15em[ {\boldsymbol{e}}(k) ]\kern-0.15em] ^q} = [ {|e_1|^q{{\mathrm{sgn}}} ({e_1}(k))}\quad{|e_2|^q{{\mathrm{sgn}}} ({e_2}(k))\quad \cdots }\quad{|e_n|^q{{\mathrm{sgn}}} ({e_n}(k))} ]^{\text{T}} .
注1 若 {e_n}(k) > 0 ,则 { [\kern-0.15em[ {e_n}(k) ]\kern-0.15em] ^q}{\text{ = }}e_n^q(k) .
为削弱切换项带来的抖振问题,采用式(17)代替符号函数[7],设计整数阶离散时间切换函数为
{\boldsymbol{s}}(k) = \Delta {\boldsymbol{e}}(k) + {l_1}{\boldsymbol{e}}(k) + {l_2}{ [\kern-0.15em[ {\boldsymbol{e}}(k) ]\kern-0.15em] ^q}, (18) 式中: {\boldsymbol{s}}(k) \in {\mathbb{R}^6} ,为滑模状态, {l_1} 、 {l_2} 为控制器参数, 0 < {l_1}h < 1 , {l_2} > 0 .
定义2 G-L形式的分数阶算子定义为
{\Delta ^a}e(k) = \frac{1}{{{h^a}}}\sum\limits_{j = 0}^k {{{( - 1)}^j}\left( {\begin{array}{*{20}{c}} a \\ j \end{array}} \right)e(k - j)} , (19) 式中: a \in \mathbb{R} , 0 < a < 1 ,h \in {\mathbb{R}^ + },二项式系数 \left( {\begin{array}{*{20}{c}} a \\ j \end{array}} \right) 如式(20).
\left( {\begin{array}{*{20}{c}} a \\ j \end{array}} \right) = \left\{ {\begin{array}{*{20}{l}} {1,}\quad {j = 0,} \\ {\dfrac{{a(a - 1) \cdots (a - j + 1)}}{{j!}},}&{j = 1,2,3, \cdots } \end{array}} \right. (20) 对于 {\boldsymbol{e}}(k) \in {\mathbb{R}^n} , {\Delta ^a}{\boldsymbol{e}}(k) = [ {{\Delta ^a}{e_1}(k)}\quad{{\Delta ^a}{e_2}(k) \quad \cdots }\quad {{\Delta ^a}{e_n}(k)} ]^{\text{T}} , {\Delta ^a}{\boldsymbol{e}}(k) = \Delta ({\Delta ^{a - 1}}{\boldsymbol{e}}(k)) .
注2 在实际系统中,设 L 为累加次数,简化计算式为
{\Delta }^{a}{\boldsymbol{e}}(k)=\frac{1}{{h}^{a}}{\displaystyle \sum _{j=0}^{L}{(-1)}^{j}\left(\begin{array}{c}a\\ j\end{array}\right){\boldsymbol{e}}(k-j)\text{,}\quad L\in {\mathbb{Z}}^{ + }}. (21) 为减小跟踪误差,可以利用分数阶微积分的全局记忆性,将式(18)改进为
{\boldsymbol{s}}(k) = \Delta {\boldsymbol{e}}(k) + {l_1}{\boldsymbol{e}}(k) + {l_2}{\Delta ^{a - 1}}{ [\kern-0.15em[ {\boldsymbol{e}}(k) ]\kern-0.15em] ^q}. (22) 由式(8)可推出
\begin{split} & {\boldsymbol{s}}(k + 1) = {\boldsymbol{s}}(k) + h\Delta {\boldsymbol{s}}(k) = \\ &\quad {\boldsymbol{s}}(k) + h({\Delta ^2}{\boldsymbol{e}}(k) + {l_1}\Delta {\boldsymbol{e}}(k) + {l_2}{\Delta ^a}{ [\kern-0.15em[ {\boldsymbol{e}}(k) ]\kern-0.15em] ^q}) = \\ &\quad {\boldsymbol{s}}(k) + h({{\boldsymbol{J}}^{ - 1}}{\boldsymbol{u}}(k) + {\boldsymbol{d}}(k) - {\Delta ^2}{{\boldsymbol{x}}_{\text{d}}}(k) + \\ &\quad {l_1}\Delta {\boldsymbol{e}}(k) + {l_2}{\Delta ^a}{ [\kern-0.15em[ {\boldsymbol{e}}(k) ]\kern-0.15em] ^q}). \end{split} (23) 提出一种由指数到达律和等效控制律组成的滑模控制律. 其中,指数到达律表达式为
{{\boldsymbol{u}}_{{\text{re}}}}(k) = - {\boldsymbol{J}}({k_1}{\boldsymbol{s}}(k) + {k_2}{ [\kern-0.15em[ {\boldsymbol{s}}(k) ]\kern-0.15em] ^b}), (24) 式中: {k_1} 、 {k_2} 、 b 为控制器参数, 0 < {k_1}h < 1 , {k_2} > 0 , b 为正奇数之比, 0 < b < 1 .
令 {\boldsymbol{s}}(k + 1) - {\boldsymbol{s}}(k) = {\boldsymbol{0}} ,并用扰动估计值 {\boldsymbol{\hat d}}(k) 代替 {\boldsymbol{d}}(k) ,可得等效控制律为
{{\boldsymbol{u}}_{{\text{eq}}}}(k) = {\boldsymbol{J}}({\Delta ^2}{{\boldsymbol{x}}_{\text{d}}}(k) - {l_1}\Delta {\boldsymbol{e}}(k) - {l_2}{\Delta ^a}{ [\kern-0.15em[ {\boldsymbol{e}}(k) ]\kern-0.15em] ^q} - {\boldsymbol{\hat d}}(k)). (25) 因此,总的滑模控制律为
{\boldsymbol{u}}(k) = {{\boldsymbol{u}}_{{\text{re}}}}(k) + {{\boldsymbol{u}}_{{\text{eq}}}}(k). (26) 将式(26)代入式(23)中,有
{\boldsymbol{s}}(k + 1) = {\boldsymbol{s}}(k) - {k_1}h{\boldsymbol{s}}(k) - {k_2}h{ [\kern-0.15em[ {\boldsymbol{s}}(k) ]\kern-0.15em] ^b} + h{\boldsymbol{\tilde d}}(k). (27) 磁悬浮转台闭环系统控制框图如图3所示.
3. 稳定性分析
3.1 系统稳定性
磁悬浮转台闭环系统稳定性可由定理2构造.
定理2 对于式(8),若选择滑模状态 {\boldsymbol{s}}(k) 的动态方程为式(22),选择干扰观测器为式(9),那么可以采用式(26)描述的分数阶滑模控制律来保证以下2个条件:
1) 当系统状态远离滑模面时,每一个滑模变量 {s_n}(k) ,将在 {k^ * } 步内进入区域 \varLambda ,如式(28)、(29).
\left\{\begin{array}{l} {\boldsymbol{\varLambda}} =\Bigg({s}_{n}(k):|{s}_{n}(k)|\leqslant \delta =\Bigg.\\ \quad\left.\psi (b)\mathrm{max}\left\{{\left(\dfrac{{d}^{\ast }}{{k}_{2}}\right)}^{\frac{1}{b}},{\left(\dfrac{{k}_{2}h}{1-{k}_{1}h}\right)}^{\frac{1}{1-b}}\right\}\right)\text{,}\\ \psi (b)=1 + {b}^{\frac{b}{1-b}}-{b}^{\frac{1}{1-b}}, \end{array}\right. (28) \left\{\begin{array}{l} {k}^{\ast }=\left\lfloor \dfrac{{s}_{n}^{2}(0)-{\delta }^{2}}{{\mu }^{2}}\right\rfloor + 1\text{,}\\ \mu ={k}_{1}h\delta + ({\psi }^{b}(b)-1){d}^{\ast }h, \end{array}\right. (29) 式中: \mu 为中间变量, \delta 为 {s_n}(k) 的上界, \psi ({\text{•}}) 为自定义函数.
2) 滑模变量 {s_n}(k) 进入区域{\rm{\varLambda}}之后将始终保持在该区域中,即若 |{s_n}(k)| \leqslant \delta ,则 |{s_n}(k + 1)| \leqslant \delta .
证明 证明分为到达段和滑动段两部分.
步骤1 (到达段) 首先证明当{s_n}(k) \notin \varLambda ,即 |{s_n}(k)| > \delta 时,滑模状态将在 {k^ * } 步内进入区域{\rm{\varLambda}}.
在离散域构造Lyapunov函数为
V(k) = s_n^2(k). (30) \begin{split} & V(k + 1) - V(k) = s_n^2(k + 1) - s_n^2(k) = \\ &\quad - ({k_1}h{s_n}(k) + {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} - h{{\tilde d}_n}(k)) \cdot \\ &\quad (2s_n^2(k) - {k_1}hs_n^2(k) - {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} + h{{\tilde d}_n}(k)). \end{split} (31) 引理1[7] 若 \psi (\zeta ) = 1 + {\zeta ^{\frac{\zeta }{{1 - \zeta }}}} - {\zeta ^{\frac{1}{{1 - \zeta }}}} , 0 < \zeta < 1 ,则有 1 < \psi (\zeta ) < 2 ,且对于任意 0 \leqslant \beta \leqslant 1 有
\beta \psi (\zeta ) - {\beta ^\zeta }{\psi ^\zeta }(\zeta ) + \psi (\zeta ) - 1 \geqslant 0. (32) 当 |{s_n}(k)| > \delta 时,考虑以下2种情况:
情况1 {s_n}(k) < - \delta .
由 {s_n}(k) < - \delta = - \psi (b)\max \left\{ {{{\left( {\dfrac{{{d^ * }}}{{{k_2}}}} \right)}^{\frac{1}{b}}},{\text{ }}{{\left( {\dfrac{{{k_2}h}}{{1 - {k_1}h}}} \right)}^{\frac{1}{{1 - b}}}}} \right\} 得
\left\{\begin{array}{l} { [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} < - {\delta ^b} \leqslant - {\psi ^b}(b)\dfrac{{{d^ * }}}{{{k_2}}} < 0\\ { [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^{1 - b}} > {\delta ^{1 - b}} \geqslant {\psi ^{1 - b}}(b)\left( {\dfrac{{{k_2}h}}{{1 - {k_1}h}}} \right) > 0. \end{array}\right. (33) 则有
\left\{ \begin{array}{l} {k_2}{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} < - {d^ * }{\psi ^b}(b),\\ (1 - {k_1}h){s_n}(k) < {\psi ^{1 - b}}(b){k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b}. \end{array} \right. (34) 由于 - {\tilde d_n}(k) < {d^ * } ,由式(34)可得
\begin{split} & {k_1}h{s_n}(k) + {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} - h{{\tilde d}_n}(k)< \\ &\quad - {k_1}h\delta - {d^ * }h{\psi ^b}(b) + {d^ * }h =\\ &\quad - ({k_1}h\delta + ({\psi ^b}(b) - 1){d^ * }h) = - \mu . \end{split} (35) 由引理1可知, 1 < \psi (b) < 2 ,所以 \mu 为正数.
由于 {\psi ^{1 - b}}(b) > 1 , { [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} < 0 ,由式(34)可得
\begin{split} & {s_n}(k) < {k_1}h{s_n}(k) + {\psi ^{1 - b}}(b){k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} <\\ &\quad {k_1}h{s_n}(k) + {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b}. \end{split} (36) 由于 {\tilde d_n}(k) < {d^ * } ,根据式(36)、(34)可推导出
\begin{split} & 2{s_n}(k) - {k_1}h{s_n}(k) - {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} + h{{\tilde d}_n}(k)< \\ &\quad 2({k_1}h{s_n}(k) + {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b}) - {k_1}h{s_n}(k) -\\ &\quad {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} + {d^ * }h = {k_1}h{s_n}(k) + {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} + \\ &\quad {d^ * }h < - {k_1}h\delta - {d^ * }h{\psi ^b}(b) + {d^ * }h = \\ &\quad - ({k_1}h\delta + ({\psi ^b}(b) - 1){d^ * }h) = - \mu . \end{split} (37) 因此,有 V(k + 1) - V(k) < - {\mu ^2} < 0 .
情况2 当 {s_n}(k) > \delta 时,同理可证明 V(k + 1) - V(k) < - {\mu ^2} < 0 .
综上所述,当 {s_n}(k) < - \delta < 0 或 {s_n}(k) > \delta > 0 时,即 {s_n}(k) \notin \varLambda 时,均有 V(k + 1) - V(k) < - {\mu ^2} < 0 成立.
记滑模变量的初始值为 {s_n}(0) ,由式(30)得
s_n^2(k) < s_n^2(k - 1) - {\mu ^2} < \cdots < s_n^2(0) - k{\mu ^2}. (38) 令 s_n^2(0) - k{\mu ^2} = {\delta ^2} ,解得 k = (s_n^2(0) - {\delta ^2})/{\mu ^2} . 即系统状态远离滑模面时可在 {k^ * } 步内进入区域{\rm{\varLambda}}.
步骤2 (滑动段) 证明系统状态一旦进入区域{\rm{\varLambda}},之后将始终保持在区域{\rm{\varLambda}},即若 |{s_n}(k)| \leqslant \delta ,则始终有 |{s_n}(k + 1)| \leqslant \delta 成立.
推论1 记 \omega = \max \left\{ {{{\left( {\dfrac{{{d^ * }}}{{{k_2}}}} \right)}^{\frac{1}{b}}},{\text{ }}{{\left( {\dfrac{{{k_2}h}}{{1 - {k_1}h}}} \right)}^{\frac{1}{{1 - b}}}}} \right\} ,则有
{d^ * }h \leqslant {k_2}h{\omega ^b} \leqslant (1 - {k_1}h)\omega . (39) 证明 当 \omega = {\left( {\dfrac{{{d^ * }}}{{{k_2}}}} \right)^{\frac{1}{b}}} \geqslant {\left( {\dfrac{{{k_2}h}}{{1 - {k_1}h}}} \right)^{\frac{1}{{1 - b}}}} 时,有 (1 - {k_1}h){\omega ^{1 - b}} \geqslant {k_2}h ,故 (1 - {k_1}h)\omega \geqslant {k_2}h{\omega ^b} = {d^ * }h ;
当 \omega = {\left( {\dfrac{{{k_2}h}}{{1 - {k_1}h}}} \right)^{\frac{1}{{1 - b}}}} \geqslant {\left( {\dfrac{{{d^ * }}}{{{k_2}}}} \right)^{\frac{1}{b}}} 时,有 {k_2}{\omega ^b} \geqslant {d^ * } , (1 - {k_1}h){\omega ^{1 - b}} = {k_2}h ,故 (1 - {k_1}h)\omega = {k_2}h{\omega ^b} \geqslant {d^ * }h . 原式成立,推论1证毕.
设 {s_n}(k) = \gamma \delta = \gamma \omega \psi (b) ,其中 - 1 \leqslant \gamma \leqslant 1 ,故 \delta > 0 ,由式(27)可得
\begin{split} & {s_n}(k + 1) = {s_n}(k) - {k_1}h{s_n}(k) - {k_2}h{ [\kern-0.15em[ {s_n}(k) ]\kern-0.15em] ^b} + \\ &\quad h{{\tilde d}_n}(k) = (1 - {k_1}h)\omega \gamma \psi (b) - {k_2}h|\omega \gamma \psi (b){|^b}{{\mathrm{sgn}}} (\gamma) + \\ &\quad h{{\tilde d}_n}(k) = (1 - {k_1}h)\omega |\gamma |\psi (b){{\mathrm{sgn}}} (\gamma ) - {k_2}h{\omega ^b}|\gamma {|^b} \times \\ &\quad {\psi ^b}(b){{\mathrm{sgn}}} (\gamma ) + h{{\tilde d}_n}(k). \end{split} (40) 当 |{s_n}(k)| < \delta 时,考虑以下2种情况:
情况1 当 - 1 \leqslant \gamma < 0 时,由式(40)得
\begin{split} &{s_n}(k + 1) = - (1 - {k_1}h)\omega {\text{|}}\gamma |\psi (b) +\\ &\quad{k_2}h{\omega ^b}{\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) + h{\tilde d_n}(k). \end{split} (41) 由于 {\tilde d_n}(k) < {d^ * } ,所以
\begin{split} &{s_n}(k + 1) < - (1 - {k_1}h)\omega {\text{|}}\gamma |\psi (b) + \\ &\quad{k_2}h{\omega ^b}{\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) + {d^ * }h. \end{split} (42) 由推论1可知,{d^ * }h \leqslant {k_2}h{\omega ^b} \leqslant (1 - {k_1}h)\omega ,则有
\begin{split} & {s_n}(k + 1) < - (1 - {k_1}h)\omega {\text{|}}\gamma |\psi (b) + \\ &\quad (1 - {k_1}h)\omega {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) +(1 - {k_1}h)\omega = \\ &\quad (1 - {k_1}h)\omega ( - {\text{|}}\gamma |\psi (b) + {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) + 1). \end{split} (43) 因为 0 < {\text{|}}\gamma | \leqslant 1 ,由引理1可得
{s_n}(k + 1) < (1 - {k_1}h)\omega \psi (b) < \omega \psi (b) = \delta . (44) 另外,由于 {\tilde d_n}(k) > - {d^ * } ,由式(41)可得
\begin{split} & {s_n}(k + 1) > - (1 - {k_1}h)\omega {\text{|}}\gamma {\text{|}}\psi (b) +\\ &\quad {k_2}h{\omega ^b}{\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) - {d^ * }h. \end{split} (45) 由推论1可知 - {d^ * }h \geqslant - {k_2}h{\omega ^b},则
\begin{split} & {s_n}(k + 1) > - (1 - {k_1}h)\omega {\text{|}}\gamma {\text{|}}\psi (b) + {k_2}h{\omega ^b}{\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) - \\ &\quad {k_2}h{\omega ^b} = - (1 - {k_1}h)\omega {\text{|}}\gamma {\text{|}}\psi (b) +\\ &\quad {k_2}h{\omega ^b}({\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) - {\text{1)}}{\text{.}} \end{split} (46) 若 \gamma \psi (b) \leqslant - 1 ,则 {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) - {\text{1 > 0}} ,故
{s_n}(k + 1) > - (1 - {k_1}h)\omega {\text{|}}\gamma {\text{|}}\psi (b) > - \omega \psi (b) = - \delta . (47) 若 - 1 < \gamma \psi (b) \leqslant 0 ,则 {\text{1}} - {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b){\text{ > 0}} ,由推论1可知 - {k_2}h{\omega ^b} \geqslant - (1 - {k_1}h)\omega ,由式(46)可推出
\begin{split} & {s_n}(k + 1) > - (1 - {k_1}h)\omega {\text{|}}\gamma {\text{|}}\psi (b) - {k_2}h{\omega ^b}(1 - {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b)) \geqslant \\ &\quad - (1 - {k_1}h)\omega {\text{|}}\gamma {\text{|}}\psi (b) - (1 - {k_1}h)\omega (1 - {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b)) =\\ &\quad - (1 - {k_1}h)\omega ({\text{|}}\gamma {\text{|}}\psi (b) - {\text{|}}\gamma {{\text{|}}^b}{\psi ^b}(b) + 1). \\[-41pt] \end{split} (48) 因为 0 < {\text{|}}\gamma | \leqslant 1 ,由引理1可得
{s_n}(k + 1) > - (1 - {k_1}h)\omega \psi (b) \geqslant - \omega \psi (b) = - \delta . (49) 综上,当 - 1 \leqslant \gamma < 0 时, |{s_n}(k + 1)| \leqslant \delta 成立.
情况2 当 0 \leqslant \gamma < 1 时,同理可证明有 |s(k + 1)| \leqslant \delta .
综上所述,当 |{s_n}(k)| \leqslant \delta 时, |{s_n}(k + 1)| \leqslant \delta 恒成立,即一旦系统状态进入区域{\rm{\varLambda}},则之后任意时刻都将始终保持在区域{\rm{\varLambda}}内.
根据步骤1和步骤2的证明结果,定理2得证.
3.2 跟踪误差界
引理2[10] 考虑有界序列{e_n}(k) ,若\max |{e_n}(k - 1)| \leqslant \rho ,\rho \in {\mathbb{R}^ + },则其分数阶微积分也有界,即
\frac{1}{{h}^{a-1}}{\displaystyle \sum _{j=1}^{L}{(-1)}^{j}}\left(\begin{array}{c}a-1\\ j\end{array}\right){e}_{n}(k-j)\leqslant \frac{(K-1)\rho }{{h}^{a-1}}\text{,}\;\;\;L\in {\mathbb{Z}}^{ + }, (50) 式中:K = \dfrac{{\varGamma (L + 2 - a)}}{{\varGamma (2 - a)\varGamma (L + 1)}} \geqslant 1,\varGamma (x) = \displaystyle\int_0^\infty {{\tau ^{x - 1}}e_n^{ - \tau }{\text{d}}\tau } .
根据定理2,当滑模变量 {s_n}(k) 进入界内后,有
{s_n}(k) = \Delta {e_n}(k) + {l_1}{e_n}(k) + {l_2}{\Delta ^{a - 1}}{ [\kern-0.15em[ {e_n}(k) ]\kern-0.15em] ^q} \leqslant \delta , (51) \begin{split} &\Delta {e_n}(k) + {l_1}{e_n}(k) + \\ &\quad \frac{{{l_2}}}{{{h^{a - 1}}}}\sum\limits_{j = 0}^L {{{( - 1)}^j}\left( {\begin{array}{*{20}{c}} {a - 1} \\ j \end{array}} \right)} { [\kern-0.15em[ {e_n}(k - j) ]\kern-0.15em] ^q} \leqslant \delta , \end{split} (52) 所以
\begin{split} & \Delta {e_n}(k) + {l_1}{e_n}(k) + \frac{{{l_2}}}{{{h^{\alpha - 1}}}}{ [\kern-0.15em[ {e_n}(k) ]\kern-0.15em] ^q}+ \\ &\quad \frac{{{l_2}}}{{{h^{a - 1}}}}\sum\limits_{j = 1}^L {{{( - 1)}^j}\left( {\begin{array}{*{20}{c}} {a - 1} \\ j \end{array}} \right)} { [\kern-0.15em[ {e_n}(k - j) ]\kern-0.15em] ^q} \leqslant \delta . \end{split} (53) 由引理2可得
\Delta {e_n}(k) + {l_1}{e_n}(k) + {l_2}{h^{1 - a}}{ [\kern-0.15em[ {e_n}(k) ]\kern-0.15em] ^q} \leqslant \delta + {l_2}\frac{{(K - 1)\rho }}{{{h^{a - 1}}}}, (54) 式中:\rho \geqslant \max { [\kern-0.15em[ {e_n}(k - 1) ]\kern-0.15em] ^q} .
由定理2可得,跟踪误差 {e_n}(k) 最终有上界,即
|{e}_{n}(k)|\leqslant \psi (q)\mathrm{max}\left\{{\left(\frac{\vartheta {h}^{a-1}}{{l}_{2}}\right)}^{\frac{1}{q}}\text{,}{\left(\frac{{l}_{2}{h}^{2-a}}{1-{l}_{1}h}\right)}^{\frac{1}{1-q}}\right\}, (55) 式中: \vartheta = \delta + {l_2}\dfrac{{(K - 1)\rho }}{{{h^{a - 1}}}} .
因此,当滑模变量进入区域{\rm{\varLambda}}后,相应的跟踪误差将在跟踪误差界内运动. 当 {\left( {\dfrac{{\vartheta {h^{a - 1}}}}{{{l_2}}}} \right)^{\frac{1}{q}}} \leqslant {\left( {\dfrac{{{l_2}{h^{2 - a}}}}{{1 - {l_1}h}}} \right)^{\frac{1}{{1 - q}}}} ,即 {l_2} \geqslant \dfrac{{{\vartheta ^{1 - q}}{{(1 - {l_1}h)}^q}}}{{{h^{1 + q - a}}}} 时,跟踪误差界由后者决定. 由于 0 < h \ll 1{\text{ s}} ,所以
{\left( {\frac{{{l_2}{h^{2 - a}}}}{{1 - {l_1}h}}} \right)^{\frac{1}{{1 - q}}}} < {\left( {\frac{{{l_2}}}{{1 - {l_1}h}}} \right)^{\frac{1}{{1 - q}}}}. (56) 此时,相较于整数阶滑模控制,分数阶滑模控制方法的跟踪误差界更小.
4. 磁悬浮转台控制系统验证
磁悬浮转台实验平台如图4所示,系统参数如表1所示. 动子由圆周型Halbach永磁阵列、背板以及测量目标物构成,定子由8个线圈、支承架和基座构成,6个位移传感器可实时测量动子的六自由度位移. 控制器采用NI的PXIe-8880和多功能IO模块PXIe-7856R. 首先,由PXIe-7856R对6路传感器信号进行A/D转换,经位移解算后数据被送入PXIe-8880控制器运算后输出控制信号,由PXIe-7856R模块产生模拟电压,通过PA12A进行功率放大并转换为驱动电流后送入8个线圈.
表 1 磁悬浮转台系统参数Table 1. Systematic parameters of maglev rotary table参数 数值 永磁体长(lm)、宽(wm)、高(hm)/mm 30、8、8 φm/(°) 7.5 r/mm 80 Br/T 1.2 线圈尺寸(lc、wc、hc、rc)/mm 60、10、10、10 线圈匝数 N/匝 300 m/kg 2.37 {J_\alpha }({J_\beta })、{J_\gamma }/(kg·m2) 9.37×10−3 、1.87×10−2 采样间隔/ms 1 重力加速度/(m·s−2) 9.8 将所提方法与PID和整数阶离散时间滑模控制(DTSMC)的运动控制效果进行对比,取a = 0.5,N = 10. 其中,PID参数kp = 50.7,ki = 1,kd = 107.4,且各自由度的增益KPID如表2所示;DTSMC和所提方法参数中取b = 3/5,q = 3/5,各自由度的其余参数如表2所示.
表 2 控制器参数Table 2. Controller parameters自由度 PID DTSMC 和所提方法 KPID k1 k2 l1 l2 x,y,z 0.124 1 19.3 8.84 21.4 α,β 4.90 × 10−4 3.95 × 10−3 0.0763 0.0348 0.0846 γ 9.78 × 10−4 7.89 × 10−3 0.152 0.0694 0.169 各实验中动子的初始位移为(0, 0, 2 mm, 0, 0, 0)T. 设总采样点数为 M ,跟踪误差的均方根和最大值分别为
\left\{\begin{array}{l} {e}_{\text{RMS}}=\sqrt{\dfrac{1}{M}{\displaystyle \sum _{k=0}^{M-1}{e}_{n}^{2}(k)}}\text{, }\\ {e}_{\text{MAX}}=\underset{k=0,\cdots ,M-1}{\mathrm{max}}\text{ }{e}_{n}(k). \end{array}\right. 4.1 控制算法耗时
在1 s运行时间内,3种控制方法耗时的平均值分别为0.087 、0.220、0.244 s,满足采样周期1.000 ms的要求. 相对于DTSMC算法,分数阶微积分的引入使得所提方法的耗时增加了约10%,可通过改变式(21)中的L值来调节算法耗时.
4.2 阶跃信号跟踪
以1 s为时间间隔,在6个自由度上依次跟踪幅度为5.0 mm、5.0 mm、0.3 mm、0.87 mrad(0.05°)、0.87 mrad(0.05°)、5.24 mrad(0.3°)的阶跃信号,验证磁悬浮转台的阶跃响应性能, 跟踪效果如图5所示.
由图5可知,磁悬浮转台在所提出的分数阶滑模控制策略下,在保证调节时间短的条件下,超调量更小,综合来看,具有最佳的阶跃响应效果.
4.3 三角波跟踪
在x、y、γ自由度对频率为1 Hz,峰值分别为5 mm、5 mm、5.24 mrad(0.3°)的三角波进行跟踪,验证工作台的往复平移和旋转性能. 跟踪效果和跟踪误差ex、 ey和 ez的波形如图6所示, {e_{{\text{RMS}}}} 和 {e_{{\text{MAX}}}} 如表3所示.
根据表3数据:磁悬浮转台在所提出的分数阶滑模控制策略下,对三角波的跟踪误差的均方根和最大值均最小; 与整数阶滑模控制方法相比,采用所提方法,x、y、γ自由度的跟踪误差均方根分别减小了12.8%、16.8%和23.7%,最大值分别减小了9.26%、13.00%和33.20%.
表 3 三角波跟踪误差Table 3. Tracking errors of triangular wave控制方法 x/mm y/mm γ/mrad eRMS eMAX eRMS eMAX eRMS eMAX PID 0.0148 0.0966 0.0155 0.0977 0.0706 0.212 DTSMC 0.0117 0.0875 0.0149 0.0953 0.0410 0.170 所提方法 0.0102 0.0794 0.0124 0.0829 0.0313 0.113 4.4 圆轨迹跟踪
令磁悬浮转台在x-y平面跟踪半径为5 mm的圆轨迹,跟踪效果如图7所示,对应跟踪误差的均方根 {e_{{\text{RMS}}}} 和最大值 {e_{{\text{MAX}}}} 如表4所示.
由表4数据可知:与整数阶滑模控制方法相比,采用所提方法,x、y自由度的跟踪误差均方值分别减小了6.39%和12.40%,最大值分别减小了9.90%和12.10%,磁悬浮转台在所提出的分数阶滑模控制策略下对圆轨迹的跟踪效果最优.
表 4 圆轨迹跟踪误差Table 4. Tracking errors of circular trajectorymm 控制方法 x y eRMS eMAX eRMS eMAX PID 0.0109 0.0686 0.0226 0.163 DTSMC 0.0101 0.0657 0.0251 0.158 所提方法 0.0091 0.0615 0.0220 0.138 5. 结 论
1) 为抑制系统非线性、强耦合性和不确定性对磁悬浮转台运动精度的影响,设计了基于扰动补偿的离散时间分数阶滑模控制策略以提高跟踪精度,并进行了实验验证.
2) 针对包含系统误差、六自由度间耦合项和外界干扰的集总扰动,设计了非线性干扰观测器进行扰动补偿,证明了估计误差有界且可调节到任意小.
3) 利用分数幂函数对抖振的抑制能力和分数阶微积分的全局记忆性设计运动控制器,证明了闭环系统稳定性和有限时间收敛特性.
4) 开展了磁悬浮转台系统六自由度轨迹跟踪实验,分别对阶跃信号、三角波和圆轨迹进行跟踪,系统在所提出的分数阶滑模控制策略下的跟踪误差最小,验证了所提控制策略的优越性.
-
表 1 磁悬浮转台系统参数
Table 1. Systematic parameters of maglev rotary table
参数 数值 永磁体长(lm)、宽(wm)、高(hm)/mm 30、8、8 φm/(°) 7.5 r/mm 80 Br/T 1.2 线圈尺寸(lc、wc、hc、rc)/mm 60、10、10、10 线圈匝数 N/匝 300 m/kg 2.37 {J_\alpha }({J_\beta })、{J_\gamma }/(kg·m2) 9.37×10−3 、1.87×10−2 采样间隔/ms 1 重力加速度/(m·s−2) 9.8 表 2 控制器参数
Table 2. Controller parameters
自由度 PID DTSMC 和所提方法 KPID k1 k2 l1 l2 x,y,z 0.124 1 19.3 8.84 21.4 α,β 4.90 × 10−4 3.95 × 10−3 0.0763 0.0348 0.0846 γ 9.78 × 10−4 7.89 × 10−3 0.152 0.0694 0.169 表 3 三角波跟踪误差
Table 3. Tracking errors of triangular wave
控制方法 x/mm y/mm γ/mrad eRMS eMAX eRMS eMAX eRMS eMAX PID 0.0148 0.0966 0.0155 0.0977 0.0706 0.212 DTSMC 0.0117 0.0875 0.0149 0.0953 0.0410 0.170 所提方法 0.0102 0.0794 0.0124 0.0829 0.0313 0.113 表 4 圆轨迹跟踪误差
Table 4. Tracking errors of circular trajectory
mm 控制方法 x y eRMS eMAX eRMS eMAX PID 0.0109 0.0686 0.0226 0.163 DTSMC 0.0101 0.0657 0.0251 0.158 所提方法 0.0091 0.0615 0.0220 0.138 -
[1] XU X Z, ZHENG C L, XU F Q. A real-time numerical decoupling method for multi-DoF magnetic levitation rotary table[J]. Applied Sciences, 2019, 9(16): 3263.1-3263.16. [2] ZHU H Y, TEO T J, PANG C K. Flexure-based magnetically levitated dual-stage system for high-bandwidth positioning[J]. IEEE Transactions on Industrial Informatics, 2019, 15(8): 4665-4675. doi: 10.1109/TII.2019.2890951 [3] LAHDO M, STRÖHLA T, KOVALEV S. Design and implementation of an new 6-DoF magnetic levitation positioning system[J]. IEEE Transactions on Magnetics, 2019, 55(12): 8107407.1-8107407.7. [4] SILVA-RIVAS J C, KIM W J. Multivariable control and optimization of a compact 6-DOF precision positioner with hybrid H2/H∞ and digital filtering[J]. IEEE Transactions on Control Systems Technology, 2013, 21(5): 1641-1651. doi: 10.1109/TCST.2012.2215035 [5] LI D F, GUTIERREZ H. Observer-based sliding mode control of a 6-DOF precision maglev positioning stage[C]//2008 34th Annual Conference of IEEE Industrial Electronics. Orlando: IEEE, 2008: 2562-2567. [6] 魏静波,罗浩,关子津. 基于干扰观测器的磁悬浮球系统全局快速终端滑模控制[J]. 西南交通大学学报,2023,58(4): 836-844.WEI Jingbo, LUO Hao, GUAN Zijin. Global fast terminal sliding mode control for maglev ball system based on disturbance observer[J]. Journal of Southwest Jiaotong University, 2023, 58(4): 836-844. [7] DU H B, YU X H, LI S H. Dynamical behaviors of discrete-time fast terminal sliding mode control systems[M]//Recent Advances in Sliding Modes: From Control to Intelligent Mechatronics. Cham: Springer, 2015: 77-97. [8] ZHANG L, ZHANG Z, HUANG L. Hybrid non-linear differentiator design for a permanent-electro magnetic suspension maglev system[J]. IET Signal Processing, 2012, 6(6): 559-567. doi: 10.1049/iet-spr.2011.0264 [9] NGUYEN S D, LAM B D, NGO V H. Fractional-order sliding-mode controller for semi-active vehicle MRD suspensions[J]. Nonlinear Dynamics, 2020, 101(2): 795-821. doi: 10.1007/s11071-020-05818-w [10] SUN G H, MA Z Q, YU J Y. Discrete-time fractional order terminal sliding mode tracking control for linear motor[J]. IEEE Transactions on Industrial Electronics, 2018, 65(4): 3386-3394. doi: 10.1109/TIE.2017.2748045 [11] KUANG Z A, GAO H J, TOMIZUKA M. Precise linear-motor synchronization control via cross-coupled second-order discrete-time fractional-order sliding mode[J]. IEEE/ASME Transactions on Mechatronics, 2021, 26(1): 358-368. [12] HU C X, WANG Z, ZHU Y, et al. Performance-oriented precision LARC tracking motion control of a magnetically levitated planar motor with comparative experiments[J]. IEEE Transactions on Industrial Electronics, 2016, 63(9): 5763-5773. doi: 10.1109/TIE.2016.2538743 [13] JANSEN J W, VAN LIEROP C M M, LOMONOVA E A, et al. Magnetically levitated planar actuator with moving magnets[J]. IEEE Transactions on Industry Applications, 2008, 44(4): 1108-1115. doi: 10.1109/TIA.2008.926065 [14] DYCK M, LU X D, ALTINTAS Y. Magnetically levitated rotary table with six degrees of freedom[J]. IEEE/ASME Transactions on Mechatronics, 2017, 22(1): 530-540. doi: 10.1109/TMECH.2016.2621108 [15] YANG H L, DENG F, HE Y, et al. Robust nonlinear model predictive control for reference tracking of dynamic positioning ships based on nonlinear disturbance observer[J]. Ocean Engineering, 2020, 215: 107885.1-107885.7. [16] LU X, XU F Q, XU X Z, et al. Directed-driven 8-phase magnetically levitated rotary table based on an analytical-numerical model[J]. IEEE Access, 2020, 8: 31159-31170. doi: 10.1109/ACCESS.2020.2973223 -