CN108051189A - 一种旋转机械故障特征提取方法及装置 - Google Patents
一种旋转机械故障特征提取方法及装置 Download PDFInfo
- Publication number
- CN108051189A CN108051189A CN201711156210.XA CN201711156210A CN108051189A CN 108051189 A CN108051189 A CN 108051189A CN 201711156210 A CN201711156210 A CN 201711156210A CN 108051189 A CN108051189 A CN 108051189A
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- signal
- complex
- scale
- 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
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
- G01M13/00—Testing of machine parts
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Complex Calculations (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
本发明涉及一种旋转机械故障特征提取方法及装置,包括:提出一种二元特征尺度分解方法——复局部特征尺度分解(Complex Local characteristic‑scale decomposition,CLCD);采集同一截面水平及垂直方向的振动信号,将其组成一个复数信号;通过CLCD将该信号自适应分解成多个复内禀尺度分量,用全矢谱技术融合各分量的实部和虚部特征。本发明通过提出的CLCD方法直接对双通道的信号进行处理,确保了各个通道信号有相同的分解尺度,便于信息融合。并且运用了双线性变换,所以该方法自适应性强,故障特征提取全面、速度快、效率高,具有良好的技术效果。
Description
技术领域
本发明涉及一种旋转机械故障特征提取方法及装置,属于机械振动技术领域。
背景技术
近年来,旋转机械故障诊断成为研究的热点,如何从故障振动信号中准确提取振动信号特征则是研究旋转机械故障的关键。旋转机械在不同时刻、不同方向呈现的振动强度不同,基于单通道信息的故障诊断方法会割裂转子不同方向振动信号的特征信息,不能全面反映故障特征。考虑到稳态时转子同一截面互相垂直方向上各谐波的振荡轨迹为椭圆,人们提出了基于同源信息技术的故障特征提取方法,如全频谱、全息谱和全矢谱。同源信息技术能把两个正交通道的振动信号很好的融合,更全面、准确地反映故障特征。但同源信息技术以傅立叶变换为基础而提取信息,适用于平稳信号的分析。为了分析非线性、非平稳的机械故障振动信号,与同源信息技术相结合的各种方法相继被提出。
如经验模态分解(EMD)和局部均值分解(LMD)方法,这两种方法均可将任一复杂信号自适应地分解为若干瞬时频率具有物理意义的相互独立的单分量信号,适用于分析非线性、非平稳的信号,但是它们处理双通道的二元振动信号时却采用一元信号处理方法,效率低,并会出现分解结果数目不统一等问题,而且EMD还存在过包络及欠包络问题、LMD也存在着信号突变问题。
复局部均值分解(CLMD)和全矢谱相结合的CLMD全矢包络技术虽然解决了分解尺度不统一的问题,但是在平滑局部均值函数和包络估计函数时采用了移动平均的方法,并且通过解调才能得到纯调频函数,因此仍存在着计算量大问题。
二元经验模态分解(BEMD)虽然可以实现二元信号的自适应分解,并有效提取旋转机械的故障特征,但作为EMD的扩展,BEMD仍存在模态混叠、端点效应和计算量大、运行时间长等缺点。
2012年,程军圣等人提出了一种新的自适应信号分解方法——局部特征尺度分解(Local characteristic-scale decomposition,LCD),由于采用线性变换的形式对相邻两个极值点之间的信号进行处理,和EMD和LMD相比,LCD在时频局部化、运算速度和包络拟合精度等方面表现出良好的优越性。
本发明在LCD的基础上,提出了一种二元特征尺度分解方法——CLCD,并将其用于旋转机械故障特征的提取中。将CLCD方法与全矢谱技术相结合,直接对双通道的信号进行处理,确保了各个通道信号有相同的分解尺度,便于信息融合。并且由于CLCD方法运用了双线性变换,所以该方法具有较强的自适应性,故障特征提取全面、速度快、效率高,具有良好的技术效果。
发明内容
本发明的目的是提供一种旋转机械故障特征提取方法及装置,以解决目前旋转机械故障识别过程中存在的运行时间长、效率低、采用一维信号处理存在的分解尺度不统一等问题。
为解决上述技术问题,本发明提供了一种旋转机械故障特征提取方法,包括以下步骤:
步骤1,提出一种二元特征尺度分解方法——复局部特征尺度分解(CLCD);
步骤2,采集旋转机械同一截面水平方向的振动信号x(t)及垂直方向的振动信号y(t),将其组成一个复数信号z(t)=x(t)+jy(t),其中,ij=i2=j2=-1。
步骤3,通过CLCD将该复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和,按照能量从高到低的顺序排列,提取前四个复内禀尺度分量;
步骤4,分别对各复内禀尺度分量进行复数傅立叶变换,复数傅立叶变换后得到的实部傅立叶谱表示水平方向振动信号的频谱,复数傅立叶变换后得到的虚部傅立叶谱表示垂直方向振动信号的频谱;
步骤5,用全矢谱技术融合各复内禀尺度分量的实部和虚部特征,得到相应的全矢谱。
步骤6,根据得到的全矢谱谱线特征提取旋转机械故障特征,得到的全矢谱谱线特征即是旋转机械故障特征。
进一步的,步骤1提出了一种二元特征尺度分解方法——复局部特征尺度分解(CLCD),具体实现步骤为:
步骤1-1,将复数信号z(t)分别投影到0方向和π/2方向,
z0(t)=Re(e-j0·z(t))
zπ/2(t)=Re(e-jπ/2·z(t));
步骤1-2,利用局部特征尺度分解分别计算z0(t)和zπ/2(t)的基线信号;
步骤1-3,将步骤1-2中得到的z0(t)和zπ/2(t)的基线信号分别从z0(t)和zπ/2(t)中分离出来,得到剩余信号I0(i,k)(t)和I(π/2)(i,k)(t),
I0(i,k)(t)=z0(t)-p0(i,k)(t)
I(π/2)(i,k)(t)=zπ/2(t)-p(π/2)(i,k)(t)
其中,p0(i,k)(t)为z0(t)的基线信号,p(π/2)(i,k)(t)为zπ/2(t)的基线信号,i为ISC分量的个数,ISC为内禀尺度分量,k为剩余信号满足ISC条件的最大迭代次数。
步骤1-4,根据步骤1-3中得到的实轴方向的剩余信号I0(i,k)(t)和虚轴方向的剩余信号I(π/2)(i,k)(t),计算复数信号的剩余信号I(i,k)(t),
I(i,k)(t)=e-j0·I0(i,k)(t)+e-jπ/2·I(π/2)(i,k)(t)
步骤1-5,根据步骤1-4得到的复数信号的剩余信号I(i,k)(t),可得到复ISC分量ci(t),
ci(t)=I(i,k)(t)
步骤1-6,将所有复内禀尺度分量和最终残留信号相加,得到复数信号的复局部特征尺度分解的完整表达式为:
其中,n为最大迭代次数,r(t)为最终残留信号。
进一步的,步骤1-3中满足的ISC条件为,
其中,Xk为每一个ISC分量在整个数据段内极值点,τk为各个极值点对应的时刻(k=1,2,...,M,其中,M为极值点的个数),设定常数a∈(0,1),典型地,a=0.5,Lk为均值点,设定变动量△,当|Lk+1|≤△时迭代结束;
进一步的,步骤1-6中通过循环处理得到最终残留信号,判断循环停止的条件是残留信号呈现单调变化特性或者是一常数。
进一步的,步骤4中的复内禀尺度分量按能量从高到低的顺序排列,对前四个复内禀尺度分量分别进行复数傅立叶变换,结果为:
其中,cn为步骤3中分解得到的复内禀尺度分量。
进一步的,步骤5中得到的全矢谱Rai为:
本发明还提供了一种旋转机械故障特征提取装置,包括检测单元、分解单元和转换单元;
所述检测单元用于采集旋转机械同一截面水平方向的振动信号x(t)及垂直方向的振动信号y(t),并将其组成一个复数信号z(t)=x(t)+jy(t),其中,ij=i2=j2=-1;
所述分解单元用于将所述复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和;
所述转换单元用于将各复内禀尺度分量从时域变换到频域,得到各复内禀尺度分量水平方向频谱、垂直方向频谱,以及融合后的全矢谱,从而提取旋转机械故障特征,得到的全矢谱谱线特征即是旋转机械故障特征。
进一步的,所述分解单元运用复局部特征尺度分解将所述复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和,具体步骤包括:
(1)将复数信号z(t)分别投影到0方向和π/2方向,
z0(t)=Re(e-j0·z(t))
zπ/2(t)=Re(e-jπ/2·z(t));
(2)利用局部特征尺度分解分别计算z0(t)和zπ/2(t)的基线信号;
(3)将步骤(2)中得到的z0(t)和zπ/2(t)的基线信号分别从z0(t)和zπ/2(t)中分离出来,得到剩余信号I0(i,k)(t)和I(π/2)(i,k)(t),
I0(i,k)(t)=z0(t)-p0(i,k)(t)
I(π/2)(i,k)(t)=zπ/2(t)-p(π/2)(i,k)(t)
其中,p0(i,k)(t)为z0(t)的基线信号,p(π/2)(i,k)(t)为zπ/2(t)的基线信号,i为ISC分量的个数,ISC为内禀尺度分量,k为剩余信号满足ISC条件的最大迭代次数。
(4)根据步骤(3)中得到的实轴方向的剩余信号I0(i,k)(t)和虚轴方向的剩余信号I(π/2)(i,k)(t),计算复数信号的剩余信号I(i,k)(t),
I(i,k)(t)=e-j0·I0(i,k)(t)+e-jπ/2·I(π/2)(i,k)(t)
(5)根据步骤(4)得到的复数信号的剩余信号I(i,k)(t),可得到复ISC分量ci(t),
ci(t)=I(i,k)(t)
(6)将所有复内禀尺度分量和最终残留信号相加,得到复数信号的复局部特征尺度分解的完整表达式为:
其中,n为最大迭代次数,r(t)为最终残留信号。
进一步的,所述步骤(3)中满足的ISC条件为,
设定变动量△,当|Lk+1|≤△时迭代结束;其中,Xk为每一个ISC分量在整个数据段内极值点,τk为各个极值点对应的时刻(k=1,2,...,M,其中,M为极值点的个数),设定常数a∈(0,1),典型地,a=0.5,Lk为均值点,
进一步的,所述步骤(6)中通过循环处理得到最终残留信号,判断循环停止的条件是残留信号呈现单调变化特性或者是一常数。
进一步的,所述转换单元将各复内禀尺度分量按能量从高到低的顺序排列,对前四个复内禀尺度分量分别进行复数傅立叶变换,结果为:
其中,cn为所述分解单元将所述复数信号自适应分解得到的复内禀尺度分量。
复数傅立叶变换后得到的实部傅立叶谱表示水平方向振动信号的频谱,复数傅立叶变换后得到的虚部傅立叶谱表示垂直方向振动信号的频谱;
进一步的,用全矢谱技术融合各复内禀尺度分量的实部和虚部特征,得到相应的全矢谱Rai为:
本发明的有益效果是:本发明首先提出一种二元特征尺度分解方法——复局部特征尺度分解(CLCD),并将它应用于旋转机械故障特征的提取中;采集同一截面上水平及垂直方向的振动信号,并将两个通道的振动信号组成一个复数信号;通过复局部特征尺度分解将复数信号自适应地分解成多个复内禀尺度分量与一个最终残留信号之和,按照能量从高到低的顺序排列,提取前四个复内禀尺度分量;分别对各复内禀尺度分量进行复数傅立叶变换,复数傅立叶变换后得到的实部傅立叶谱表示水平方向振动信号的频谱,复数傅立叶变换后得到的虚部傅立叶谱表示垂直方向振动信号的频谱;用全矢谱技术融合各复内禀尺度分量的实部和虚部特征,得到相应的全矢谱;根据得到的全矢谱谱线特征提取旋转机械故障特征。本发明通过提出的CLCD方法及提取装置直接对双通道的信号进行处理,确保了各个通道信号有相同的分解尺度,便于信息融合。并且运用了双线性变换,所以该方法自适应性强,故障特征提取全面、速度快、效率高,具有良好的技术效果。
附图说明
图1-a为水平方向振动信号的时域图;
图1-b为垂直方向振动信号的时域图;
图1-c为复数信号z(t)平面图;
图1-d为复数信号z(t)三维图;
图1-e为复数信号z(t)实部傅立叶谱;
图1-f为复数信号z(t)虚部傅立叶谱;
图2为各复内禀尺度分量的实部和虚部时域波形图;
图3-a为各复内禀尺度分量水平方向振动信号的频谱图;
图3-b为各复内禀尺度分量垂直方向振动信号的频谱图;
图3-c为各复内禀尺度分量的全矢频谱图。
具体实施方式
下面结合附图以及具体的实施例对本发明的具体实施方式作进一步详细说明。
本发明针对现有技术的不足,提出了一种二元特征尺度分解方法——复局部特征尺度分解(CLCD),将其与全矢谱技术相结合,直接对双通道的信号进行处理,确保了各个通道信号有相同的分解尺度,便于信息融合。由于CLCD运用了双线性变换,所以该方法自适应性强,故障特征提取全面、速度快、效率高,具有良好的技术效果。该方法的具体实施过程如下。
第一步,提出一种二元特征尺度分解方法——复局部特征尺度分解(CLCD),具体实现步骤为:
1)将复数信号z(t)投影到0方向,
z0(t)=Re(e-j0·z(t)) (1)
2)利用局部特征尺度分解计算z0(t)的基线信号;
3)将z0(t)的基线信号从z0(t)中分离出来,得到剩余信号I0(i,k)(t),
I0(i,k)(t)=z0(t)-p0(i,k)(t) (2)
其中,p0(i,k)(t)为z0(t)的基线信号,i为ISC分量的个数,ISC为内禀尺度分量,k为剩余信号满足ISC条件的最大迭代次数。
4)剩余信号满足的ISC条件为,
其中,Xk为每一个ISC分量在整个数据段内极值点,τk为各个极值点对应的时刻(k=1,2,...,M,其中,M为极值点的个数),设定常数a∈(0,1),典型地,a=0.5,Lk为均值点,设定变动量△,当|Lk+1|≤△时迭代结束;
5)采用同样方法,将复数信号z(t)投影到π/2方向,
zπ/2(t)=Re(e-jπ/2·z(t)) (3)
6)利用局部特征尺度分解计算zπ/2(t)的基线信号p(π/2)(i,k)(t);
7)将zπ/2(t)的基线信号p(π/2)(i,k)(t)从zπ/2(t)中分离出来,得到剩余信号I(π/2)(i,k)(t),
I(π/2)(i,k)(t)=zπ/2(t)-p(π/2)(i,k)(t) (4)
8)复数信号的剩余信号I(i,k)(t)可通过下式求取,
I(i,k)(t)=e-j0·I0(i,k)(t)+e-jπ/2·I(π/2)(i,k)(t) (5)
9)复ISC分量ci(t)为,
ci(t)=I(i,k)(t) (6)
10)残留信号ri(t)可表示为,
ri(t)=z(t)-ci(t) (7)
11)判断残留信号ri(t)是否呈现单调变化特性或者是一常数,通过循环处理直至残留信号呈现单调变化特性或者是一常数,最终残留信号用r(t)表示。将所有复内禀尺度分量和最终残留信号相加,得到复数信号的复局部特征尺度分解的完整表达式为:
12)步骤2)和步骤6)采用局部特征尺度分解计算基线信号的分解过程如下(以z0(t)的局部特征尺度分解过程为例说明)。
A.确定信号z0(t)的所有极值点Ni及对应的时刻τi(i=1,2,3,…,n),并构造时间序列。由于该序列不包含端点的数值,因此采用两端各增加一个极值点的方法进行延拓:
B.根据式(9)和(10)计算基线信号控制点的纵坐标Pi,其中,i=1,2,…,n;
C.求原始信号的第i个(i=1,2,…,n-1)区间的基线信号段Pi(t)。采用式(11)所示的分段线性的方法:
其中,t∈(τi,τi+1]。
D.依次连接所求的各个基线信号段,得到基线信号p0(t)。
第二步,采集同一截面水平方向及垂直方向上的振动信号,本实施例中利用电涡流传感器检测转子外表面同一截面水平方向的振动信号x(t)、垂直方向的振动信号y(t)(采样频率为2048Hz,采样时长为0.5s,转子转速1703转/s),水平方向振动信号的时域图如图1-a所示,垂直方向振动信号的时域图如图1-b所示,从波形图上可知,水平方向振动信号和垂直方向振动信号的时域波形图相似,但是振动能量却不相同。
第三步,将两个方向的振动信号直接组成一个复数信号z(t),z(t)=x(t)+jy(t),其中,ij=i2=j2=-1,如图1-c所示为该复数信号的平面图,图1-d是该复数信号的三维图。采用复数傅立叶变换对复数信号z(t)进行变换,得到实部傅立叶谱(如图1-e所示)和虚部傅立叶谱(如图1-f所示);实部傅立叶谱表示水平方向振动信号的频谱,虚部傅立叶谱表示垂直方向振动信号的频谱;由图1-e和图1-f可知,旋转机械故障出现了转子不平衡故障,而且信号在水平方向的频谱图和垂直方向的频谱图无论是幅值大小,还是结构分布均不相同,这说明利用单方向的频谱图提取旋转机械的故障特征不够全面,易产生误判,因此需要采用综合分析的方法。
第四步,通过本发明提出的复局部特征尺度分解将该复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和,按照能量从高到低的顺序排列,提取前四个复内禀尺度分量,并将它们投影到水平和垂直方向,得到提取的各复内禀尺度分量的实部和虚部时域波形图,如图2所示。
第五步,分别对提取的各复内禀尺度分量进行复数傅立叶变换,复数傅立叶变换后得到的实部傅立叶谱表示水平方向振动信号的频谱,即提取的各复内禀尺度分量的实部特征,如图3-a所示;复数傅立叶变换后得到的虚部傅立叶谱表示垂直方向振动信号的频谱,即提取的各复内禀尺度分量的虚部特征,如图3-b所示。由图3-a和图3-b可以看出,同源振动信号x(t)和y(t)的振动能量不同,利用CLCD自适应分解得到的水平方向各振动分量和垂直方向各振动分量的频谱结构、能量分布均呈现不同特征。第一个复内禀尺度分量水平方向分量cx1在5倍频处振动能量最强,其垂直方向分量cy1在3倍频和5倍频处呈现了一定强度的振动能量;第二个复内禀尺度分量水平方向分量cx2和垂直方向分量cy2,无论是频谱结构还是能量分布均呈现出很大的不同;第三个复内禀尺度分量水平方向分量cx3在特征频率处振动能量最大,虽然垂直方向分量cy3和水平方向分量cx3的频谱结构相似,但其振动能量明显较小;第四个复内禀尺度分量水平方向分量cx4在二分频处的振动能量约为垂直方向分量cy4的一倍。
第六步,为了提高旋转机械故障特征提取的准确性,将全矢谱技术引入到CLCD方法中,融合两个方向的信号后再做频谱分析,实现CLCD分解的同源信息融合,图3-c为融合后得到的四个复内禀尺度分量的全矢频谱图,得到的全矢谱谱线特征即是旋转机械故障特征。对比图3-a至图3-c可知,全矢谱能有效地融合提取的各复内禀尺度分量的实部特征和虚部特征,频谱图更能真实地提取旋转机械的故障特征。
本发明提出了一种二元特征尺度分解方法——复局部特征尺度分解(CLCD),并将该方法与全矢谱技术相结合,直接对双通道的信号进行处理,确保了各个通道信号有相同的分解尺度,便于信息融合。并且由于CLCD方法运用了双线性变换,所以该方法具有较强的自适应性,故障特征提取全面、速度快、效率高,具有良好的技术效果。
以上所述,仅为本发明的具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,可轻易想到变化或替换,都应涵盖在本发明的保护范围之内。因此,本发明的保护范围应所述以权利要求的保护范围为准。
Claims (10)
1.一种旋转机械故障特征提取方法,其特征在于,包括以下步骤:
步骤1,提出一种二元特征尺度分解方法——复局部特征尺度分解(CLCD);
步骤2,采集旋转机械同一截面水平方向的振动信号x(t)及垂直方向的振动信号y(t),将其组成一个复数信号z(t)=x(t)+jy(t),其中,ij=i2=j2=-1;
步骤3,通过CLCD将该复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和,按照能量从高到低的顺序排列,提取前四个复内禀尺度分量;
步骤4,分别对各复内禀尺度分量进行复数傅立叶变换,复数傅立叶变换后得到的实部傅立叶谱表示水平方向振动信号的频谱,复数傅立叶变换后得到的虚部傅立叶谱表示垂直方向振动信号的频谱;
步骤5,用全矢谱技术融合各复内禀尺度分量的实部和虚部特征,得到相应的全矢谱;
步骤6,根据得到的全矢谱谱线特征提取旋转机械故障特征。
2.根据权利要求1所述的一种旋转机械故障特征提取方法,其特征在于,步骤1提出了一种二元特征尺度分解方法——复局部特征尺度分解(CLCD),具体实现步骤为:
步骤1-1,将复数信号z(t)分别投影到0方向和π/2方向,
z0(t)=Re(e-j0·z(t))
zπ/2(t)=Re(e-jπ/2·z(t));
步骤1-2,利用局部特征尺度分解分别计算z0(t)和zπ/2(t)的基线信号;
步骤1-3,将步骤1-2中得到的z0(t)和zπ/2(t)的基线信号分别从z0(t)和zπ/2(t)中分离出来,得到剩余信号I0(i,k)(t)和I(π/2)(i,k)(t),
I0(i,k)(t)=z0(t)-p0(i,k)(t)
I(π/2)(i,k)(t)=zπ/2(t)-p(π/2)(i,k)(t)
其中,p0(i,k)(t)为z0(t)的基线信号,p(π/2)(i,k)(t)为zπ/2(t)的基线信号,i为ISC分量的个数,ISC为内禀尺度分量,k为剩余信号满足ISC条件的最大迭代次数;
步骤1-4,根据步骤1-3中得到的实轴方向的剩余信号I0(i,k)(t)和虚轴方向的剩余信号I(π/2)(i,k)(t),计算复数信号的剩余信号I(i,k)(t),
I(i,k)(t)=e-j0·I0(i,k)(t)+e-jπ/2·I(π/2)(i,k)(t)
步骤1-5,根据步骤1-4得到的复数信号的剩余信号I(i,k)(t),可得到复ISC分量ci(t),
ci(t)=I(i,k)(t)
步骤1-6,将所有复内禀尺度分量和最终残留信号相加,得到复数信号的复局部特征尺度分解的完整表达式为:
<mrow>
<mi>z</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>n</mi>
</munderover>
<msub>
<mi>c</mi>
<mi>i</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mi>r</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
其中,n为最大迭代次数,r(t)为最终残留信号。
3.根据权利要求2所述的一种旋转机械故障特征提取方法,其特征在于,所述步骤1-3中满足的ISC条件为,
<mrow>
<msub>
<mi>L</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<mi>a</mi>
<mo>&lsqb;</mo>
<msub>
<mi>X</mi>
<mi>k</mi>
</msub>
<mo>+</mo>
<mfrac>
<mrow>
<msub>
<mi>&tau;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>&tau;</mi>
<mi>k</mi>
</msub>
</mrow>
<mrow>
<msub>
<mi>&tau;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>2</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>&tau;</mi>
<mi>k</mi>
</msub>
</mrow>
</mfrac>
<mrow>
<mo>(</mo>
<msub>
<mi>X</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>2</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>X</mi>
<mi>k</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<mo>+</mo>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>-</mo>
<mi>a</mi>
<mo>)</mo>
</mrow>
<msub>
<mi>X</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
其中,Xk为每一个ISC分量在整个数据段内极值点,τk为各个极值点对应的时刻(k=1,2,...,M,其中,M为极值点的个数),设定常数a∈(0,1),典型地,a=0.5,Lk为均值点,设定变动量△,当|Lk+1|≤△时迭代结束;
所述步骤1-6中通过循环处理得到最终残留信号,判断循环停止的条件是残留信号呈现单调变化特性或者是一常数。
4.根据权利要求1所述的一种旋转机械故障特征提取方法,其特征在于,所述步骤4中的复内禀尺度分量按能量从高到低的顺序排列,对前四个复内禀尺度分量分别进行复数傅立叶变换,结果为:
<mrow>
<msub>
<mi>Z</mi>
<mi>k</mi>
</msub>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>n</mi>
<mo>=</mo>
<mn>0</mn>
</mrow>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</munderover>
<msub>
<mi>c</mi>
<mi>n</mi>
</msub>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mi>j</mi>
<mn>2</mn>
<mi>&pi;</mi>
<mi>n</mi>
<mi>k</mi>
<mo>/</mo>
<mi>N</mi>
</mrow>
</msup>
<mo>,</mo>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>=</mo>
<mn>0</mn>
<mo>,</mo>
<mn>1</mn>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
其中,cn为步骤3中分解得到的复内禀尺度分量。
5.根据权利要求4所述的一种旋转机械故障特征提取方法,其特征在于,所述步骤5中得到的全矢谱Rai为:
<mrow>
<msub>
<mi>R</mi>
<mrow>
<mi>a</mi>
<mi>i</mi>
</mrow>
</msub>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mrow>
<mn>2</mn>
<mi>N</mi>
</mrow>
</mfrac>
<mo>&lsqb;</mo>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mi>i</mi>
</msub>
<mo>|</mo>
<mo>+</mo>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mi>i</mi>
</mrow>
</msub>
<mo>|</mo>
<mo>&rsqb;</mo>
<mo>,</mo>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
<mo>,</mo>
<mn>2</mn>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<mi>N</mi>
<mo>/</mo>
<mn>2</mn>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mo>.</mo>
</mrow>
6.一种旋转机械故障特征提取装置,其特征在于,包括检测单元、分解单元和转换单元;
所述检测单元用于采集旋转机械同一截面水平方向的振动信号x(t)及垂直方向的振动信号y(t),并将其组成一个复数信号z(t)=x(t)+jy(t),其中,ij=i2=j2=-1;
所述分解单元用于将所述复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和;
所述转换单元用于将各复内禀尺度分量从时域变换到频域,得到各复内禀尺度分量水平方向频谱、垂直方向频谱,以及融合后的全矢谱,从而提取旋转机械故障特征。
7.根据权利要求6所述的一种旋转机械故障特征提取装置,其特征在于,所述分解单元运用复局部特征尺度分解将所述复数信号自适应分解成多个复内禀尺度分量与一个最终残留信号之和;具体步骤包括:
(1)将复数信号z(t)分别投影到0方向和π/2方向,
z0(t)=Re(e-j0·z(t))
zπ/2(t)=Re(e-jπ/2·z(t));
(2)利用局部特征尺度分解分别计算z0(t)和zπ/2(t)的基线信号;
(3)将步骤(2)中得到的z0(t)和zπ/2(t)的基线信号分别从z0(t)和zπ/2(t)中分离出来,得到剩余信号I0(i,k)(t)和I(π/2)(i,k)(t),
I0(i,k)(t)=z0(t)-p0(i,k)(t)
I(π/2)(i,k)(t)=zπ/2(t)-p(π/2)(i,k)(t)
其中,p0(i,k)(t)为z0(t)的基线信号,p(π/2)(i,k)(t)为zπ/2(t)的基线信号,i为ISC分量的个数,ISC为内禀尺度分量,k为剩余信号满足ISC条件的最大迭代次数;
(4)根据步骤(3)中得到的实轴方向的剩余信号I0(i,k)(t)和虚轴方向的剩余信号I(π/2)(i,k)(t),计算复数信号的剩余信号I(i,k)(t),
I(i,k)(t)=e-j0·I0(i,k)(t)+e-jπ/2·I(π/2)(i,k)(t)
(5)根据步骤(4)得到的复数信号的剩余信号I(i,k)(t),可得到复ISC分量ci(t),ci(t)=I(i,k)(t)
(6)将所有复内禀尺度分量和最终残留信号相加,得到复数信号的复局部特征尺度分解的完整表达式为:
<mrow>
<mi>z</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>n</mi>
</munderover>
<msub>
<mi>c</mi>
<mi>i</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mi>r</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
其中,n为最大迭代次数,r(t)为最终残留信号。
8.根据权利要求7所述的一种旋转机械故障特征提取装置,其特征在于,所述步骤(3)中满足的ISC条件为,
<mrow>
<msub>
<mi>L</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<mi>a</mi>
<mo>&lsqb;</mo>
<msub>
<mi>X</mi>
<mi>k</mi>
</msub>
<mo>+</mo>
<mfrac>
<mrow>
<msub>
<mi>&tau;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>&tau;</mi>
<mi>k</mi>
</msub>
</mrow>
<mrow>
<msub>
<mi>&tau;</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>2</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>&tau;</mi>
<mi>k</mi>
</msub>
</mrow>
</mfrac>
<mrow>
<mo>(</mo>
<msub>
<mi>X</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>2</mn>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>X</mi>
<mi>k</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<mo>+</mo>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>-</mo>
<mi>a</mi>
<mo>)</mo>
</mrow>
<msub>
<mi>X</mi>
<mrow>
<mi>k</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
其中,Xk为每一个ISC分量在整个数据段内极值点,τk为各个极值点对应的时刻(k=1,2,...,M,其中,M为极值点的个数),设定常数a∈(0,1),典型地,a=0.5,Lk为均值点,设定变动量△,当|Lk+1|≤△时迭代结束;
所述步骤(6)中通过循环处理得到最终残留信号,判断循环停止的条件是残留信号呈现单调变化特性或者是一常数。
9.根据权利要求6-8中任一项所述的一种旋转机械故障特征提取装置,其特征在于,所述转换单元将各复内禀尺度分量按能量从高到低的顺序排列,对前四个复内禀尺度分量分别进行复数傅立叶变换,结果为:
<mrow>
<msub>
<mi>Z</mi>
<mi>k</mi>
</msub>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>n</mi>
<mo>=</mo>
<mn>0</mn>
</mrow>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</munderover>
<msub>
<mi>c</mi>
<mi>n</mi>
</msub>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mi>j</mi>
<mn>2</mn>
<mi>&pi;</mi>
<mi>n</mi>
<mi>k</mi>
<mo>/</mo>
<mi>N</mi>
</mrow>
</msup>
<mo>,</mo>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>=</mo>
<mn>0</mn>
<mo>,</mo>
<mn>1</mn>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
其中,cn为所述分解单元将所述复数信号自适应分解得到的复内禀尺度分量;
复数傅立叶变换后得到的实部傅立叶谱表示水平方向振动信号的频谱,复数傅立叶变换后得到的虚部傅立叶谱表示垂直方向振动信号的频谱。
10.根据权利要求9所述的一种旋转机械故障特征提取装置,其特征在于,用全矢谱技术融合各复内禀尺度分量的实部和虚部特征,得到相应的全矢谱Rai为:
<mrow>
<msub>
<mi>R</mi>
<mrow>
<mi>a</mi>
<mi>i</mi>
</mrow>
</msub>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mrow>
<mn>2</mn>
<mi>N</mi>
</mrow>
</mfrac>
<mo>&lsqb;</mo>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mi>i</mi>
</msub>
<mo>|</mo>
<mo>+</mo>
<mo>|</mo>
<msub>
<mi>Z</mi>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mi>i</mi>
</mrow>
</msub>
<mo>|</mo>
<mo>&rsqb;</mo>
<mo>,</mo>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
<mo>,</mo>
<mn>2</mn>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<mi>N</mi>
<mo>/</mo>
<mn>2</mn>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mo>.</mo>
</mrow>
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201711156210.XA CN108051189B (zh) | 2017-11-20 | 2017-11-20 | 一种旋转机械故障特征提取方法及装置 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201711156210.XA CN108051189B (zh) | 2017-11-20 | 2017-11-20 | 一种旋转机械故障特征提取方法及装置 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN108051189A true CN108051189A (zh) | 2018-05-18 |
CN108051189B CN108051189B (zh) | 2019-12-31 |
Family
ID=62118911
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201711156210.XA Active CN108051189B (zh) | 2017-11-20 | 2017-11-20 | 一种旋转机械故障特征提取方法及装置 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN108051189B (zh) |
Cited By (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109659957A (zh) * | 2019-01-10 | 2019-04-19 | 东北电力大学 | 基于apit-memd电力系统低频振荡模式辨识方法 |
CN111289275A (zh) * | 2019-12-31 | 2020-06-16 | 北京化工大学 | 旋转机械稳定性辨识方法、装置、计算机设备及存储介质 |
CN111307438A (zh) * | 2020-03-11 | 2020-06-19 | 上海纵行实业有限公司 | 一种基于信息熵的旋转机械振动故障诊断方法及其系统 |
CN111982489A (zh) * | 2020-08-27 | 2020-11-24 | 江苏师范大学 | 选择集成的改进局部特征分解的微弱故障特征提取方法 |
CN112345247A (zh) * | 2019-08-09 | 2021-02-09 | 郑州工程技术学院 | 一种滚动轴承的故障诊断方法及装置 |
CN112345248A (zh) * | 2019-08-09 | 2021-02-09 | 郑州工程技术学院 | 一种滚动轴承的故障诊断方法及装置 |
CN113358282A (zh) * | 2021-06-04 | 2021-09-07 | 重庆大学 | 一种动平衡机上复合材料尾轴低速正交融合动平衡法 |
CN113640566A (zh) * | 2021-07-08 | 2021-11-12 | 国网江苏省电力有限公司电力科学研究院 | 一种foct漂移故障特征提取方法 |
CN115876511A (zh) * | 2022-12-29 | 2023-03-31 | 郑州机械研究所有限公司 | 一种旋转机械的故障确定方法、装置及电子设备 |
Citations (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH0961300A (ja) * | 1995-08-22 | 1997-03-07 | Honda Motor Co Ltd | ギアの噛合せ検査方法 |
CN101929917A (zh) * | 2010-03-24 | 2010-12-29 | 陈先利 | 一种旋转机械的故障诊断方法 |
CN103033362A (zh) * | 2012-12-31 | 2013-04-10 | 湖南大学 | 一种基于改进多变量预测模型的齿轮故障诊断方法 |
CN104165759A (zh) * | 2014-06-17 | 2014-11-26 | 中州大学 | 一种基于复局部均值分解的转子碰摩故障特征提取方法 |
CN104165686A (zh) * | 2014-06-17 | 2014-11-26 | 中州大学 | 一种基于二元经验模态分解的转子轴心轨迹提纯方法 |
CN104614178A (zh) * | 2015-03-10 | 2015-05-13 | 中电投河南电力有限公司技术信息中心 | 一种基于全矢谱的故障征兆提取方法 |
CN104634571A (zh) * | 2015-02-06 | 2015-05-20 | 北京航空航天大学 | 一种基于lcd-mf的滚动轴承故障诊断方法 |
CN106772032A (zh) * | 2016-11-29 | 2017-05-31 | 西安理工大学 | 一种水轮发电机组的故障特征提取方法 |
CN106932191A (zh) * | 2017-04-14 | 2017-07-07 | 中州大学 | 一种转子三维空间轴心轨迹提纯方法及装置 |
CN106989907A (zh) * | 2017-04-14 | 2017-07-28 | 中州大学 | 一种转子故障特征提取方法及装置 |
-
2017
- 2017-11-20 CN CN201711156210.XA patent/CN108051189B/zh active Active
Patent Citations (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH0961300A (ja) * | 1995-08-22 | 1997-03-07 | Honda Motor Co Ltd | ギアの噛合せ検査方法 |
CN101929917A (zh) * | 2010-03-24 | 2010-12-29 | 陈先利 | 一种旋转机械的故障诊断方法 |
CN103033362A (zh) * | 2012-12-31 | 2013-04-10 | 湖南大学 | 一种基于改进多变量预测模型的齿轮故障诊断方法 |
CN104165759A (zh) * | 2014-06-17 | 2014-11-26 | 中州大学 | 一种基于复局部均值分解的转子碰摩故障特征提取方法 |
CN104165686A (zh) * | 2014-06-17 | 2014-11-26 | 中州大学 | 一种基于二元经验模态分解的转子轴心轨迹提纯方法 |
CN104634571A (zh) * | 2015-02-06 | 2015-05-20 | 北京航空航天大学 | 一种基于lcd-mf的滚动轴承故障诊断方法 |
CN104614178A (zh) * | 2015-03-10 | 2015-05-13 | 中电投河南电力有限公司技术信息中心 | 一种基于全矢谱的故障征兆提取方法 |
CN106772032A (zh) * | 2016-11-29 | 2017-05-31 | 西安理工大学 | 一种水轮发电机组的故障特征提取方法 |
CN106932191A (zh) * | 2017-04-14 | 2017-07-07 | 中州大学 | 一种转子三维空间轴心轨迹提纯方法及装置 |
CN106989907A (zh) * | 2017-04-14 | 2017-07-28 | 中州大学 | 一种转子故障特征提取方法及装置 |
Non-Patent Citations (1)
Title |
---|
黄传金 等: "全矢局部特征尺度分解及其在转子碰摩故障特征提取中的应用", 《机械强度》 * |
Cited By (16)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109659957B (zh) * | 2019-01-10 | 2022-05-17 | 东北电力大学 | 基于apit-memd电力系统低频振荡模式辨识方法 |
CN109659957A (zh) * | 2019-01-10 | 2019-04-19 | 东北电力大学 | 基于apit-memd电力系统低频振荡模式辨识方法 |
CN112345247B (zh) * | 2019-08-09 | 2022-11-25 | 郑州工程技术学院 | 一种滚动轴承的故障诊断方法及装置 |
CN112345248B (zh) * | 2019-08-09 | 2022-11-25 | 郑州工程技术学院 | 一种滚动轴承的故障诊断方法及装置 |
CN112345247A (zh) * | 2019-08-09 | 2021-02-09 | 郑州工程技术学院 | 一种滚动轴承的故障诊断方法及装置 |
CN112345248A (zh) * | 2019-08-09 | 2021-02-09 | 郑州工程技术学院 | 一种滚动轴承的故障诊断方法及装置 |
CN111289275A (zh) * | 2019-12-31 | 2020-06-16 | 北京化工大学 | 旋转机械稳定性辨识方法、装置、计算机设备及存储介质 |
CN111289275B (zh) * | 2019-12-31 | 2021-06-01 | 北京化工大学 | 旋转机械稳定性辨识方法、装置、计算机设备及存储介质 |
CN111307438A (zh) * | 2020-03-11 | 2020-06-19 | 上海纵行实业有限公司 | 一种基于信息熵的旋转机械振动故障诊断方法及其系统 |
CN111982489B (zh) * | 2020-08-27 | 2022-05-06 | 江苏师范大学 | 选择集成的改进局部特征分解的微弱故障特征提取方法 |
CN111982489A (zh) * | 2020-08-27 | 2020-11-24 | 江苏师范大学 | 选择集成的改进局部特征分解的微弱故障特征提取方法 |
CN113358282A (zh) * | 2021-06-04 | 2021-09-07 | 重庆大学 | 一种动平衡机上复合材料尾轴低速正交融合动平衡法 |
CN113358282B (zh) * | 2021-06-04 | 2023-04-21 | 重庆大学 | 一种动平衡机上复合材料尾轴低速正交融合动平衡方法 |
CN113640566A (zh) * | 2021-07-08 | 2021-11-12 | 国网江苏省电力有限公司电力科学研究院 | 一种foct漂移故障特征提取方法 |
CN113640566B (zh) * | 2021-07-08 | 2024-04-26 | 国网江苏省电力有限公司电力科学研究院 | 一种foct漂移故障特征提取方法 |
CN115876511A (zh) * | 2022-12-29 | 2023-03-31 | 郑州机械研究所有限公司 | 一种旋转机械的故障确定方法、装置及电子设备 |
Also Published As
Publication number | Publication date |
---|---|
CN108051189B (zh) | 2019-12-31 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN108051189A (zh) | 一种旋转机械故障特征提取方法及装置 | |
Chen et al. | High-accuracy fault feature extraction for rolling bearings under time-varying speed conditions using an iterative envelope-tracking filter | |
CN101806832B (zh) | 一种低频率信号的频率测量方法 | |
Wang et al. | Bearing fault diagnosis under time-varying rotational speed via the fault characteristic order (FCO) index based demodulation and the stepwise resampling in the fault phase angle (FPA) domain | |
Yang et al. | Fault feature extraction based on combination of envelope order tracking and cICA for rolling element bearings | |
CN105738696A (zh) | 全相位时移相位差频率估计方法及装置 | |
CN101403774B (zh) | 一种基于非同步采样的谐波分析方法 | |
CN106483374A (zh) | 一种基于Nuttall双窗全相位FFT的谐波间谐波检测方法 | |
Ye et al. | An improved empirical mode decomposition based on adaptive weighted rational quartic spline for rolling bearing fault diagnosis | |
CN102288821B (zh) | 三相电路相位差的测量方法、测量装置 | |
Dashtdar et al. | Fault location in the transmission network using artificial neural network | |
CN109459131A (zh) | 一种旋转机械多通道振动信号的时频特征提取方法及装置 | |
CN101509945A (zh) | 正负序电量实时检测的方法 | |
CN102221639A (zh) | 正负序电流实时检测的方法 | |
CN102809687B (zh) | 一种交流电频率的数字化测量方法 | |
CN113468760B (zh) | 基于字典学习的电机微弱故障检测方法及系统 | |
CN104034412A (zh) | 一种基于分数阶全息原理的旋转机械故障特征提取方法 | |
Kabla et al. | Bearing fault diagnosis using Hilbert-Huang transform (HHT) and support vector machine (SVM) | |
CN110907864B (zh) | 电机定子绕组的故障检测方法、装置、设备及存储介质 | |
CN104698325A (zh) | 一种电力系统负阻尼机理低频振荡和强迫振荡的判别方法 | |
Xu et al. | Match-extracting chirplet transform with application to bearing fault diagnosis | |
Yi et al. | Time-varying fault feature extraction of rolling bearing via time–frequency sparsity | |
CN103513133B (zh) | Mcr型svc装置动态响应时间的测量装置及方法 | |
Granjon et al. | Estimation of geometric properties of three-component signals for system monitoring | |
Lv et al. | High-order synchroextracting chirplet transform for accurate instantaneous frequency estimation and its application in fault diagnosis of rotary machinery |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |