皮纳卫星具有成本低、研制周期短、风险小但技术集成度高等优点, 受到了广泛关注.特别地, 立方体皮纳卫星具有标准化、模块化、通用性强等优点, 能够快速响应任务需要, 适应空间应用发展趋势^[1-2], 实用化立方星技术入选《科学》杂志“2014年全球十大科学突破”^[3].随着皮纳卫星承担任务需求难度的增大(如：对地成像), 姿态控制要求也逐步提高^[4-6], 被动式或者纯磁力矩控制难以满足要求.反作用飞轮具有控制精度高、响应速度快、寿命长等优点^[7-8], 其性能直接影响了姿态控制精度, 因此有必要开展实用的高精度微型飞轮研究.

飞轮输出形式为角动量, 本质上属于速度控制系统, 因此速度测量估计和控制策略直接影响了飞轮的控制输出精度.受到微型飞轮体积(≤36 mm×36 mm×30 mm)、功耗(≤500 mW)和微处理器计算能力(主频为32 MHz)等限制, 常规飞轮提高性能的方法难以应用.光电编码器、旋转变压器以及测速电机等高精度测速传感器受到光源寿命、系统复杂性、体积庞大等限制不能应用.朱娜等^[9]通过将飞轮霍尔元件(HALL)数量增加到27个, 来提高角位置分辨率, 从而提高测速精度.但是, 由于微型飞轮空间极为有限, 这种方法难以应用.陈冬等^[10]采用定常数增益的卡尔曼滤波算法对控制力矩陀螺进行角位置和速度实时估计, 并在高速DSP TMS320F2812(主频180 MHz)上进行验证, 运算时间为12.4 μs.但微型飞轮微处理器计算能力很差, 在控制周期内很难完成实时估计.Cheon等^[11]采用均值滤波提高飞轮速度测量精度, 但这不利于低速运行.Hoevenaars等^{[5, 12]}采用两种PID控制策略减小了过零区域非线性状态的影响, 但没进行实物测试验证算法的有效性, 只完成了数学仿真.Zheng等^[13-14]采用模型参考自适应控制算法, 并在硬件平台进行验证, 计算量大, 微型飞轮处理器难以满足实时性要求.Sinclair等^[6]仅提出了高速功率控制和低速速度控制双闭环控制方案, 硬件电路复杂, 且未给出实验验证结果.

针对微型飞轮系统尺寸空间、功耗以及处理器的计算能力等方面的局限性, 本文采用超低功耗微处理器作为控制器、开关式HALL测速、无刷直流电机驱动、PWM速率控制输出的硬件总体方案.本文提出一种将变分频T法与带速度切换迟滞环的滑动均值滤波相结合的速度估计方法, 以及离散PID控制算法, 满足飞轮的体积、质量、功耗、控制精度等总体要求.

1 总体方案

飞轮系统一般可采用速率控制模式和力矩控制模式.在力矩控制模式下, 采用电流反馈, 仅能实现电磁力矩闭环, 不能克服飞轮内部力矩干扰, 不能降低飞轮执行指令的误差.特别是当采用滚珠轴承支撑转子时, 由于摩擦力矩干扰大, 会造成飞轮执行指令产生很大误差.当采用速度控制模式时, 以转速为反馈量实现转速闭环控制.由于转速是输出力矩的积分, 转速测量值中含有干扰力矩的信息, 包含了飞轮的动力学特性.转速闭环可将飞轮系统的摩擦力矩等干扰包括在反馈环内, 从而抑制系统内部干扰, 提高飞轮系统的控制精度.

摩擦干扰力矩会对卫星控制精度造成较大的影响, 必须对摩擦干扰力矩进行抑制.为了获得较高的控制精度, 本文采用速度控制方案, 如图 1所示.考虑功耗、总线接口、体积、计算速度等多方面因素, 采用超低功耗、主频32 MHz微处理器作为控制器, 三相无刷永磁直流电机作为执行机构.通过中断的方式对电机内部集成的HALL进行采样, 实现换相控制和速度测量, 与转速指令进行比较, 经控制算法计算后, 通过PWM输出驱动MOSFET管桥电路, 实现转速闭环，实时监控电流和温度, 实现对飞轮的过流、过温保护.

2 基于HALL的速度估计方法

基于开关式HALL进行速度估计的基本原理：无刷直流电机定子绕组采用“Y”型连接, 并嵌入3个安装位置按120°分布的HALL.转子磁极经过HALL附近时, 其信号状态发生改变.易知，转子每转过60°电角度, 其中一个HALL就会改变状态.因此, 在一个完整的电周期内, 3路HALL信号输出相位差为60°、占空比为50%的脉冲信号.这3路信号异或后将电周期六等分为θ=60°的HALL位置区间, 完成一圈完整的转子机械转动周期要重复的电周期数取决于转子磁极对数P.HALL输出信号和电周期、机械周期的关系如图 2所示.

由于HALL信号表征的是离散化的角位置增量θ, 对θ进行微分运算, 可以估计出在该HALL位置区间内的平均速度值：

$ \omega = \frac{\theta }{{tP}}. $

(1)

根据式(1)可知, 电机内HALL安装位置偏差和处理器对HALL信号的时间间隔t测量误差均会造成速度估计值ω误差.本文采用以下2种方法提高测速精度.

2.1 变分频T法测速

对于脉冲信号, 常用测速法有T法、M法和M/T法^[15].以上方法各有优缺点, 其中M/T法理论上可以得到优于前两者的测速精度, 但实现方法复杂, 实现比较困难.T法测速是利用高频时钟脉冲, 通过测量2个相邻输入脉冲之间的时间来确定被测速度.本文采取适用于低性能处理器的变分频T法进行测速.

对HALL位置区间的时间测量方法, 如图 3所示, 3路HALL信号经过异或逻辑后, 触发处理器的定时器f_T中断, 根据捕获的计数值N即可求得在HALL位置区间的信号时间间隔：

$ t = \frac{{{p_{{\rm{pre}}}}N}}{{{f_0}}}. $

(2)

由于T法在极端情况下对时间的测量会出现±1个时钟脉冲周期的误差, 定时器的频率f_T越高, 误差越小.但计数器最大计数值N_max一定, 如果f_T过高, 会导致飞轮低速转动时计数值的溢出.可使用FPGA或高速DSP提高系统主频f₀和N_max来解决.考虑飞轮硬件性能指标, f₀=32 MHz, N_max=65 536, n_max=15 000 r/min, 需要调节定时器分频值p_pre防止计数器溢出.为提高测量准确度, 本文采用基于中断捕获计数值N动态调整分频值p_pre的T法测速.由于分频值p_pre决定着定时器的频率f_T, 通过实时改变f_T值可以提高测速精度.分频值p_pre预装载在处理器寄存器, 可在工作时被改变, 而新分频值p_pre仅在下一个脉冲中断捕获事件时起作用, 不会干扰当前测量结果.定时器的分频值p_pre根据实时转速值进行调整.当转速较高时, 为了获得较高的分辨率, 减小分频值p_pre, 可以适当地增加相邻捕获脉冲间的脉冲数, 即增大计数器值, 如图 4(a)所示；当速度过低, 计数溢出时, 增大分频值p_pre, 降低计数器值, 如图 4(b)所示.根据以上分析可知, 飞轮转速计算公式为

$ n = 60\frac{\omega }{{2{\rm{ \mathsf{ π} }}}} = \frac{{30\theta }}{{{\rm{ \mathsf{ π} }}tP}} = \frac{{30{f_0}\theta }}{{{\rm{ \mathsf{ π} }}{p_{{\rm{pre}}}}PN}}. $

(3)

2.2 滑动均值滤波算法

在飞轮电机定子中安装的HALL, 不可避免地会出现与理想安装位置的偏差, 导致HALL信号所对应的电角度不一致, 进而对转速的测量会造成一定的影响.当多相相轴均出现安装偏差时, 情况更加复杂.

为减小机械安装偏差, 可以采用在高速下标定补偿的方法, 但标定过程对飞轮恒速控制精度要求很高^[16].本文采用滑动均值滤波方法来解决.假定在第i位置区间内的速度ω_i为最近的m个HALL位置区间内的平均速度, 则存在下式：

$ \begin{array}{l} {\omega _i} = \frac{{{\theta _{i - m + 1}} + {\theta _{i - m + 2}} + \cdots + {\theta _i}}}{{\Delta {t_{i - m + 1}} + \Delta {t_{i - m + 2}} + \cdots + \Delta {t_i}}} = \\ \;\;\;\;\;\;\;\frac{{m \times {{60}^ \circ }}}{{\Delta {t_{i - m + 1}} + \Delta {t_{i - m + 2}} + \cdots + \Delta {t_i}}}. \end{array} $

(4)

式中：θ_i为分区i的电角度, Δt_i为分区i所需的时间, m为滤波长度.

由式(4)可知, 当m=6时, 相当于电机转子旋转了一个完整的电周期, HALL安装偏差导致信号对应的电角度偏差之和为0, 从而最大限度地减小安装偏差导致的测速误差.用于本HALL位置区间计算的速度实际上是前m个区间的平均速度, 且m值越大, 滤波越明显, 这会导致速度测量的延迟较大, 实时性较差, 动态性能较低, 因此, 该方法适用于较高速运行的状态, 不适用于低速区和速度变化较大的情况.在较高转速时选择m=6, 以追求稳态测速精度, 在起动和低速时选择m=1, 以追求动态响应速度.为避免这2种方式在临界转速附近频繁来回切换, 导致飞轮响应出现变慢甚至停转的现象, 在控制算法上设置速度切换迟滞环, 如图 5所示.

3 离散PID控制算法

根据牛顿定律和动量矩定理, 推导飞轮运动学模型^[5]：

$ {U_{\rm{c}}} = \frac{{JR}}{{{k_{\rm{m}}}}}\dot \omega + \left( {\frac{{R{C_{\rm{v}}}}}{{{k_{\rm{m}}}}} + {k_{\rm{e}}}} \right)\omega + \frac{{R{C_{\rm{0}}}}}{{{k_{\rm{m}}}}}. $

(5)

式中：J为转动惯量，R为电机相电阻，k_m为电机转矩常数，k_e为电机反电动势常数，C₀为电机静态摩擦转矩，C_v为电机动态摩擦转矩系数，U_c为控制器输出PWM等效控制电压，ω为飞轮的旋转速度.

飞轮系统要属于速度调节系统, 对阶跃转速指令的响应无稳态误差, 系统应为Ⅰ型系统, 因此应选取PI控制器：

$ {G_{\rm{c}}}\left( s \right) = \frac{{{k_{\rm{p}}}s + {k_i}}}{s}. $

(6)

式中：k_p为PID控制比例控制项，k_i为PID控制积分控制项.

综上, 飞轮系统闭环传递函数为

$ \mathit{\Phi }\left( s \right) = \frac{{{k_{\rm{m}}}{k_{\rm{p}}}s + {k_{\rm{m}}}{k_{\rm{i}}}}}{{JR{s^2} + \left( {{k_{\rm{m}}}{k_{\rm{e}}} + R{C_{\rm{v}}} + {k_{\rm{m}}}{k_{\rm{p}}}} \right)s + {k_{\rm{m}}}{k_{\rm{i}}}}}. $

(7)

为k_p和k_i, 可根据卫星总体对飞轮稳态和动态性能指标要求来计算, 其中最大速度响应时超调量σ=1%, 调节时间t_s=30 s.卫星姿控系统更关注动态响应特性, 飞轮系统应属于欠阻尼二阶系统, 根据二阶系统标准方程式可知：

