CN105510023A - 基于散度指标的变工况风电行星齿轮箱故障诊断方法 - Google Patents

基于散度指标的变工况风电行星齿轮箱故障诊断方法 Download PDF

Info

Publication number
CN105510023A
CN105510023A CN201510831633.1A CN201510831633A CN105510023A CN 105510023 A CN105510023 A CN 105510023A CN 201510831633 A CN201510831633 A CN 201510831633A CN 105510023 A CN105510023 A CN 105510023A
Authority
CN
China
Prior art keywords
fault
frequency
gear
divergence
rank
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201510831633.1A
Other languages
English (en)
Other versions
CN105510023B (zh
Inventor
吴冠宇
王方胜
滕海刚
陈国伟
史昌明
肖冰
王瑞
安国军
卢博伦
陈禹
赵海
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
State Grid Inner Mongolia East Power Integrated Energy Service Co ltd
State Grid Corp of China SGCC
Electric Power Research Institute of State Grid Eastern Inner Mongolia Power Co Ltd
Original Assignee
State Grid East Inner Mongolia Electric Power Energy-Saving Service Co Ltd
State Grid Corp of China SGCC
Electric Power Research Institute of State Grid Eastern Inner Mongolia Power Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by State Grid East Inner Mongolia Electric Power Energy-Saving Service Co Ltd, State Grid Corp of China SGCC, Electric Power Research Institute of State Grid Eastern Inner Mongolia Power Co Ltd filed Critical State Grid East Inner Mongolia Electric Power Energy-Saving Service Co Ltd
Priority to CN201510831633.1A priority Critical patent/CN105510023B/zh
Publication of CN105510023A publication Critical patent/CN105510023A/zh
Application granted granted Critical
Publication of CN105510023B publication Critical patent/CN105510023B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Landscapes

  • Testing Of Devices, Machine Parts, Or Other Structures Thereof (AREA)

Abstract

本发明属于旋转机械故障诊断技术领域,尤其涉及一种基于散度指标的变工况风电行星齿轮箱故障诊断方法,特别适用于变工况风电行星齿轮箱的故障诊断领域。本发明操作步骤如下:根据阶比重采样技术,将变工况风电行星齿轮箱传感器所采集的振动信号进行预处理,将非线性、非平稳的时域信号转化为具有平稳性的角域信号;行星齿轮箱不同于传统定轴齿轮箱,针对其结构特点及诊断的难度,将行星齿轮箱的故障分级进行诊断;提取故障特征集合;故障诊断参数;实验验证。本发明可避免振动信号非平稳的特点,有效清晰的识别故障特征阶比;J-散度和KL-散度均可表征行星齿轮箱故障发生位置及类型;在发生故障时,可相互辅助进行故障诊断,具有较高的灵敏度。

Description

基于散度指标的变工况风电行星齿轮箱故障诊断方法
技术领域
本发明属于旋转机械故障诊断技术领域,尤其涉及一种基于散度指标的变工况风电行星齿轮箱故障诊断方法,特别适用于变工况风电行星齿轮箱的故障诊断领域。
背景技术
行星齿轮箱广泛的应用于风力发电机组中,在实际运行中,行星齿轮箱不仅承受动态的重载负荷的影响,同时运行工况频繁的变化,使得其在运行过程中极其容易发生故障,其中太阳轮、行星轮、齿圈等关键零部件容易出现损伤故障对齿轮箱的影响非常巨大。因此,对其开展检测诊断对于保障风力发电机组安全高效稳定运行至关重要。但是,目前提出的时间同步平均、包络解调、倒谱、小波变换、Hilbert-Huang变换等多种诊断方法主要针对传统的定轴齿轮箱,行星齿轮箱的独特结构和运动特点使得其振动信号比传统的定轴齿轮箱更为复杂,故障诊断的难度较大,不可照搬直接应用。
行星齿轮箱不同于各个齿轮以其固定的中心轴旋转的定轴传动齿轮箱。行星齿轮传动系统由太阳轮、多个行星轮、内齿圈和行星架等结构组成。通常内齿圈固定不动,太阳轮绕自身的中心轴旋转,而几个行星轮不仅绕各自的中心轴自转,而且围绕太阳轮的中心轴公转,同时与太阳轮和内齿圈啮合,其齿轮运动的典型复合运动,使其振动响应比定轴传动齿轮箱更为复杂,因此,对于行星齿轮进行相应的故障诊断将具有更大的难度。
发明内容
针对现有技术中存在的不足,本发明提供一种基于散度指标的变工况风电行星齿轮箱故障诊断方法,是针对变工况下风电行星齿轮箱,以散度指标为特征参数的故障诊断方法,确切地说是一种可以同时实现行星齿轮箱故障模式的识别以及故障严重程度量化的故障诊断方法。目的是通过对振动信号进行阶比重采样,避免变工况对振动信号造成不平稳的影响,使得经过频谱分析后的频谱图不受变工况的影响。通过故障机理分析,确定行星齿轮箱的故障特征集合,计算故障样本与正常标准样本之间的散度值,通过观察不同故障特征集合所对应的散度值的变化情况,确定行星齿轮箱的故障模式以及故障的严重程度,本方法具有较高的灵敏度,可以避免查看复杂的频谱图形,减轻风电运维人员的工作难度。
为了达到上述发明目的,本发明的技术方案是这样实现的:
基于散度指标的变工况风电行星齿轮箱故障诊断方法,具体操作步骤如下:
(1)根据阶比重采样技术,将变工况风电行星齿轮箱传感器所采集的振动信号进行预处理,将非线性、非平稳的时域信号转化为具有平稳性的角域信号;是基于线性插值方法的非平稳振动时域信号的阶比重构技术,将等时间间隔采样的非平稳振动时域信号转化为具有平稳特性的角域振动信号,保证行星齿轮箱振动角域信号的整周期性;EMD经验模态分解方法根据信号的局部时变特性,自适应的将任意一个复杂信号分解为一系列分量,通过相关系数法则对信号进行重构,剔除原始信号中的干扰成分;
(2)行星齿轮箱不同于传统定轴齿轮箱,针对其结构的特点及诊断的难度,将行星齿轮箱的故障分级进行诊断;行星轮系的故障分为两类:分布故障和局部故障;对行星轮系的分布故障和局部故障的特征频率的进行分析计算,形成一个频率集合,并在阶比重采样的技术下,频率转化为阶比,相应的故障特征阶比不会随工况的变化而变化,形成固定的故障特征集合;
(3)提取故障特征集合;以行星齿轮箱为研究对象,按齿轮级数把行星齿轮箱分为三级:一级行星轮系、二级行星轮系和平行级;并按故障模式总体分为分布故障和局部故障两类,最终把行星齿轮箱的故障特征集合分为5个子集合,此时把平行级齿轮故障集合归结为一个子集合,由此实现对行星齿轮箱的分级分类诊断;
(4)故障诊断参数;由J-散度和KL-散度两个散度值的计算过程可以看出,两个散度值可以计算两个样本之间的差异程度;根据行星齿轮箱在处于正常状态和故障状态时其故障特征阶比所对应的幅值会发生变化,计算步骤(3)中得到的5个子集合中故障特征阶比所对应的幅值之间的散度值变化,即可实现对行星齿轮箱的故障诊断;可以说明散度值可以充分作为行星齿轮箱故障诊断的特征参数;
(5)实验验证;以行星齿轮箱处于正常运行状态时的振动数据为正常标准样本,计算不同级数、不同故障模式下的散度指标,通过观察各故障特征集合所对应的散度指标值的变化情况,实现对行星齿轮箱故障模式以及严重程度的识别;即利用J-散度和KL-散度,通过对风电行星齿轮箱处于不同状态下的故障特征集合所对应的幅值进行计算,一次性实现了风电行星齿轮箱故障模式的识别以及故障严重程度的量化,避免故障诊断分析过程中的重复操作;通过对行星齿轮箱不同运行状态下的振动数据按步骤(1)、(3)所述计算散度指标,发现散度指标J-散度和KL-散度可以作为复杂结构行星齿轮箱的故障诊断参数,并最终总结出针对行星齿轮箱的故障诊断流程。
根据步骤(1)所述的根据阶比重采样技术,将变工况风电行星齿轮箱传感器所采集的振动信号进行预处理;风电机组行星齿轮箱处于变转速、变工况的工作环境下,其采集的振动信号为非平稳信号,如直接进行频谱分析,很难得到清晰的频谱图,这对齿轮箱的故障诊断产生很大的困难,为了得到清晰正确的频谱图,采用阶比重采样技术对振动信号进行角域重采样,得到的频谱图中的阶比固定不变,便于对振动信号的分析;
阶比重采样技术的核心在于获得相对参考轴的恒定角增量采样数据,因此需要能准确获得阶次采样的时刻及相应的基准转速,即实现阶次跟踪;常见的阶次跟踪方法有硬件阶次跟踪法、计算阶次跟踪法和基于瞬时频率估计的阶次跟踪法;本发明采用计算阶次跟踪法,实现振动信号的重采样计算;
实际的齿轮箱振动信号一般情况下都含有多种干扰成分,这就使得其故障特征的提取变得比较困难;经验模态分解(EmpiricalModeDecomposition,EMD)可以根据信号的局部时变特性,自适应的将任意一个复杂信号分解为一系列分量,通过相关系数法则对信号进行重构,剔除原始信号中的干扰成分;
当齿轮发生故障时,其振动信号都具有调制特征,从信号中提取调制信息,并分析其强度和频次就可以判断故障的部位和损伤程度;信号包络谱,可反映周期性的冲击及其剧烈程度。
根据步骤(2)所述,行星齿轮箱不同于传统定轴齿轮箱,针对其结构的特点及诊断的难度,将行星齿轮箱的故障分级进行诊断;实现行星轮故障特征频率计算;风电机组齿轮增速箱结构多样,传动比大,为减小齿轮箱的尺寸,一般为行星齿轮结构,本发明对某一风电机组行星齿轮箱进行分析,其风电行星齿轮箱由两级行星轮、一级平行齿轮构成;
风电机组行星齿轮箱的两级行星轮系和平行级齿轮结构;
对于行星轮系和平行级齿轮,其故障可分为局部故障和分布式故障;在单级行星齿轮箱中,太阳轮–行星轮和行星轮–齿圈两种啮合副的啮合频率相同;通常齿圈固定不动,太阳轮、行星轮和行星架旋转,在这种情况下,啮合频率:
fm=fcZr=(fs (r)-fc)Zs(1);
式中:Zr和Zs分别为齿圈和太阳轮的齿数;fm为啮合频率;fc为行星架的旋转频率;fs (r)为太阳轮的绝对旋转频率;
太星轮局部故障特征频率为:
式中:fm为啮合频率;Zs为太阳轮齿数;N为行星轮数量,fs为太阳轮局部故障特征频率;行星轮局部故障特征频率为:
式中:fm为啮合频率;Zp为行星轮齿数;fp为行星轮局部故障特征频率;齿圈局部故障特征频率为:
式中:fm为啮合频率;fr为齿圈局部故障特征频率;N为行星轮数量,Zr为齿圈齿数。
行星齿轮箱中各种齿轮的分布式故障特征频率等于齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率。已知行星齿轮箱的啮合频率fm和某个齿轮的齿数Zg,则该齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率:
fg=fm/Zg(5);
fg该齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率,Zg为某个齿轮的齿数;则太阳轮、行星轮和齿圈分布式故障的特征频率分别为:
fs'=fm/Zs(6);
fp'=fm/Zp(7);
fr'=fm/Zr(8);
式中:fm为啮合频率;fs'、fp'、fr'太阳轮、行星轮和齿圈分布式故障的特征频率;Zs为太阳轮齿数;Zp为行星轮的齿数;Zr为齿圈齿数;
以与主轴相连接的一级行星轮系的行星架为参考转速,对行星齿轮箱中各级的各个齿轮的局部故障和分布故障的特征阶比进行计算。
根据步骤(3)所述的行星齿轮箱故障特征量的提取;行星轮系局部故障,可分为太阳轮局部故障、行星轮局部故障和内齿圈局部故障;
对于太阳轮局部故障振动信号,在包络谱中,峰值出现在太阳轮的局部故障特征频率fs、太阳轮的绝对旋转频率fs (r)、以及它们的组合fs±fs (r)等位置处;若考虑太阳轮局部故障特征频率的倍频和太阳轮绝对旋转频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在太阳轮局部故障特征频率及其倍频nfs、太阳轮的绝对旋转频率及其倍频mfs (r)、及其组合nfs±mfs (r)等位置处;
对于行星轮局部故障,在包络谱中,峰值出现在行星轮局部故障特征频率fp、行星架的旋转频率fc、以及它们的fp±fc组合等位置处,若考虑行星轮局部故障特征频率的倍频和行星架旋转频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在行星轮局部故障特征频率及其倍频nfp、行星架的旋转频率及其倍频mfc、以及它们的组合nfp±mfc等位置处;
对于齿圈局部故障,在包络谱中,峰值出现在齿圈局部故障特征频率fr位置处,若考虑齿圈局部故障特征频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在齿圈局部故障特征频率及其倍频nfr位置处;
行星轮系发生分布故障,在包络谱中,峰值出现在齿轮分布式故障特征频率fg、行星轮通过频率Nfc、及其组合fg±Nfc位置处;若考虑齿轮分布式故障特征频率的倍频和行星轮通过频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在齿轮分布式故障特征频率及其倍频nfg(n为正整数)、行星轮通过频率及其倍频mNfc(m为正整数)、以及它们的组合nfg±mNfc位置处;
定轴齿轮发生包括齿根部有较大裂纹、局部齿面磨损、轮齿折断、局部齿形误差局部故障时,其振动信号波形是以齿轮旋转频率为周期的冲击脉冲,在频率域表现为包含旋转频率的各次谐波mfr(m=1,2,···)、各阶啮合频率nfm(n=1,2,···)以及以故障齿轮的旋转频率为间隔的边频nfm±mfr(n,m=1,2···);
定轴齿轮发生分布故障时,其频域特征表现为啮合频率及其谐波分量nfm(n=1,2,···)在频谱图上的位置保持不变,但其幅值大小发生改变,而且高次谐波幅值相对增大较多;分析时,要分析3个以上谐波幅值的变化,从而从频谱上检测出这种特征;
根据上面所述,在阶比谱和阶比包络谱中选取能表征故障类型的特征频率处的幅值为特征量,表示为行星轮系局部故障行星轮系分布式故障平行级啮合齿轮故障n=1,2,3,其中平行级啮合齿轮故障归为一个特征向量,此向量从阶比图中提取,行星轮系的特征向量从包络谱中提取。
根据步骤(4)所述的故障诊断参数;即故障模式及严重程度的识别;J-散度(J-divergence)是一种谱距离,作为一种状态识别的指标,能很好地反映两个信号的相似程度,克服时域分析中的相位问题,很显然,同一信号的J-散度为零;
式(9)中,S为样本信号的幅值谱;τ为待检信号的幅值谱;J(s,τ)为两者之间的J-散度,N为信号幅值谱中幅值的个数,i为信号幅值谱中幅值的序列。
J-散度作为一种状态识别的指标,它能很好地反映两个信号的相似程度,克服时域分析中的相位问题,很显然,同一信号的J-散度为零;
KL-散度(Kullback-Lciblerdivergence,KLD),用来度量分布P和Q之间的差异性,典型情况下,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;
对离散分布,P和Q的KL-散度定义为:
式(10)中,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;Dkl为数据P和Q的KL-散度,i为数据的序列号;P(i)为序列号i所对应的P数据点;Q(i)为序列号i所对应的Q数据点。
也有人称其为KL距离,但是它并不是严格的距离概念,其不满足三角不等式;所以,把它变为对称形式:
Dkls(P||Q)=[Dkl(P||Q)+Dkl(Q||P)]/2(11);
式(11)中,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;D'kl数据点P、Q之间的对称散度值。
当齿轮箱正常时,本文所选择的故障特征阶比的幅值为0或很小,而当发生故障时,其所对应的故障阶比的幅值发生很大的变化,这使得采用采用散度指标的算法得以实现,对上述所选特征量进行J-散度和KL-散度计算,计算方法:首先收集风电行星齿轮箱的正常状态的样本,表示正常情况下的标准样本,分别表示一级行星轮系、二级行星轮系和平行级齿轮故障样本,对每一样本进行阶比重采样、EMD重构信号以及进行阶比谱、Hilbert包络谱分析,寻找对应的故障特征阶比,形成对应故障特征阶比的幅值集合,最后计算待检样本和标准样本的故障特征集合之间的J-散度和KL-散度,定位故障位置以及故障模式,对行星齿轮箱的故障实现完全的诊断;
在时域信号进行FFT变换时的分辨率的原理推出,角域信号进行FFT的分辨率为2π/θ,其中θ为角域信号的长度;信号为matlab仿真信号,信号中忽略齿轮箱中轴承以及各齿轮相互之间的影响;
由此,仿真行星齿轮箱各级处于正常工况时,其振动信号模型为:
x1(t)=A·[1-cos(2π·3fc1·t)]cos(2π·fm1·t+θ1)
(12);
+B·[1-cos(2π·3fc2·t)]·cos(2π·fm2·t+θ2)+C·cos(2π·fm3·t+θ3)
行星齿轮箱一级行星轮系太阳轮发生局部故障时,其振动信号模型为:
行星齿轮箱一级行星轮系太阳轮发生分布故障时,其振动信号模型为:
式(12)~(14)中::x1(t)、x2(t)、x3(t)为行星齿轮箱处于正常一级行星轮系太阳轮发生局部故障、分布故障时的振动信号序列;t为时间序列;θ1、θ2、θ3、φ、为初始相位;fm1、fm2、fm3为各级啮合频率;fc1、fc2为一、二级行星架的旋转频率;为一级太阳轮的绝对旋转频率;fs1、fs1'为一级太阳轮发生局部故障和分布故障时的特征频率;A、B、C为无量纲常数,行星齿轮箱各个状态时值会有所不同;各个振动信号采用频率为8192HZ。
根据步骤(5)所述的实验验证;是在步骤1、步骤4的基础上,以行星齿轮箱的正常、行星齿轮箱一级行星轮系太阳轮局部故障、行星齿轮箱一级行星轮系太阳轮分布故障时的振动信号为例来进行分析。
所述的阶比重采样方法有效避免振动信号非平稳特点引发的频率冗杂、难以分析的特点,通过EMD分解消噪后重构后的阶比谱和包络谱清晰,容易找到对应的特征阶次,避免噪声对有效信号的干扰。
所述的行星齿轮箱一级行星轮系太阳轮发生局部故障采集的振动信号经阶比重采样和EMD分解重构后的阶比谱;阶比谱中一级行星轮系啮合频率附近出现边频带,与正常时对比发生很明显的变化,同时对比正常与故障的阶比包络谱,故障时的调制成分发生很大的变化,与正常标准样本进行散度指标的计算,正常时故障特征阶比集合所对应的故障阶比的幅值很小,甚至为0,而发生故障时其故障特征阶比的幅值变化很大,通过与正常样本进行散度指标的计算,确定散度指标越大,发生故障的可能性越大,确定故障发生的粗略位置以及故障的严重程度。
所述的行星齿轮箱一级行星轮系太阳轮,发生分布故障时所采集的振动信号经阶比重采样和EMD分解重构后进行频谱分析和包络谱分析;阶比谱中一级行星轮系啮合频率附近出现边频带,与正常时对比发生很明显的变化,同时对比分布故障与局部故障的阶比谱不同,观察阶比包络谱,故障信号中存在复杂的调频信息,计算故障特征向量对应的幅值的散度指标,可以实现对分布故障进行有效的诊断,同时也可以有效的识别出局部故障和分布故障;
以风电行星齿轮箱为研究对象,对行星齿轮箱进行拆分为一级行星轮系、二级行星轮系、平行级齿轮三部分,同时行星轮系故障分为局部故障、分布式故障,暂不考虑轴承的影响,同时相比较基于散度指标的轴承故障诊断,其所采用的对样本之间进行散度指标的计算,来观察之间的相似性,而行星齿轮箱结构复杂,缺少特定的故障样本数据,对故障位置、故障模式的确定具有更大的难度,决定采用故障样本只与正常标准样进行散度指标计算,输入散度指标的故障特征向量分别为一级行星轮系A11、A12,二级行轮系A21、A22,平行级齿轮A3,通过计算观察某一故障特征向量较其他发生很大变化来确定故障发生的位置以及故障模式,对行星齿轮箱正常标准样本N与行星齿轮箱一级行星轮系太阳轮局部故障、分布故障样本的故障特征向量进行计算,对比分析,找出行星齿轮箱故障发生的哪一级以及哪一种故障模式,即此样本故障发生位置和模式,再通过对A11进行进一步计算,完全推出故障发生的位置,实现对故障的完全诊断;
由其他位置发生故障时各个故障样本与标准样本的J-散度、KL-散度值的变化可以发现,故障样本与正常标准样本的散度指标计算值越大,就可以确定此样本为在此故障特征向量所包含的故障模式集合中的一类,这样可以基本锁定齿轮箱发生故障的位置和类型,通过对此故障特征集合进行细致区别计算,完全诊断出故障位置和模式。
所述的散度指标,经过计算可能够诊断出行星齿轮箱的故障模式及其位置,这使得风电机组运维人员可以完全避免查看复杂的频谱图,通过观察一些指标的变化可以诊断出齿轮箱的运行状态,同时散度指标的大小也可以衡量故障的严重程度。
本发明的优点及有益效果如下:
本发明提供一种针对变工况下风电行星齿轮箱,以散度指标为特征参数的故障诊断方法,通过对风电机组行星齿轮箱正常以及发生局部故障、分布式故障的振动信号的分析结果表明:基于阶比谱、阶比包络谱的分析,可以避免振动信号非平稳的特点,有效清晰的识别故障特征阶比;J-散度和KL-散度均可表征行星齿轮箱故障发生位置及类型;在发生故障时,KL-散度、J-散度可以相互辅助进行故障诊断,具有较高的灵敏度,基于散度指标的故障诊断方法可以避免查看复杂的频谱图形,此故障诊断方法对风电场风力发电机行星齿轮箱齿轮箱的故障诊断具有一定的参考价值。
下面结合附图和具体实施例对本发明作进一步详细的说明。
附图说明
图1为本发明中行星齿轮箱结构简略图;
图2为本发明中行星齿轮箱正常时的角域信号以及阶比谱、阶比包络谱;
图3为本发明中行星齿轮箱一级太阳轮局部故障角域图以及阶比谱、阶比包络谱;
图4为本发明中行星齿轮箱一级太阳轮分布故障角域图以及阶比谱、阶比包络谱;
图5基于阶比分析、EMD和散度指标的行星齿轮箱故障诊断方法流程图
具体实施方式
本发明提供一种针对变工况下风电行星齿轮箱,以散度指标为特征参数的故障诊断方法,适用于变工况下的风电机组星型齿轮箱的故障诊断方法,具体操作步骤如下:
(1)根据阶比重采样技术,将处于变工况的风电机组行星齿轮箱传感器所采集的振动信号进行预处理,将非线性、非平稳的时域信号转化为具有平稳性的角域信号,避免了使用硬件方式实现等角度采样的昂贵成本。本步骤是基于线性插值方法的非平稳振动时域信号的阶比重构技术,将等时间间隔采样的非平稳振动时域信号转化为具有平稳特性的角域振动信号,保证了行星齿轮箱振动角域信号的整周期性。EMD经验模态分解方法根据信号的局部时变特性,自适应的将任意一个复杂信号分解为一系列分量,通过相关系数法则对信号进行重构,剔除原始信号中的干扰成分。
(2)行星齿轮箱不同于传统定轴齿轮箱,针对其结构的特点及诊断的难度,将行星齿轮箱的故障分级进行诊断。行星轮系的故障分为两类:分布故障和局部故障。对行星轮系的分布故障和局部故障的特征频率的进行分析计算,形成一个频率集合,并在阶比重采样的技术下,频率转化为阶比,相应的故障特征阶比不会随工况的变化而变化,形成固定的故障特征集合。
(3)提取故障特征集合。以行星齿轮箱为研究对象,按齿轮级数把行星齿轮箱分为三级:一级行星轮系、二级行星轮系和平行级。并按故障模式总体分为分布故障和局部故障两类,最终把行星齿轮箱的故障特征集合分为5个子集合,此时把平行级齿轮故障集合归结为一个子集合,由此实现对行星齿轮箱的分级分类诊断。
(4)故障诊断参数。由J-散度和KL-散度两个散度值的计算过程可以看出,两个散度值可以计算两个样本之间的差异程度。根据行星齿轮箱在处于正常状态和故障状态时其故障特征阶比所对应的幅值会发生变化,计算步骤(3)中得到的5个子集合中故障特征阶比所对应的幅值之间的散度值变化,即可实现对行星齿轮箱的故障诊断。可以说明散度值可以充分作为行星齿轮箱故障诊断的特征参数。
(5)实验验证。以行星齿轮箱处于正常运行状态时的振动数据为正常标准样本,计算不同级数、不同故障模式下的散度指标,通过观察各故障特征集合所对应的散度指标值的变化情况,实现对行星齿轮箱故障模式以及严重程度的识别。即利用J-散度和KL-散度,通过对风电行星齿轮箱处于不同状态下的故障特征集合所对应的幅值进行计算,一次性实现了风电行星齿轮箱故障模式的识别以及故障严重程度的量化,避免故障诊断分析过程中的重复操作。通过对行星齿轮箱不同运行状态下的振动数据按步骤(1)、(3)所述计算散度指标,发现散度指标J-散度和KL-散度可以作为复杂结构行星齿轮箱的故障诊断参数。并最终总结出针对行星齿轮箱的故障诊断流程,如图5所示。
下面针对本发明方法的具体操作步骤做一详细的分析和解释说明如下,:
步骤1,齿行星轮箱振动信号预处理。风电机组行星齿轮箱处于变转速、变工况的工作环境下,其采集的振动信号为非平稳信号,如直接进行频谱分析,很难得到清晰的频谱图,这对齿轮箱的故障诊断产生很大的困难,为了得到清晰正确的频谱图,采用阶比重采样技术对振动信号进行角域重采样,得到的频谱图中的阶比固定不变,便于对振动信号的分析。
阶比重采样技术的核心在于获得相对参考轴的恒定角增量采样数据,因此需要能准确获得阶次采样的时刻及相应的基准转速,即实现阶次跟踪。常见的阶次跟踪方法有硬件阶次跟踪法、计算阶次跟踪法和基于瞬时频率估计的阶次跟踪法等。本发明采用计算阶次跟踪法,实现振动信号的重采样计算。
实际的齿轮箱振动信号一般情况下都含有多种干扰成分,这就使得其故障特征的提取变得比较困难。经验模态分解(EmpiricalModeDecomposition,EMD)可以根据信号的局部时变特性,自适应的将任意一个复杂信号分解为一系列分量,通过相关系数法则对信号进行重构,剔除原始信号中的干扰成分。
当齿轮发生故障时,其振动信号都具有调制特征,从信号中提取调制信息,并分析其强度和频次就可以判断故障的部位和损伤程度。信号包络谱,可反映周期性的冲击及其剧烈程度。
步骤2,故障特征频率计算。风电机组齿轮增速箱结构多样,传动比大,为减小齿轮箱的尺寸,一般为行星齿轮结构,本发明对某一风电机组行星齿轮箱进行分析,其行星齿轮箱由两级行星轮、一级平行齿轮构成,其结构如图1所示。
风电机组行星齿轮箱的两级行星轮系和平行级齿轮结构,其结构参数如表1所示。
风电行星齿轮箱一般分为行星轮系、平行级齿轮,对于行星轮系和平行级齿轮,其故障可分为局部故障和分布式故障。在单级行星齿轮箱中,太阳轮–行星轮和行星轮–齿圈两种啮合副的啮合频率相同。通常齿圈固定不动,太阳轮、行星轮和行星架旋转,在这种情况下,啮合频率:
fm=fcZr=(fs (r)-fc)Zs(1);
式中:Zr和Zs分别为齿圈和太阳轮的齿数;fm为啮合频率;fc为行星架的旋转频率;fs (r)为太阳轮的绝对旋转频率。
太阳轮局部故障特征频率为:
式中:fm为啮合频率;Zs为太阳轮齿数;N为行星轮数量,fs为太阳轮局部故障特征频率;行星轮局部故障特征频率为:
式中:fm为啮合频率;Zp为行星轮齿数;fp为行星轮局部故障特征频率;齿圈局部故障特征频率为:
式中:fm为啮合频率;fr为齿圈局部故障特征频率;N为行星轮数量,Zr为齿圈齿数。
行星齿轮箱中各种齿轮的分布式故障特征频率等于齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率。已知行星齿轮箱的啮合频率fm和某个齿轮的齿数Zg,则该齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率:
fg=fm/Zg(5);
fg该齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率,Zg为某个齿轮的齿数;则太阳轮、行星轮和齿圈分布式故障的特征频率分别为:
fs'=fm/Zs(6);
fp'=fm/Zp(7);
fr'=fm/Zr(8);
式中:fm为啮合频率;fs'、fp'、fr'太阳轮、行星轮和齿圈分布式故障的特征频率;Zs为太阳轮齿数;Zp为行星轮的齿数;Zr为齿圈齿数。
以与主轴相连接的一级行星轮系的行星架为参考转速,对行星齿轮箱中各级的各个齿轮的局部故障和分布故障的特征阶比进行计算,如表2所示。
步骤3,行星齿轮箱故障特征量的提取。行星轮系局部故障,可分为太阳轮局部故障、行星轮局部故障和内齿圈局部故障。对于太阳轮局部故障振动信号,在包络谱中,峰值出现在太阳轮的局部故障特征频率fs、太阳轮的绝对旋转频率fs (r)、以及它们的组合fs±fs (r)等位置处。若考虑太阳轮局部故障特征频率的倍频和太阳轮绝对旋转频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在太阳轮局部故障特征频率及其倍频nfs、太阳轮的绝对旋转频率及其倍频mfs (r)、及其组合nfs±mfs (r)等位置处。对于行星轮局部故障,在包络谱中,峰值出现在行星轮局部故障特征频率fp、行星架的旋转频率fc、以及它们的fp±fc组合等位置处,若考虑行星轮局部故障特征频率的倍频和行星架旋转频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在行星轮局部故障特征频率及其倍频nfp、行星架的旋转频率及其倍频mfc、以及它们的组合nfp±mfc等位置处。对于齿圈局部故障,在包络谱中,峰值出现在齿圈局部故障特征频率fr位置处,若考虑齿圈局部故障特征频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在齿圈局部故障特征频率及其倍频nfr位置处。
行星轮系发生分布故障,在包络谱中,峰值出现在齿轮分布式故障特征频率fg、行星轮通过频率Nfc、及其组合fg±Nfc等位置处。若考虑齿轮分布式故障特征频率的倍频和行星轮通过频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在齿轮分布式故障特征频率及其倍频nfg(n为正整数)、行星轮通过频率及其倍频mNfc(m为正整数)、以及它们的组合nfg±mNfc等位置处。
定轴齿轮发生包括齿根部有较大裂纹、局部齿面磨损、轮齿折断、局部齿形误差等局部故障时,其振动信号波形是以齿轮旋转频率为周期的冲击脉冲,在频率域表现为包含旋转频率的各次谐波mfr(m=1,2,···)、各阶啮合频率nfm(n=1,2,···)以及以故障齿轮的旋转频率为间隔的边频nfm±mfr(n,m=1,2···)等。
定轴齿轮发生分布故障时,其频域特征表现为啮合频率及其谐波分量nfm(n=1,2,···)在频谱图上的位置保持不变,但其幅值大小发生改变,而且高次谐波幅值相对增大较多。分析时,要分析3个以上谐波幅值的变化,从而从频谱上检测出这种特征。
根据上面所述,在阶比谱和阶比包络谱中选取能表征故障类型的特征频率处的幅值为特征量,表示为行星轮系局部故障行星轮系分布式故障平行级啮合齿轮故障n=1,2,3,其中平行级啮合齿轮故障归为一个特征向量,此向量从阶比图中提取,行星轮系的特征向量从包络谱中提取。
步骤4,故障模式及严重程度的识别。J-散度(J-divergence)是一种谱距离,可以作为一种状态识别的指标,它能很好地反映两个信号的相似程度,克服时域分析中的相位问题,很显然,同一信号的J-散度为零。
式(9)中,S为样本信号的幅值谱;τ为待检信号的幅值谱;J(s,τ)为两者之间的J-散度,N为信号幅值谱中幅值的个数,i为信号幅值谱中幅值的序列。
J-散度作为一种状态识别的指标,它能很好地反映两个信号的相似程度,克服时域分析中的相位问题,很显然,同一信号的J-散度为零。
KL-散度(Kullback-Lciblerdivergence,KLD),用来度量分布P和Q之间的差异性,典型情况下,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布。
对离散分布,P和Q的KL-散度定义为:
式(10)中,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;Dkl为数据P和Q的KL-散度,i为数据的序列号;P(i)为序列号i所对应的P数据点;Q(i)为序列号i所对应的Q数据点。
也有人称其为KL距离,但是它并不是严格的距离概念,其不满足三角不等式。所以,把它变为对称形式:
D'kl(P||Q)=[Dkl(P||Q)+Dkl(Q||P)]/2(11);
式(11)中,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;D'kl数据点P、Q之间的对称散度值。
当齿轮箱正常时,本文所选择的故障特征阶比的幅值为0或很小,而当发生故障时,其所对应的故障阶比的幅值发生很大的变化,这使得采用采用散度指标的算法得以实现,对上述所选特征量进行J-散度和KL-散度计算,计算方法:首先收集风电行星齿轮箱的正常状态的样本,表示正常情况下的标准样本,分别表示一级行星轮系、二级行星轮系和平行级齿轮故障样本,对每一样本进行阶比重采样、EMD重构信号以及进行阶比谱、Hilbert包络谱分析,寻找对应的故障特征阶比,形成对应故障特征阶比的幅值集合,最后计算待检样本和标准样本的故障特征集合之间的J-散度和KL-散度,定位故障位置以及故障模式,对行星齿轮箱的故障实现完全的诊断。
本发明在时域信号进行FFT变换时的分辨率的原理推出,角域信号进行FFT的分辨率为2π/θ,其中θ为角域信号的长度。本发明的信号为matlab仿真信号,信号中忽略齿轮箱中轴承以及各齿轮之间的影响,假设各齿轮箱传动级之间的振动影响不存在,对本发明的方法进行仿真验证针对行星轮系处于正常、分布故障、局部故障时的振动信号模型参考文献,这里就不再赘述。
由此,可以仿真行星齿轮箱各级处于正常工况时,其振动信号模型为:
x1(t)=A·[1-cos(2π·3fc1·t)]cos(2π·fm1·t+θ1)
(12);
+B·[1-cos(2π·3fc2·t)]·cos(2π·fm2·t+θ2)+C·cos(2π·fm3·t+θ3)
行星齿轮箱一级行星轮系太阳轮发生局部故障时,其振动信号模型为:
行星齿轮箱一级行星轮系太阳轮发生分布故障时,其振动信号模型为:
式(12)~(14)中:x1(t)、x2(t)、x3(t)为行星齿轮箱处于正常一级行星轮系太阳轮发生局部故障、分布故障时的振动信号序列;t为时间序列;θ1、θ2、θ3、φ、为初始相位;fm1、fm2、fm3为各级啮合频率;fc1、fc2为一、二级行星架的旋转频率;为一级太阳轮的绝对旋转频率;fs1、fs1'为一级太阳轮发生局部故障和分布故障时的特征频率;A、B、C为无量纲常数,行星齿轮箱各个状态时值会有所不同,这里就不再详述。各个振动信号采用频率为8192HZ。
步骤5,在步骤1、步骤4的基础上,以行星齿轮箱的正常、行星齿轮箱一级行星轮系太阳轮局部故障、行星齿轮箱一级行星轮系太阳轮分布故障时的振动信号为例来进行分析。
图2是行星齿轮箱正常振动信号经过阶比角域重采样以及EMD重构后的平稳角域信号图,其信号长度为,进行阶比谱、阶比包络谱分析,可以看出,采用阶比重采样技术可以有效避免振动信号非平稳特点引发的频率冗杂、难以分析的特点,通过EMD分解消噪后重构后的阶比谱和包络谱清晰,容易找到对应的特征阶次,避免噪声对有效信号的干扰。
图3为行星齿轮箱一级行星轮系太阳轮发生局部故障采集的振动信号经阶比重采样和EMD分解重构后的阶比谱。可以看出,阶比谱中一级行星轮系啮合频率附近出现边频带,与正常时对比发生很明显的变化,同时对比正常与故障的阶比包络谱,故障时的调制成分发生很大的变化,与正常标准样本进行散度指标的计算,正常时故障特征阶比集合所对应的故障阶比的幅值很小,甚至为0,而发生故障时其故障特征阶比的幅值变化很大,通过与正常样本进行散度指标的计算,可以确定,散度指标越大,发生故障的可能性越大,可以确定故障发生的粗略位置以及故障的严重程度。
图4为行星齿轮箱一级行星轮系太阳轮发生分布故障时所采集的振动信号,该振动信号经阶比重采样和EMD分解重构后进行频谱分析和包络谱分析。由图可以看出,阶比谱中一级行星轮系啮合频率附近出现边频带,与正常时对比发生很明显的变化,同时对比分布故障与局部故障的阶比谱可以发现有很大的不同,观察阶比包络谱可以看出,故障信号中存在复杂的调频信息,计算故障特征向量对应的幅值的散度指标,可以实现对分布故障进行有效的诊断,同时也可以有效的识别出局部故障和分布故障。
以风电行星齿轮箱为研究对象,对行星齿轮箱进行拆分为一级行星轮系、二级行星轮系、平行级齿轮三部分,同时行星轮系故障分为局部故障、分布式故障,暂不考虑轴承的影响,同时相比较基于散度指标的轴承故障诊断,其所采用的对样本之间进行散度指标的计算,来观察之间的相似性,而行星齿轮箱结构复杂,缺少特定的故障样本数据,对故障位置、故障模式的确定具有更大的难度,决定采用故障样本只与正常标准样进行散度指标计算,输入散度指标的故障特征向量分别为一级行星轮系A11、A12,二级行轮系A21、A22,平行级齿轮A3,通过计算观察某一故障特征向量较其他发生很大变化来确定故障发生的位置以及故障模式,本发明对行星齿轮箱正常标准样本N与行星齿轮箱一级行星轮系太阳轮局部故障、分布故障样本的故障特征向量进行计算,对比分析,找出行星齿轮箱故障发生的哪一级以及哪一种故障模式,即此样本故障发生位置和模式,再通过对A11进行进一步计算,完全推出故障发生的位置,实现对故障的完全诊断,计算散度指标如表3样本S1、S7所示,表3是风电机行星齿轮箱各类故障模式下的故障样本与正常状态下的标准样本之间的J-散度和KL-散度值。
同时表3也罗列了其他位置发生故障时各个故障样本与标准样本的J-散度、KL-散度值的变化,可以发现故障样本与正常标准样本的散度指标计算值越大,可以确定此样本为在此故障特征向量所包含的故障模式集合,这样可以基本锁定齿轮箱发生故障的位置和类型,通过对此故障特征集合进行细致区别计算,可完全诊断出故障位置和模式,样本S1、S7在上述的基础上,样本S1、S7散度指标的细化计算为表4所示。
由表4可以看出,经过散度指标的计算可以完全的诊断出行星齿轮箱的故障模式及其位置,这使得风电机组运维人员可以完全避免查看复杂的频谱图,通过观察一些指标的变化可以诊断出齿轮箱的运行状态,诊断变得更加简单,同时散度指标的大小也可以衡量故障的严重程度。具体的行星齿轮箱故障中诊断流程如图5所示。
表1风电行星齿轮箱结构参数。
表2行星齿轮箱故障特征频率计算。
表3风电机行星齿轮箱故障样本与标准样本之间的J-散度和KL-散度。
表4样本S1、S7散度指标的细化计算。

Claims (10)

1.基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:具体操作步骤如下:
(1)根据阶比重采样技术,将变工况风电行星齿轮箱传感器所采集的振动信号进行预处理,将非线性、非平稳的时域信号转化为具有平稳性的角域信号;是基于线性插值方法的非平稳振动时域信号的阶比重构技术,将等时间间隔采样的非平稳振动时域信号转化为具有平稳特性的角域振动信号,保证行星齿轮箱振动角域信号的整周期性;EMD经验模态分解方法根据信号的局部时变特性,自适应的将任意一个复杂信号分解为一系列分量,通过相关系数法则对信号进行重构,剔除原始信号中的干扰成分;
(2)行星齿轮箱不同于传统定轴齿轮箱,针对其结构的特点及诊断的难度,将行星齿轮箱的故障分级进行诊断;行星轮系的故障分为两类:分布故障和局部故障;对行星轮系的分布故障和局部故障的特征频率的进行分析计算,形成一个频率集合,并在阶比重采样的技术下,频率转化为阶比,相应的故障特征阶比不会随工况的变化而变化,形成固定的故障特征集合;
(3)提取故障特征集合;以行星齿轮箱为研究对象,按齿轮级数把行星齿轮箱分为三级:一级行星轮系、二级行星轮系和平行级;并按故障模式总体分为分布故障和局部故障两类,最终把行星齿轮箱的故障特征集合分为5个子集合,此时把平行级齿轮故障集合归结为一个子集合,由此实现对行星齿轮箱的分级分类诊断;
(4)故障诊断参数;由J-散度和KL-散度两个散度值的计算过程可以看出,两个散度值可以计算两个样本之间的差异程度;根据行星齿轮箱在处于正常状态和故障状态时其故障特征阶比所对应的幅值会发生变化,计算步骤(3)中得到的5个子集合中故障特征阶比所对应的幅值之间的散度值变化,即可实现对行星齿轮箱的故障诊断;可以说明散度值可以充分作为行星齿轮箱故障诊断的特征参数;
(5)实验验证;以行星齿轮箱处于正常运行状态时的振动数据为正常标准样本,计算不同级数、不同故障模式下的散度指标,通过观察各故障特征集合所对应的散度指标值的变化情况,实现对行星齿轮箱故障模式以及严重程度的识别;即利用J-散度和KL-散度,通过对风电行星齿轮箱处于不同状态下的故障特征集合所对应的幅值进行计算,一次性实现了风电行星齿轮箱故障模式的识别以及故障严重程度的量化,避免故障诊断分析过程中的重复操作;通过对行星齿轮箱不同运行状态下的振动数据按步骤(1)、(3)所述计算散度指标,发现散度指标J-散度和KL-散度可以作为复杂结构行星齿轮箱的故障诊断参数,并最终总结出针对行星齿轮箱的故障诊断流程。
2.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:根据步骤(1)所述的根据阶比重采样技术,将变工况风电行星齿轮箱传感器所采集的振动信号进行预处理;风电机组行星齿轮箱处于变转速、变工况的工作环境下,其采集的振动信号为非平稳信号,如直接进行频谱分析,很难得到清晰的频谱图,这对齿轮箱的故障诊断产生很大的困难,为了得到清晰正确的频谱图,采用阶比重采样技术对振动信号进行角域重采样,得到的频谱图中的阶比固定不变,便于对振动信号的分析;
阶比重采样技术的核心在于获得相对参考轴的恒定角增量采样数据,因此需要能准确获得阶次采样的时刻及相应的基准转速,即实现阶次跟踪;常见的阶次跟踪方法有硬件阶次跟踪法、计算阶次跟踪法和基于瞬时频率估计的阶次跟踪法;本发明采用计算阶次跟踪法,实现振动信号的重采样计算;
实际的齿轮箱振动信号一般情况下都含有多种干扰成分,这就使得其故障特征的提取变得比较困难;经验模态分解(EmpiricalModeDecomposition,EMD)可以根据信号的局部时变特性,自适应的将任意一个复杂信号分解为一系列分量,通过相关系数法则对信号进行重构,剔除原始信号中的干扰成分;
当齿轮发生故障时,其振动信号都具有调制特征,从信号中提取调制信息,并分析其强度和频次就可以判断故障的部位和损伤程度;信号包络谱,可反映周期性的冲击及其剧烈程度。
3.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:根据步骤(2)所述,行星齿轮箱不同于传统定轴齿轮箱,针对其结构的特点及诊断的难度,将行星齿轮箱的故障分级进行诊断;实现行星轮故障特征频率计算;风电机组齿轮增速箱结构多样,传动比大,为减小齿轮箱的尺寸,一般为行星齿轮结构,本发明对某一风电机组行星齿轮箱进行分析,其风电行星齿轮箱由两级行星轮、一级平行齿轮构成;
风电机组行星齿轮箱的两级行星轮系和平行级齿轮结构;
对于行星轮系和平行级齿轮,其故障可分为局部故障和分布式故障;在单级行星齿轮箱中,太阳轮–行星轮和行星轮–齿圈两种啮合副的啮合频率相同;通常齿圈固定不动,太阳轮、行星轮和行星架旋转,在这种情况下,啮合频率:
fm=fcZr=(fs (r)-fc)Zs(1);
式中:Zr和Zs分别为齿圈和太阳轮的齿数;fm为啮合频率;fc为行星架的旋转频率;fs (r)为太阳轮的绝对旋转频率;
太星轮局部故障特征频率为:
f s = f m Z s N - - - ( 2 ) ;
式中:fm为啮合频率;Zs为太阳轮齿数;N为行星轮数量,fs为太阳轮局部故障特征频率;行星轮局部故障特征频率为:
f p = 2 f m Z p - - - ( 3 ) ;
式中:fm为啮合频率;Zp为行星轮齿数;fp为行星轮局部故障特征频率;齿圈局部故障特征频率为:
f r = f m Z r N - - - ( 4 ) ;
式中:fm为啮合频率;fr为齿圈局部故障特征频率;N为行星轮数量,Zr为齿圈齿数;
行星齿轮箱中各种齿轮的分布式故障特征频率等于齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率;已知行星齿轮箱的啮合频率fm和某个齿轮的齿数Zg,则该齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率:
fg=fm/Zg(5);
fg该齿轮相对行星架(太阳轮和齿圈故障)或齿圈(行星轮故障)的旋转频率,Zg为某个齿轮的齿数;则太阳轮、行星轮和齿圈分布式故障的特征频率分别为:
fs'=fm/Zs(6);
fp'=fm/Zp(7);
fr'=fm/Zr(8);
式中:fm为啮合频率;fs'、fp'、fr'太阳轮、行星轮和齿圈分布式故障的特征频率;Zs为太阳轮齿数;Zp为行星轮的齿数;Zr为齿圈齿数;
以与主轴相连接的一级行星轮系的行星架为参考转速,对行星齿轮箱中各级的各个齿轮的局部故障和分布故障的特征阶比进行计算。
4.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:根据步骤(3)所述的行星齿轮箱故障特征量的提取;行星轮系局部故障,可分为太阳轮局部故障、行星轮局部故障和内齿圈局部故障;
对于太阳轮局部故障振动信号,在包络谱中,峰值出现在太阳轮的局部故障特征频率fs、太阳轮的绝对旋转频率fs (r)、以及它们的组合fs±fs (r)等位置处;若考虑太阳轮局部故障特征频率的倍频和太阳轮绝对旋转频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在太阳轮局部故障特征频率及其倍频nfs、太阳轮的绝对旋转频率及其倍频mfs (r)、及其组合nfs±mfs (r)等位置处;
对于行星轮局部故障,在包络谱中,峰值出现在行星轮局部故障特征频率fp、行星架的旋转频率fc、以及它们的fp±fc组合等位置处,若考虑行星轮局部故障特征频率的倍频和行星架旋转频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在行星轮局部故障特征频率及其倍频nfp、行星架的旋转频率及其倍频mfc、以及它们的组合nfp±mfc等位置处;
对于齿圈局部故障,在包络谱中,峰值出现在齿圈局部故障特征频率fr位置处,若考虑齿圈局部故障特征频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在齿圈局部故障特征频率及其倍频nfr位置处;
行星轮系发生分布故障,在包络谱中,峰值出现在齿轮分布式故障特征频率fg、行星轮通过频率Nfc、及其组合fg±Nfc位置处;若考虑齿轮分布式故障特征频率的倍频和行星轮通过频率的倍频作为调制频率的情形,则在包络谱中,峰值将出现在齿轮分布式故障特征频率及其倍频nfg(n为正整数)、行星轮通过频率及其倍频mNfc(m为正整数)、以及它们的组合nfg±mNfc位置处;
定轴齿轮发生包括齿根部有较大裂纹、局部齿面磨损、轮齿折断、局部齿形误差局部故障时,其振动信号波形是以齿轮旋转频率为周期的冲击脉冲,在频率域表现为包含旋转频率的各次谐波mfr(m=1,2,···)、各阶啮合频率nfm(n=1,2,···)以及以故障齿轮的旋转频率为间隔的边频nfm±mfr(n,m=1,2···);
定轴齿轮发生分布故障时,其频域特征表现为啮合频率及其谐波分量nfm(n=1,2,···)在频谱图上的位置保持不变,但其幅值大小发生改变,而且高次谐波幅值相对增大较多;分析时,要分析3个以上谐波幅值的变化,从而从频谱上检测出这种特征;
根据上面所述,在阶比谱和阶比包络谱中选取能表征故障类型的特征频率处的幅值为特征量,表示为行星轮系局部故障 A 1 = [ A nf s , A mf s ( r ) , A nf s ± mf s ( r ) , A nf p , A mf c , A nf p ± mf c , A nf r ] ; 行星轮系分布式故障 A 2 = [ A nf g , A mNf c ] ; 平行级啮合齿轮故障 A 3 = [ A mf r , A nf m , A nf m ± mf r ] ; m,n=1,2,3,其中平行级啮合齿轮故障归为一个特征向量,此向量从阶比图中提取,行星轮系的特征向量从包络谱中提取。
5.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:根据步骤(4)所述的故障诊断参数;即故障模式及严重程度的识别;J-散度(J-divergence)是一种谱距离,作为一种状态识别的指标,能很好地反映两个信号的相似程度,克服时域分析中的相位问题,很显然,同一信号的J-散度为零;
J ( s , τ ) = 1 2 N Σ i = 0 N - 1 [ S ( i ) τ ( i ) + τ ( i ) S ( i ) ] - 1 - - - ( 9 ) ;
式(9)中,S为样本信号的幅值谱;τ为待检信号的幅值谱;J(s,τ)为两者之间的J-散度,N为信号幅值谱中幅值的个数,i为信号幅值谱中幅值的序列;
J-散度作为一种状态识别的指标,它能很好地反映两个信号的相似程度,克服时域分析中的相位问题,很显然,同一信号的J-散度为零;
KL-散度(Kullback-Lciblerdivergence,KLD),用来度量分布P和Q之间的差异性,典型情况下,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;
对离散分布,P和Q的KL-散度定义为:
D k l ( P | | Q ) = Σ i P ( i ) lg P ( i ) Q ( i ) - - - ( 10 ) ;
式(10)中,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;Dkl为数据P和Q的KL-散度,i为数据的序列号;P(i)为序列号i所对应的P数据点;Q(i)为序列号i所对应的Q数据点;
也有人称其为KL距离,但是它并不是严格的距离概念,其不满足三角不等式;所以,把它变为对称形式:
Dkls(P||Q)=[Dkl(P||Q)+Dkl(Q||P)]/2(11);
式(11)中,P表示数据点真实分布,Q表示数据的理论分布、模型分布或P的近似分布;D'kl数据点P、Q之间的对称散度值;
当齿轮箱正常时,本文所选择的故障特征阶比的幅值为0或很小,而当发生故障时,其所对应的故障阶比的幅值发生很大的变化,这使得采用采用散度指标的算法得以实现,对上述所选特征量进行J-散度和KL-散度计算,计算方法:首先收集风电行星齿轮箱的正常状态的样本,表示正常情况下的标准样本,分别表示一级行星轮系、二级行星轮系和平行级齿轮故障样本,对每一样本进行阶比重采样、EMD重构信号以及进行阶比谱、Hilbert包络谱分析,寻找对应的故障特征阶比,形成对应故障特征阶比的幅值集合,最后计算待检样本和标准样本的故障特征集合之间的J-散度和KL-散度,定位故障位置以及故障模式,对行星齿轮箱的故障实现完全的诊断;
在时域信号进行FFT变换时的分辨率的原理推出,角域信号进行FFT的分辨率为2π/θ,其中θ为角域信号的长度;信号为matlab仿真信号,信号中忽略齿轮箱中轴承以及各齿轮相互之间的影响;
由此,仿真行星齿轮箱各级处于正常工况时,其振动信号模型为:
x1(t)=A·[1-cos(2π·3fc1·t)]cos(2π·fm1·t+θ1)
(12);
+B·[1-cos(2π·3fc2·t)]·cos(2π·fm2·t+θ2)+C·cos(2π·fm3·t+θ3)
行星齿轮箱一级行星轮系太阳轮发生局部故障时,其振动信号模型为:
行星齿轮箱一级行星轮系太阳轮发生分布故障时,其振动信号模型为:
式(12)~(14)中::x1(t)、x2(t)、x3(t)为行星齿轮箱处于正常一级行星轮系太阳轮发生局部故障、分布故障时的振动信号序列;t为时间序列;θ1、θ2、θ3、φ、为初始相位;fm1、fm2、fm3为各级啮合频率;fc1、fc2为一、二级行星架的旋转频率;为一级太阳轮的绝对旋转频率;fs1、fs1'为一级太阳轮发生局部故障和分布故障时的特征频率;A、B、C为无量纲常数,行星齿轮箱各个状态时值会有所不同;各个振动信号采用频率为8192HZ。
6.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:根据步骤(5)所述的实验验证;是在步骤1、步骤4的基础上,以行星齿轮箱的正常、行星齿轮箱一级行星轮系太阳轮局部故障、行星齿轮箱一级行星轮系太阳轮分布故障时的振动信号为例来进行分析。
7.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:所述的阶比重采样方法有效避免振动信号非平稳特点引发的频率冗杂、难以分析的特点,通过EMD分解消噪后重构后的阶比谱和包络谱清晰,容易找到对应的特征阶次,避免噪声对有效信号的干扰。
8.根据权利要求5所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法其特征是:所述的行星齿轮箱一级行星轮系太阳轮发生局部故障采集的振动信号经阶比重采样和EMD分解重构后的阶比谱;阶比谱中一级行星轮系啮合频率附近出现边频带,与正常时对比发生很明显的变化,同时对比正常与故障的阶比包络谱,故障时的调制成分发生很大的变化,与正常标准样本进行散度指标的计算,正常时故障特征阶比集合所对应的故障阶比的幅值很小,甚至为0,而发生故障时其故障特征阶比的幅值变化很大,通过与正常样本进行散度指标的计算,确定散度指标越大,发生故障的可能性越大,确定故障发生的粗略位置以及故障的严重程度。
9.根据权利要求5所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:所述的行星齿轮箱一级行星轮系太阳轮,发生分布故障时所采集的振动信号经阶比重采样和EMD分解重构后进行频谱分析和包络谱分析;阶比谱中一级行星轮系啮合频率附近出现边频带,与正常时对比发生很明显的变化,同时对比分布故障与局部故障的阶比谱不同,观察阶比包络谱,故障信号中存在复杂的调频信息,计算故障特征向量对应的幅值的散度指标,可以实现对分布故障进行有效的诊断,同时也可以有效的识别出局部故障和分布故障;
以风电行星齿轮箱为研究对象,对行星齿轮箱进行拆分为一级行星轮系、二级行星轮系、平行级齿轮三部分,同时行星轮系故障分为局部故障、分布式故障,暂不考虑轴承的影响,同时相比较基于散度指标的轴承故障诊断,其所采用的对样本之间进行散度指标的计算,来观察之间的相似性,而行星齿轮箱结构复杂,缺少特定的故障样本数据,对故障位置、故障模式的确定具有更大的难度,决定采用故障样本只与正常标准样进行散度指标计算,输入散度指标的故障特征向量分别为一级行星轮系A11、A12,二级行轮系A21、A22,平行级齿轮A3,通过计算观察某一故障特征向量较其他发生很大变化来确定故障发生的位置以及故障模式,对行星齿轮箱正常标准样本N与行星齿轮箱一级行星轮系太阳轮局部故障、分布故障样本的故障特征向量进行计算,对比分析,找出行星齿轮箱故障发生的哪一级以及哪一种故障模式,即此样本故障发生位置和模式,再通过对A11进行进一步计算,完全推出故障发生的位置,实现对故障的完全诊断;
由其他位置发生故障时各个故障样本与标准样本的J-散度、KL-散度值的变化可以发现,故障样本与正常标准样本的散度指标计算值越大,就可以确定此样本为在此故障特征向量所包含的故障模式集合中的一类,这样可以基本锁定齿轮箱发生故障的位置和类型,通过对此故障特征集合进行细致区别计算,完全诊断出故障位置和模式。
10.根据权利要求1所述的基于散度指标的变工况风电行星齿轮箱故障诊断方法,其特征是:所述的散度指标,经过计算可能够诊断出行星齿轮箱的故障模式及其位置,这使得风电机组运维人员可以完全避免查看复杂的频谱图,通过观察一些指标的变化可以诊断出齿轮箱的运行状态,同时散度指标的大小也可以衡量故障的严重程度。
CN201510831633.1A 2015-11-24 2015-11-24 基于散度指标的变工况风电行星齿轮箱故障诊断方法 Active CN105510023B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201510831633.1A CN105510023B (zh) 2015-11-24 2015-11-24 基于散度指标的变工况风电行星齿轮箱故障诊断方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201510831633.1A CN105510023B (zh) 2015-11-24 2015-11-24 基于散度指标的变工况风电行星齿轮箱故障诊断方法

Publications (2)

Publication Number Publication Date
CN105510023A true CN105510023A (zh) 2016-04-20
CN105510023B CN105510023B (zh) 2019-10-11

Family

ID=55718166

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201510831633.1A Active CN105510023B (zh) 2015-11-24 2015-11-24 基于散度指标的变工况风电行星齿轮箱故障诊断方法

Country Status (1)

Country Link
CN (1) CN105510023B (zh)

Cited By (22)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107192524A (zh) * 2017-04-05 2017-09-22 天津大学 一种考虑强谐波干扰的风电结构工作模态参数识别方法
CN107727394A (zh) * 2017-12-01 2018-02-23 华能国际电力股份有限公司 一种获取风力机齿轮箱齿数的方法
CN108151869A (zh) * 2017-11-27 2018-06-12 广州航新航空科技股份有限公司 一种机械振动特征指标提取方法、系统及装置
CN108181105A (zh) * 2017-11-28 2018-06-19 杭州安脉盛智能技术有限公司 基于逻辑回归和j散度的滚动轴承故障预诊方法及系统
CN108226297A (zh) * 2018-01-15 2018-06-29 国网江苏省电力公司检修分公司特高压交直流运检中心 一种基于光纤光栅的真空管波纹管表面裂纹检测方法
CN108613806A (zh) * 2018-07-23 2018-10-02 潍柴动力股份有限公司 一种齿轮检测方法及装置
CN108982135A (zh) * 2017-06-02 2018-12-11 上海金艺检测技术有限公司 热轧立辊轧机运行状态的在线监测方法
CN110068461A (zh) * 2019-05-29 2019-07-30 温州职业技术学院 齿轮齿条一体化传动系统多工况测试设备
CN110132578A (zh) * 2019-06-01 2019-08-16 吉林大学 一种齿轮系统复合故障特征提取方法及故障试验装置
CN110398362A (zh) * 2018-04-19 2019-11-01 中国科学院沈阳自动化研究所 一种机器人rv减速器故障诊断和定位方法
CN110686830A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机活塞环状态检测方法
CN110686768A (zh) * 2019-10-17 2020-01-14 昆明理工大学 一种改进的旋转机械非平稳振动信号计算阶比分析方法
CN110686892A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机弹性传动齿轮状态检测方法
CN110686890A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机气阀状态检测方法
CN110686879A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机气缸套状态检测方法
CN111222247A (zh) * 2020-01-13 2020-06-02 北京化工大学 一种旋转机械早期故障预警方法
CN113723732A (zh) * 2020-05-25 2021-11-30 中国石油化工股份有限公司 用于离心泵的状态判定方法及系统
CN113761675A (zh) * 2021-07-23 2021-12-07 东北大学 基于边频分布规律的行星轮轮齿裂纹故障特征判定方法
CN114235388A (zh) * 2021-12-15 2022-03-25 盛瑞传动股份有限公司 变速箱故障检测方法、装置、设备及存储介质
CN115165340A (zh) * 2022-06-23 2022-10-11 上海交通大学 复杂传动系统的振动信号频域特征提取方法和系统
CN115468645A (zh) * 2022-08-26 2022-12-13 国网湖北省电力有限公司黄冈供电公司 一种基于振动信号分段时频图谱优选的有载分接开关故障智能诊断方法
WO2024138791A1 (zh) * 2022-12-27 2024-07-04 江苏徐工国重实验室科技有限公司 变速箱故障的检测方法和装置

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE3311618A1 (de) * 1983-03-30 1984-10-04 Zahnrad- und Getriebefabrik Siegfried F. Tandler, 2800 Bremen Vorrichtung zur guetebestimmung von zahnradgetrieben
CN103308152A (zh) * 2013-06-06 2013-09-18 沈阳大学 基于瞬时频率估计的旋转机械振动信号角域重采样方法
CN103411774A (zh) * 2013-07-17 2013-11-27 华北电力大学 波动工况下的风电机组在线预警方法
CN104392082A (zh) * 2014-07-10 2015-03-04 中山火炬职业技术学院 一种基于振动监测的风力发电机组齿轮箱早期故障诊断方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE3311618A1 (de) * 1983-03-30 1984-10-04 Zahnrad- und Getriebefabrik Siegfried F. Tandler, 2800 Bremen Vorrichtung zur guetebestimmung von zahnradgetrieben
CN103308152A (zh) * 2013-06-06 2013-09-18 沈阳大学 基于瞬时频率估计的旋转机械振动信号角域重采样方法
CN103411774A (zh) * 2013-07-17 2013-11-27 华北电力大学 波动工况下的风电机组在线预警方法
CN104392082A (zh) * 2014-07-10 2015-03-04 中山火炬职业技术学院 一种基于振动监测的风力发电机组齿轮箱早期故障诊断方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
吴冠宇 等: "基于阶比分析和散度指标的风电机组行星齿轮箱分级故障诊断分析方法", 《风力发电》 *

Cited By (29)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107192524A (zh) * 2017-04-05 2017-09-22 天津大学 一种考虑强谐波干扰的风电结构工作模态参数识别方法
CN108982135A (zh) * 2017-06-02 2018-12-11 上海金艺检测技术有限公司 热轧立辊轧机运行状态的在线监测方法
CN108151869B (zh) * 2017-11-27 2020-07-03 广州航新航空科技股份有限公司 一种机械振动特征指标提取方法、系统及装置
CN108151869A (zh) * 2017-11-27 2018-06-12 广州航新航空科技股份有限公司 一种机械振动特征指标提取方法、系统及装置
CN108181105A (zh) * 2017-11-28 2018-06-19 杭州安脉盛智能技术有限公司 基于逻辑回归和j散度的滚动轴承故障预诊方法及系统
CN107727394A (zh) * 2017-12-01 2018-02-23 华能国际电力股份有限公司 一种获取风力机齿轮箱齿数的方法
CN108226297A (zh) * 2018-01-15 2018-06-29 国网江苏省电力公司检修分公司特高压交直流运检中心 一种基于光纤光栅的真空管波纹管表面裂纹检测方法
CN110398362A (zh) * 2018-04-19 2019-11-01 中国科学院沈阳自动化研究所 一种机器人rv减速器故障诊断和定位方法
CN108613806A (zh) * 2018-07-23 2018-10-02 潍柴动力股份有限公司 一种齿轮检测方法及装置
CN110068461A (zh) * 2019-05-29 2019-07-30 温州职业技术学院 齿轮齿条一体化传动系统多工况测试设备
CN110068461B (zh) * 2019-05-29 2020-05-05 温州职业技术学院 齿轮齿条一体化传动系统多工况测试设备
CN110132578A (zh) * 2019-06-01 2019-08-16 吉林大学 一种齿轮系统复合故障特征提取方法及故障试验装置
CN110132578B (zh) * 2019-06-01 2019-12-24 吉林大学 一种齿轮系统复合故障特征提取方法及故障试验装置
CN110686768A (zh) * 2019-10-17 2020-01-14 昆明理工大学 一种改进的旋转机械非平稳振动信号计算阶比分析方法
CN110686890A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机气阀状态检测方法
CN110686892B (zh) * 2019-10-23 2021-07-13 中船动力有限公司 在线柴油机弹性传动齿轮状态检测方法
CN110686892A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机弹性传动齿轮状态检测方法
CN110686879A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机气缸套状态检测方法
CN110686830A (zh) * 2019-10-23 2020-01-14 中船动力有限公司 在线柴油机活塞环状态检测方法
CN110686879B (zh) * 2019-10-23 2021-07-13 中船动力有限公司 在线柴油机气缸套状态检测方法
CN110686830B (zh) * 2019-10-23 2021-07-13 中船动力有限公司 在线柴油机活塞环状态检测方法
CN111222247A (zh) * 2020-01-13 2020-06-02 北京化工大学 一种旋转机械早期故障预警方法
CN113723732A (zh) * 2020-05-25 2021-11-30 中国石油化工股份有限公司 用于离心泵的状态判定方法及系统
CN113761675A (zh) * 2021-07-23 2021-12-07 东北大学 基于边频分布规律的行星轮轮齿裂纹故障特征判定方法
CN113761675B (zh) * 2021-07-23 2023-09-22 东北大学 基于边频分布规律的行星轮轮齿裂纹故障特征判定方法
CN114235388A (zh) * 2021-12-15 2022-03-25 盛瑞传动股份有限公司 变速箱故障检测方法、装置、设备及存储介质
CN115165340A (zh) * 2022-06-23 2022-10-11 上海交通大学 复杂传动系统的振动信号频域特征提取方法和系统
CN115468645A (zh) * 2022-08-26 2022-12-13 国网湖北省电力有限公司黄冈供电公司 一种基于振动信号分段时频图谱优选的有载分接开关故障智能诊断方法
WO2024138791A1 (zh) * 2022-12-27 2024-07-04 江苏徐工国重实验室科技有限公司 变速箱故障的检测方法和装置

Also Published As

Publication number Publication date
CN105510023B (zh) 2019-10-11

Similar Documents

Publication Publication Date Title
CN105510023A (zh) 基于散度指标的变工况风电行星齿轮箱故障诊断方法
Teng et al. Compound faults diagnosis and analysis for a wind turbine gearbox via a novel vibration model and empirical wavelet transform
Han et al. Fault feature extraction of low speed roller bearing based on Teager energy operator and CEEMD
Antoniadou et al. A time–frequency analysis approach for condition monitoring of a wind turbine gearbox under varying load conditions
Ma et al. Deep residual learning with demodulated time-frequency features for fault diagnosis of planetary gearbox under nonstationary running conditions
Li et al. A fault diagnosis method for planetary gearboxes under non-stationary working conditions using improved Vold-Kalman filter and multi-scale sample entropy
Wang et al. Unknown fault feature extraction of rolling bearings under variable speed conditions based on statistical complexity measures
Pan et al. Nonlinear sparse mode decomposition and its application in planetary gearbox fault diagnosis
He et al. A novel order tracking method for wind turbine planetary gearbox vibration analysis based on discrete spectrum correction technique
Li et al. A new noise-controlled second-order enhanced stochastic resonance method with its application in wind turbine drivetrain fault diagnosis
Feng et al. A phase angle based diagnostic scheme to planetary gear faults diagnostics under non-stationary operational conditions
CN105938468A (zh) 一种滚动轴承的故障诊断方法
Urbanek et al. Time–frequency approach to extraction of selected second-order cyclostationary vibration components for varying operational conditions
Wang et al. Optimal demodulation subband selection for sun gear crack fault diagnosis in planetary gearbox
Liu et al. Rotating machinery fault diagnosis under time-varying speeds: A review
Liu et al. Flexible generalized demodulation for intelligent bearing fault diagnosis under nonstationary conditions
Liu et al. Fault diagnosis of wind turbines under nonstationary conditions based on a novel tacho-less generalized demodulation
Wei et al. Variational nonlinear component decomposition for fault diagnosis of planetary gearboxes under variable speed conditions
CN101587017A (zh) 一种基于局部均值分解循环频率谱的齿轮故障诊断方法
CN105928702B (zh) 基于形态分量分析的变工况齿轮箱轴承故障诊断方法
CN105865776A (zh) 一种基于eemd和广义s变换的风电齿轮箱故障诊断方法
CN105445022A (zh) 一种基于双树复小波变换-熵特征融合的行星齿轮故障诊断方法
CN105548595A (zh) 一种提取风电齿轮箱各级轴转速检测方法
Lin et al. A review and strategy for the diagnosis of speed-varying machinery
CN103234702A (zh) 一种叶片不平衡故障诊断方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CP01 Change in the name or title of a patent holder

Address after: Room 502000, building 0108, Hohhot autonomous region, Inner Mongolia

Patentee after: ELECTRIC POWER RESEARCH INSTITUTE OF STATE GRID EASTERN INNER MONGOLIA POWER Co.,Ltd.

Patentee after: State Grid Inner Mongolia East Power Integrated Energy Service Co.,Ltd.

Patentee after: STATE GRID CORPORATION OF CHINA

Address before: Room 502000, building 0108, Hohhot autonomous region, Inner Mongolia

Patentee before: ELECTRIC POWER RESEARCH INSTITUTE OF STATE GRID EASTERN INNER MONGOLIA POWER Co.,Ltd.

Patentee before: STATE GRID INNER MONGOLIA EASTERN ELECTRIC POWER ENERGY SAVING SERVICES Co.,Ltd.

Patentee before: State Grid Corporation of China

CP01 Change in the name or title of a patent holder