200 kW/180 MJ飞轮储能系统转子有限元分析与实验分析

金 梅,颜廷鑫,张立国,李光军,王 玮,王 娜

(1.燕山大学 河北省测试计量技术与仪器重点实验室，河北 秦皇岛 066004;2.北京航空航天大学,北京 100083)

1 引 言

飞轮储能技术一般通过永磁同步电动发电一体机[1]将机械能与电能相互转化的技术,具有响应速度快、大功率充放电、环保无污染、安全高效等优点。从20纪90年代清华大学等单位研究飞轮储能开始,到目前市场上成型的250 kW/3 kWh,仅适用于小功率负载应用,尚缺少大功率、高储能量的研究和应用。飞轮储能系统工作时需要频繁的调整升降速来维持电网质量,因此工作转速范围应排除临界转速,更要避免储能单体的自振频率。刘彩霞等[2]利用有限元分析软件ANSYS,对传感器加载垂直力和剪切力进行有限元分析,得到传感器在三维力作用下的应力应变分布状态,并进行了实验验证。戴兴建等[3]对20 kW/1 kWh飞轮进行了飞轮轴系的动力学研究;史涔溦等[4]对不同材料的超高速飞轮轴系进行了分析;周元伟等[5]等运用控制变量法分析磁轴承刚度与阻尼随控制参数变化下的600 Wh轴系转子的模态频率分布位置及不平衡振动响应规律;文献[6]研究了采用AMB悬浮的刚性转子的非线性动力学,发现转子的径向和轴向振动特性受转子质量的影响;Yang S M[7]发现,当旋转频率接近共振频率时,微型磁悬浮旋转机械的径向振荡变得过大;文献[8]调节控制系统的比例系数,控制系统的刚度可调,从而获得理想的自振频率。

本文以解决电网系统载荷波动大问题为背景,研制了200 kW/180 MJ飞轮储能系统,首先建立以机械轴承为支撑的重载飞轮转子的有限元模型,进行临界转速的求解和谐波响应分析,其次对飞轮系统进行升降速全周期的试验检测,检测飞轮储能系统在充放电过程中的单体振动,最后将实验结果与分析结果进行对比实验。

2 飞轮储能系统

飞轮储能系统主要有电动/发电机、飞轮、轴承、机组外壳、控制系统[9]、刹车电阻等组成。其组成部分如图1所示。

图1 飞轮储能系统示意图Fig.1 Schematic diagram of flywheel energy storage system

飞轮储能系统初始时由电能转化为飞轮高速旋转的机械能,在电网正常运行时时以额定转速运行,当电网频率、电压甚至断电时,飞轮储能系统立刻将机械能转化为电能弥补电网功率并启动发电机[10],其运行原理如图2所示。

图2 飞轮储能控制系统图Fig.2 Diagram of flywheel energy storage control system

本文研究的飞轮转子用40CrNiMoA铸成,轴系弹性模量2.11×1011Pa,泊松比为0.31,转子质量 3 215 kg(不含磁钢),极转动惯量为668 kg·m2,额定转速7 500 r/min。总储能量180 MJ。

飞轮轴系上下端均采用滚动轴承作为径向支撑,安全可靠。轴向采用电磁轴承承担大部分重量,可以大大降低径向轴承重载摩擦损耗。永磁电动/发电机定子绕组安装在壳体上,铷铁硼永磁体安装在电机轴转子上。

3 有限元分析

对于飞轮储能系统轴系,有阻尼的线性系统的振动方程[11]如下:

(1)

将u分为相互垂直方向的x(t)和y(t)，C分为转子系统外部阻尼矩阵Cn和转子系统内部阻尼矩阵Cr,F(t)分为x和y方向的激励列阵Fx、Fy。则转子动力学微分方程[12]又可以表示为:

(2)

做拉普拉斯变换,并设所有坐标的初速度和初位移均为零,则有:

Ms2+(Cr+Cn)s+K-jωCr=0

(3)

式中ω为圆频率。求解式(3)可得:

(4)

从式(4)可以看出,s存在2个根,其中虚部代表自由旋转时的固有频率。自由旋转时,可忽略阻尼的影响,得到无阻尼模态频率为:

(5)

式中:k为刚度系数;ω1为正向进动时的固有频率;ω2为反向进动时的固有频率。根据式(5)可得表1。

表1 不同刚度下的正反进动固有频率Tab.1 Positive and negative precession natural frequencies with different stiffness Hz

3.1 静力学分析

有限元分析法[13]是将实体模型分为有限个网格模型,并通过数值分析求解,选择合适的求解模型计算的结果可以用来设计参考和理论分析。

根据飞轮转子的设计结构,可以拆装为电机轴、飞轮2部分,现将电机轴、飞轮、整体飞轮分别进行静力学分析、模态求解。

假设飞轮轴系为线性结构。因飞轮空腔内气压<20 bar,故仅考虑重力和轴承因素影响,建立实体模型后,根据陀螺效应对飞轮转子划分网格;使用COMBI214模拟支撑特性。飞轮转子的有限元模型以及生成的网络分割分别如图3和图4所示,分别生成30 862、7 787、45 417个节点,17 717、4 394、27 259 个单元。

图3 有限元模型图Fig.3 Finite element model diagram

图4 有限元网格分解图Fig.4 Finite element model diagram

通过对APDL编程进行静力学分析可分别得到飞轮轴系的静态力学分析结果,如图5所示,电机轴受力最大值及最小值分别为8.024 8×107和7.271 2×105Pa,受力部位主要为电机轴与飞轮连接处内部;飞轮受力最大值及最小值分别为4.142 2×108和906 49 Pa,受力部位主要为飞轮中心;轴系整体受力最大值及最小值分别为1.157 8×109和1.044 2×105Pa,受力部位主要为飞轮中心;同时计算得到其最大变形,x方向最大位移分别为3.522 2×10-3m、3.879 6×10-2m、15.869 m,1.990 3×10-5m、4.551 1×10-4m、1.055 2×10-3m。

图5 有限元静力分析图Fig.5 Finite element static analysis diagram

3.2 模态分析

对于转子动力学研究,阻尼是不可忽略的问题[14],为计算飞轮轴系的临界转速,以避免升降速实验中系统装置出现共振[15]而带来的不良影响,因此采用完全阻尼法对飞轮轴系进行模态分析。完全阻尼法采用Lanczos算法计算并得到复数特征值和特征向量。

现假设飞轮轴系的径向支撑是各向同性,则弯曲振动会出现在同一平面。轴承的刚度一般为 5×107～2×108N/mm,而包含轴承座的刚度计算范围在1×107～1×109N/mm,因此取刚度为1×108N/mm,运用第3.1节中的有限元模型继续进行模态分析求解。图6给出了飞轮轴系系统的前5阶模态振型云图,表2列举了飞轮运行过程中可能出现的模态及临界转速，表3所示给出了飞轮轴系系统临界转速表。

图6 飞轮轴系模态振型云图Fig.6 Cloud diagram of flywheel shafting mode

图7所示为飞轮轴系的坎贝尔图,在激励频率下与1阶固有频率相交即是共振频率，图7中三角形显示共振频率发生200 rad/s附近，正反进动稳定。分析表3可得临界转速在1 900 r/min附近,这样就需要在升降速试验中快速跨越临界转速,以防止共振。根据前5阶的振型云图可得最大的位移35.255 mm发生在第5阶模态,前3阶模态的最大值为 17.09 mm。同时,其固有频率主要为31 Hz和52 Hz对应的临界转速也相对确定。

表2 飞轮轴系模态求解列表Tab.2 The modal solution list of flywheel shafting

表3 飞轮轴系的临界转速列表Tab.3 Critical speed list of flywheel shafting r/min

图7 飞轮轴系的坎贝尔图Fig.7 Campbell diagram of flywheel shafting

3.3 谐响应谱分析

为计算出飞轮轴系在固定频率下的响应并得到位移响应所对应的频率曲线进行响应谱分析。

预先在飞轮侧面施加一个X方向的正弦力,幅值为300 N,相位角为0°,在飞轮轴系底端面施加一个Y方向的正弦力,幅值为500 N,相位角为90°。激励频率为0到100 Hz,模拟人为敲击振动响应，因每间隔1 Hz至少求解一次故求解间隔设置为200,即每隔0.5 Hz求解一次,求解方法为完全法。

分别指定电机轴、飞轮、飞轮轴系机械轴承内接触面,设置方向为x方向,频率响应结果分别如图 8～图13所示。由图8可知，电机轴的上端轴承内接触面随着激励频率的增大，其振动位移逐渐增大。在100 Hz激励下，振动位移为0.003 1 mm。仿真结果符合设计要求。

图8 电机轴上端轴承内接触面频率响应曲线Fig.8 Frequency response curve of contact surface in end bearing of motor shaft