$ \sigma = \exp \left( { - \frac{{{\rm{ \mathsf{ π} }}\xi }}{{\sqrt {1 - {\xi ^2}} }}} \right) \times 100\% ,\Delta \le \frac{{{{\rm{e}}^{ - \xi {\omega _{\rm{n}}}{t_{\rm{s}}}}}}}{{\sqrt {1 - {\xi ^2}} }}. $

(8)

式中：ξ为阻尼比，ω_n为自然频率，ε为预定门限值.Δ为稳态误差带, 取Δ=0.05.

由式(8)可分别求得ξ、ω_n.由特征方程的极点配置方法可知：

$ {k_i} = \frac{{{T_{\rm{m}}}\omega _n^2}}{K},{k_{\rm{p}}} = \frac{{2\xi \sqrt {K{k_{\rm{i}}}{T_{\rm{m}}}} }}{K} - 1. $

(9)

式中：

$ K = \frac{{{k_{\rm{m}}}}}{{{C_{\rm{v}}}R + {k_{\rm{m}}}{k_{\rm{e}}}}},{T_{\rm{m}}} = \frac{J}{{{C_{\rm{v}}} + {k_{\rm{m}}}{k_{\rm{e}}}/R}}. $

由于飞轮使用数字处理器, 假定控制周期为T_c, 将控制算法离散化可得离散PID控制算法.系统在PID控制器的作用下, 由于积分控制的引入, 使得控制输出不断累加增大, 使得控制输出进入饱和区域, 一旦系统出现反向误差, 控制输出逐渐从饱和区退出, 进入饱和区愈深, 则退出饱和区所需要的时间愈长, 导致此时控制系统反应迟钝, 控制性能下降, 调节品质恶化.为克服以上问题, 对传统的PID控制算法进行改进, 采用积分分离的PID算法.当偏差大于某个规定的门限值时, 取消积分作用；只有当偏差较小时, 引入积分控制, 以便消除静差, 提高精度：

$ \begin{array}{l} {U_{\rm{c}}}\left( k \right) = {k_{\rm{p}}}e\left( k \right) + \beta {k_i}\sum\limits_{j = 0}^k {e\left( j \right)T} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\frac{{{k_{\rm{p}}}}}{T}\left[ {e\left( k \right) - e\left( {k - 1} \right)} \right]. \end{array} $

(10)

式中：

$ \beta = \left\{ \begin{array}{l} 0,\;\;\;\;\;\left| {e\left( k \right)} \right| > \left| \varepsilon \right|;\\ 1,\;\;\;\;\;\left| {e\left( k \right)} \right| \le \left| \varepsilon \right|. \end{array} \right. $

其中, ε值的选取对克服积分饱和有重要的影响, 一般需要通过实验整定的方法获得.k_d是PID微分控制项，相当于补偿了被控对象的一个极点, 能够扩大稳定域, 改善系统的动态特性, k_d需要通过实验测试, 采用工程整定的方法获取.

4 实验测试与在轨验证

如图 6所示为所研制的飞轮实物, 飞轮控制系统相关参数如表 1所示, 由式(9)和实验整定可得，k_p=16、k_i=0.002、k_d=0.020、ε=420 r/min.

为验证系统速度估计和控制算法, 分别进行阶跃响应、正弦跟踪等实验, 以测试飞轮稳态和动态特性.为确保飞轮在空间环境下能够正常运行, 开展环境试验, 检测飞轮功能和性能变化.

飞轮时域特性和稳态特性如图 7所示.其中图 7(a)为飞轮在最大转速为9 000 r/min的阶跃指令下实验测试曲线, 飞轮动态特性如下：上升时间为25.1 s, 调节时间为28.6 s, 峰值时间为44.5 s, 超调量为3.2%.飞轮进入稳态(从50 s起)后, 速度误差统计特性如下：最大误差为1.5%, 平均误差为0.07‰, 误差方差为0.65.如图 7(b)所示为飞轮在转速为6 000 r/min的阶跃指令下实验测试曲线, 飞轮动态特性如下：上升时间为7.1 s, 调节时间为8.0 s, 峰值时间为9.3 s, 超调量为3.20%.飞轮进入稳态(从50 s起)后, 速度误差统计特性如下：最大误差为1.02%, 平均误差为0.4‰, 误差方差为0.21.通过对比可知, 飞轮在6 000 r/min转速以上时动态特性变差.

飞轮正弦跟踪、误差特性和低速过零特性曲线分别如图 8所示.正弦指令幅值为6 000 r/min, 频率为0.01 Hz.实验测试结果表明：飞轮跟踪特性总体上良好, 误差≤1.5%.跟踪误差e随着指令速度的增大缓慢增大, 在速度过零附近由于摩擦力矩突变导致误差瞬间变大, 速度跟踪出现迟滞现象.

在完成冲击谱、振动等力学试验后开展热真空试验(温度范围为-20 ℃~+50 ℃, 真空度为10^-3 Pa), 测试飞轮功能和性能, 对比控制系统在常规(标准大气压、室温)条件下与热真空试验条件下的性能变化.如表 2所示为飞轮在转速为6 000 r/min阶跃指令的实验测试对比情况.其中, t_r为上升时间, t_p为峰值时间, σ为超调量, e_a为最大误差, e_r为平均误差, σ_v为误差方差.结果表明：飞轮阶跃响应的动态过程有很小变化, 这是由于真空环境使得飞轮受到空气阻尼变小, 稳态过程误差变化非常小.

如图 9所示为飞轮在幅值为6 000 r/min, 频率为0.01 Hz的正弦实验测试中的热真空环境下不同温度时结果对比情况.可知, 飞轮速度跟踪特性良好, 跟踪误差e特性几无变化, 均小于1.5%, 在低速过零区域也未明显发生变化.

以上常规和环境试验测试结果表明：经过力学环境后飞轮结构无损害, 功能正常, 在热真空环境下能够正常起转、运行和停止, 性能无明显变化.

2016年6月25日, 该飞轮作为姿控系统执行机构, 随12U立方星“翱翔之星”, 搭载“长征七号”发射入轨.在轨运行99 d直到卫星坠入大气层, 期间各项功能正常, 成功实现卫星姿控任务.如图 10所示为地面站接收到的卫星对日跟踪过程中飞轮执行情况的遥测数据.飞轮接收控制指令后进行速度调节, 三轴由初始转速开始改变, 其中X轴只进行加速, Y、Z轴则先减速为负速度后再加速为正速度最后减速为负速度等过程.整个过程持续时间约280 s, 飞轮运行正常, 顺利完成姿控系统任务.运行结果表明：在轨空间环境下, 飞轮速度估计与控制策略可行有效, 满足星上任务需求.

所研制的飞轮兼顾了体积、质量、功耗和性能指标的设计要求, 与国外同类飞轮产品性能参数^[17-19]对比如表 3所示.

5 结论

针对皮纳卫星、立方星的姿控系统精度需求, 本文提出了一套适用于超低功耗微处理器的微型反作用飞轮测速与控制方案, 包括总体架构、速度估计方法和离散控制算法等, 并通过了地面性能测试实验、环境鉴定试验和在轨运行验证.地面实验测试结果如下：阶跃响应调节时间约为8 s, t_p≤10 s, e_a≤1.5%, e_r≤0.4‰, σ_v≈0.25；正弦跟踪误差随转速变大有所增大, 在过零区域由于低速摩擦特性引起误差突变, 但均小于1.5%.在经过力学环境后热真空试验条件下, 飞轮的动态和稳态特性基本保持不变.飞轮随12U立方星“翱翔之星”在轨运行99 d, 工作正常, 成功完成姿态控制系统任务.地面实验和在轨运行结果表明：飞轮所采用的方案有效可行, 速度估计与控制策略能够满足卫星姿态控制系统的任务需求.