图9 飞轮下端轴承内接触面频率响应曲线Fig.9 Frequency response of contact surface in lower bearing of flywheel

飞轮单体经激励后，低频振动较明显，并在 4 Hz、21 Hz、53 Hz左右出现临界点。最大振动位移为0.03 mm，体现出飞轮单体的低频共振特性。

图10 飞轮轴系上端轴承内接触面频率响应曲线Fig.10 Frequency response curve of contact surface in upper end bearing of flywheel shafting

经电机轴与飞轮单体组合而成的飞轮轴系振动图与飞轮单体相似，但临界振动频率均出现逐渐减小的趋势。这与飞轮单体占主要质量关系密切。

图11 电机轴上端轴承内接触面激励频率100 Hz对应的x方向的相位角响应曲线Fig.11 The phase angle response curve in the x direction corresponding to the excitation frequency of the inner contact surface of the bearing at the upper end of the motor shaft is 100 Hz

图12 飞轮下端轴承内接触面激励频率100 Hz对应的x方向的相位角响应曲线Fig.12 The phase Angle response curve in the x direction corresponding to the excitation frequency of the inner contact surface of the flywheel lower bearing is 100 Hz

图13 飞轮轴系上端轴承内接触面激励频率100 Hz对应的x方向的相位角响应曲线Fig.13 The phase angle response curve in the x direction corresponding to the excitation frequency of the inner contact surface of the bearing at the upper end of the flywheel shafting is 100 Hz

将固定激励频率100 Hz分别施加于电机轴、飞轮、飞轮轴系可得：飞轮单体的振动响应与激励反向;电机轴与飞伦轴系的振动响应与激励同向。

通过完全法求解谐响应,可以对电机轴、飞轮、飞轮轴系得到如下结论:

1)电机轴的变形最大值为1.621 3 mm,相对其它频点的差别不大(最大差别1.6 mm),在相位角的变化下,指定方向的最大振幅仅1.019 mm,综上说明电机轴设计尺寸合理;

2)飞轮的总变形在不同频点处的差别不大(最大变形为100 Hz时的1.7619×10-5mm),说明飞轮的频率响应不明显;在相位角变化下,指定方向的变形量与总变形量“差值”比较小,对飞轮的配合结构件影响较小,综上说明飞轮设计尺寸合理;

3)飞轮轴系的总变形及指定方向的振幅都较小,符合设计要求。

4 试验研究

为掌握飞轮储能系统的运行特点,验证仿真结果,需要通过持续的升降速试验进行实时测量。在壳体的上端和下端分别安装CT1050LC加速度传感器(图14)实时测量振动情况,并在电机、轴承等部位安装温度传感器实时监测关键部位温度参数。

图14 试验装置Fig.14 Test monomer and platform

图15所示为飞轮储能系统温度转速表，对照图15,对飞轮储能系统法进行反复的升降速试验,分别从0 r/min升速至4 500 r/min和2 500 r/min降速至 1 800 r/min,在充分的升降速试验后,储能系统内各关键部位温度保持稳定,各关键器件运行正常,进一步说明,飞轮轴系设计规格符合设计要求。

为检测飞轮储能系统轴承振动状态,将加速度传感器分别安装在轴承的垂直径向和轴向,进行9 000 r/min的升速振动检测,采样频率为2.5 KHz,而共振频率一般在5～500 Hz,生成500 Hz的频谱瀑布图,如图16所示。

图15 飞轮储能系统温度转速表Fig.15 Temperature and tachometer of flywheel energy storage system

图16 频谱瀑布图Fig.16 Spectral cascade diagram

在飞轮储能系统运行中,可以看到在1 800～2 300 r/min(图16中黑框区域)出现稍大的振动,其中在2 150 r/min振动最为强烈。

临界转速与有限元分析稍有出入,综合分析为储能系统壳体与飞轮轴系装配及轴承系数稍有不同所引起的现象,属于正常现象,有限元法仍然可以提供较理想的技术参考和理论分析。

5 结 论

通过有限元软件对飞轮轴系建模和仿真分析,以及全周期多次检测飞轮储能系统运行状态,可知仿真临界转速分析结果31 Hz与实际临界转速 2 150 r/min相近,实际转速52 Hz附近未发现明显的振动情况,综合分析为加工精度及装配问题所致。该研究结果对精确了解飞轮轴系的临界转速和运行状态以及轴承的刚度对飞轮转子稳定性的影响有非常重要的参考意义。