技术领域Technical field
本发明涉及电力系统稳定与控制技术领域,具体涉及一种基于卷积神经网络(convolutional neural network,CNN)的电力系统振荡类型的快速辨识方法。The invention relates to the technical field of power system stability and control, and in particular to a method for rapid identification of power system oscillation types based on a convolutional neural network (CNN).
背景技术Background technique
自电力系统诞生以来,振荡研究就成为了关乎电力系统动态性能与稳定性的重要研究之一。经过多年的研究,传统电力系统中低频振荡(low frequency oscillation,LFO)和次同步振荡(sub-synchronous oscillation,SSO)的特性已被充分揭示。它们的共同特性是由较大物理惯量的旋转单元所主导,尤其是大型同步发电机组。然而,由于可再生能源的快速发展,如太阳能、风能、地热能、生物质能等,大量的新能源并网发电装备通过电力电子逆变装置接入电网,这给电力系统的规划以及运行控制带来了巨大的挑战。传统控制策略的电力电子逆变器侧重于发电和电能质量的控制,具有响应速度快、低惯量等特性,大量的采用传统电力电子技术,会导致系统呈现出惯量减弱的趋势,严重危及电力系统的稳定运行。例如,2009年10月,美国德克萨斯州的双馈风力发电系统和串联补偿电网之间发生了频率约为20Hz的次同步谐振事故,造成消弧电路损坏和大规模断路故障。考虑到电力系统振荡对电力系统稳定运行的严重威胁,快速辨识出电力系统振荡类型对电网后续动态稳定性的评估及调整具有重要意义。Since the birth of power systems, oscillation research has become one of the important studies related to the dynamic performance and stability of power systems. After years of research, the characteristics of low frequency oscillation (LFO) and sub-synchronous oscillation (SSO) in traditional power systems have been fully revealed. Their common characteristic is that they are dominated by rotating units with large physical inertia, especially large synchronous generator units. However, due to the rapid development of renewable energy, such as solar energy, wind energy, geothermal energy, biomass energy, etc., a large number of new energy grid-connected power generation equipment is connected to the grid through power electronic inverter devices, which poses challenges to the planning and operation control of the power system. brought huge challenges. Power electronic inverters with traditional control strategies focus on power generation and power quality control, and have characteristics such as fast response and low inertia. Extensive use of traditional power electronics technology will cause the system to show a trend of weakening inertia, seriously endangering the power system. stable operation. For example, in October 2009, a subsynchronous resonance accident with a frequency of approximately 20 Hz occurred between the doubly-fed wind power generation system and the series compensation grid in Texas, causing damage to the arc suppression circuit and large-scale circuit breakage. Considering that power system oscillation poses a serious threat to the stable operation of the power system, quickly identifying the type of power system oscillation is of great significance to the subsequent assessment and adjustment of the dynamic stability of the power grid.
目前,电力系统振荡类型的辨识主要是模态辨识,对应的方法分为基于模型的方法和基于实测信号的方法。由于无需系统精确的模型和参数,后者在电力系统中得到了广泛应用。常见算法有基于快速傅里叶变换(fast fourier transformation,FFT)的算法、小波分析、Prony、希尔伯特-黄变换(HHT)等以及近年来发展起来的借助旋转不变技术估计信号参数(ESPRIT:estimation ofsignal parameters via rotational invariancetechnique),矩阵束方法(MP:matrix pencil algorithm),随机子空间方法(SSI:stochastic subspace identification)等。经过多年的研究与改进,这些方法可以较为准确地获取系统的模态,在抗噪性上也有一定的效果,但在辨识的数学机理上,这些算法将被测信号视作平稳随机过程,并且往往需要采集数秒钟乃至数分钟的数据。随着可再生能源并网比例的提高以及大量的电力电子设备接入电网,所产生的振荡表现出不同于传统电力系统振荡的新特性。由于可再生能源的发电机单元和AC/DC电网之间的复杂动态交互,所产生的振荡一般具有从几赫兹到几千赫兹的宽频范围,同时还具有强时变性、非线性、非平稳性以及强噪声干扰等特性。在这种背景下,传统电力系统中振荡信号近似平稳的假设难以保证,亟需发展新的振荡类型辨识方法来适应电力系统的发展。At present, the identification of power system oscillation types is mainly mode identification, and the corresponding methods are divided into model-based methods and methods based on measured signals. The latter has been widely used in power systems since accurate models and parameters of the system are not required. Common algorithms include algorithms based on fast Fourier transformation (FFT), wavelet analysis, Prony, Hilbert-Huang transform (HHT), etc., as well as the estimation of signal parameters developed in recent years with the help of rotation invariant technology ( ESPRIT: estimation of signal parameters via rotational invariancetechnique), matrix pencil method (MP: matrix pencil algorithm), stochastic subspace method (SSI: stochastic subspace identification), etc. After years of research and improvement, these methods can obtain the mode of the system more accurately and have certain effects on noise immunity. However, in terms of the mathematical mechanism of identification, these algorithms regard the measured signal as a stationary random process, and It is often necessary to collect data from seconds to minutes. As the proportion of renewable energy connected to the grid increases and a large number of power electronic equipment are connected to the grid, the resulting oscillations show new characteristics that are different from the traditional power system oscillations. Due to the complex dynamic interaction between the generator unit of renewable energy and the AC/DC grid, the oscillation generated generally has a wide frequency range from a few hertz to several thousand hertz, and also has strong time variability, nonlinearity, and non-stationarity. and strong noise interference and other characteristics. In this context, the assumption that oscillation signals in traditional power systems are approximately stationary is difficult to guarantee, and there is an urgent need to develop new oscillation type identification methods to adapt to the development of power systems.
近年来,深度学习(deep learning,DL)发展迅猛,在众多领域已经得到了成功的应用。在电力系统振荡类型辨识方面,也有学者开展了探索性的研究。有研究提出使用深度学习算法辨识电力系统振荡模式是区域振荡还是局部振荡。还有研究提出指数型衰减正弦神经网络(exponentially damped sinusoids neural network,EDSNN)的辨识方法。但这些研究,均存在对被测振荡信号非平稳性特性考虑不足的缺陷。目前,同步相量测量单元(PMU)的广域测量系统(WAMS)在电力系统中获得了广泛应用,为电力系统振荡分析提供了数据来源,让深度学习在电力系统振荡类型辨识方面的应用成为可能。同时,由于事先通过大量的数据学习,深度学习的引入还使得电力系统振荡类型的辨识具有神经网络类算法的共同优点——辨识速度迅速。而这往往是系统状态强时变的高比例可再生能源并网电力系统中迫切需要的一个特点。因此,随着可再生能源渗透率的进一步提高,电力系统中的振荡模式愈发频繁且复杂,如何借鉴深度学习算法实现电力系统振荡类型的在线快速辨识,成为了一个亟待解决的技术问题。In recent years, deep learning (DL) has developed rapidly and has been successfully applied in many fields. In terms of identification of power system oscillation types, some scholars have also carried out exploratory research. Some studies have proposed using deep learning algorithms to identify whether the power system oscillation mode is regional oscillation or local oscillation. There are also studies that propose the identification method of exponentially damped sinusoids neural network (EDSNN). However, these studies all have the defect of insufficient consideration of the non-stationary characteristics of the measured oscillation signal. At present, the wide area measurement system (WAMS) of the synchronized phasor measurement unit (PMU) has been widely used in power systems, providing a data source for power system oscillation analysis, making the application of deep learning in the identification of power system oscillation types a possible. At the same time, due to a large amount of data learning in advance, the introduction of deep learning also enables the identification of power system oscillation types to have the common advantage of neural network algorithms - rapid identification. This is often an urgently needed feature in high-proportion renewable energy grid-connected power systems with strong time-varying system conditions. Therefore, as the penetration rate of renewable energy further increases, oscillation patterns in the power system become more frequent and complex. How to use deep learning algorithms to achieve online rapid identification of power system oscillation types has become an urgent technical problem that needs to be solved.
发明内容Contents of the invention
本发明的目的在于针对高比例可再生能源并网的电力电子化系统中易出现的复杂振荡现象,提供一种分析准确,设计合理的基于CNN的电力系统振荡类型的快速辨识方法。The purpose of the present invention is to provide a fast identification method of power system oscillation types based on CNN that is accurate in analysis and reasonably designed in view of the complex oscillation phenomena that are prone to occur in power electronic systems with a high proportion of renewable energy connected to the grid.
本发明可通过下述技术方案实现:The present invention can be realized through the following technical solutions:
一种基于CNN的电力系统振荡类型的快速辨识方法,由以下步骤构成:A fast identification method of power system oscillation types based on CNN, consisting of the following steps:
步骤1:生成振荡样本数据,即根据电力系统振荡信号数学模型为CNN模型生成训练样本数据与测试样本数据。Step 1: Generate oscillation sample data, that is, generate training sample data and test sample data for the CNN model based on the mathematical model of the power system oscillation signal.
所选振荡信号数学模型为指数型衰减正弦量(exponentially dampedsinusoids,EDSs),其公式为:The selected mathematical model of the oscillation signal is exponentially damped sinusoids (EDSs), and its formula is:
式(1)中,x(tj)为tj时刻的信号;Ai为幅值;σi为衰减因子;fi为频率;为相位;j=0,1,2,…,n-1;n为采样点数;i=1,2,…,m;m为模型实际阶数;/>为第i阶信号模态加入的时间节点;η(tj)代表噪声信号。In formula (1), x(tj ) is the signal at time tj ; Ai is the amplitude; σi is the attenuation factor; fiis the frequency; is the phase; j=0,1,2,…,n-1; n is the number of sampling points; i=1,2,…,m; m is the actual order of the model;/> is the time node added to the i-th order signal mode; eta (tj ) represents the noise signal.
步骤2:采用预处理算法处理电力系统振荡样本数据以便后续训练过程中的特征提取;Step 2: Use a preprocessing algorithm to process the power system oscillation sample data for feature extraction in the subsequent training process;
预处理算法采用平铺的方式将序列数据处理为n阶方阵,其公式为:The preprocessing algorithm uses a tiling method to process the sequence data into an n-order square matrix, and its formula is:
式(2)中,X为离散的时序序列,长度为N,由于目标矩阵为n阶方阵,因此N=n2;In formula (2), X is a discrete time series with length N. Since the target matrix is an n-order square matrix, N=n2 ;
此外,为了实现对特征矩阵的边缘滤波,对上述n阶方阵实施零填充操作,即在输入矩阵的边缘使用零值进行填充,使原输入矩阵转换为m(m>n)阶方阵。In addition, in order to achieve edge filtering of the feature matrix, a zero-filling operation is performed on the above-mentioned n-order square matrix, that is, zero values are used to fill the edges of the input matrix, so that the original input matrix is converted into an m (m>n)-order square matrix.
步骤3:确定电力系统振荡类型的分类准则,根据分类准则为样本数据添加类别标签。Step 3: Determine the classification criteria for power system oscillation types, and add category labels to the sample data according to the classification criteria.
分类准则为按存在性分类。一般而言,电力系统振荡信号按不同频率可大致分为四类,分别为:低频振荡(f≤2.5Hz)、次同步振荡(2.5Hz<f≤45Hz)、超同步振荡(50Hz<f<100Hz)、中高频振荡(f≥100Hz);对于某一高阶振荡信号,其可能包含任一类型振荡模态或同时包含多种类型振荡模态,按存在性分类即分别判断该四种振荡模态于某一振荡信号中是否存在,再按存在与不存在分为两类进行辨识,因此本发明需按照四种振荡类型训练四个不同的分类网络。The classification criterion is classification by existence. Generally speaking, power system oscillation signals can be roughly divided into four categories according to different frequencies, namely: low-frequency oscillation (f≤2.5Hz), sub-synchronous oscillation (2.5Hz<f≤45Hz), and super-synchronous oscillation (50Hz<f< 100Hz), medium and high frequency oscillation (f≥100Hz); for a certain high-order oscillation signal, it may contain any type of oscillation mode or multiple types of oscillation modes at the same time, and the four oscillations are judged separately according to the existence classification Whether the mode exists in a certain oscillation signal is divided into two categories for identification based on presence and absence. Therefore, the present invention requires training four different classification networks according to four oscillation types.
步骤4:根据振荡样本数据特征以及分类需求构建CNN模型。Step 4: Build a CNN model based on the oscillation sample data characteristics and classification requirements.
本发明所搭建的CNN模型,具体结构包含输入层、卷积层、池化层、全连接层、Dropout层以及输出层;The specific structure of the CNN model built by this invention includes an input layer, a convolution layer, a pooling layer, a fully connected layer, a Dropout layer and an output layer;
其中,输入层用以输入序列数据;卷积层为网络核心部分,负责将原始数据映射到隐层的特征空间;池化层用以缩小数据尺寸,减少全连接层的参数并加快训练速度;全连接层用以完成分布式特征表示到样本标记空间的映射;Dropout层通过令数据随机失活来预防过拟合现象的发生;输出层则输出最终的分类结果。Among them, the input layer is used to input sequence data; the convolution layer is the core part of the network, responsible for mapping the original data to the feature space of the hidden layer; the pooling layer is used to reduce the data size, reduce the parameters of the fully connected layer and speed up the training; The fully connected layer is used to complete the mapping of distributed feature representation to the sample label space; the Dropout layer prevents overfitting by randomly inactivating the data; the output layer outputs the final classification result.
步骤5:将训练样本输入CNN模型,训练模型的分类能力。Step 5: Input the training samples into the CNN model to train the classification ability of the model.
CNN模型训练过程为一个循环迭代过程,主要通过前向传播算法和反向传播算法更新网络各层级状态与参数,具体包括以下步骤:The CNN model training process is a cyclic iterative process, which mainly updates the status and parameters of each level of the network through the forward propagation algorithm and the back propagation algorithm. It specifically includes the following steps:
S51:初始化各网络参数,包括各神经元的权重和偏置等;S51: Initialize each network parameter, including the weight and bias of each neuron;
S52:输入训练样本数据;S52: Input training sample data;
S53:计算CNN网络内部各层级输出;S53: Calculate the output of each level within the CNN network;
S54:计算当前序列索引预测输出及其损失函数;S54: Calculate the current sequence index prediction output and its loss function;
S55:计算所有参数基于损失函数的偏导数,通过梯度下降法更新网络所有的参数;S55: Calculate the partial derivatives of all parameters based on the loss function, and update all parameters of the network through the gradient descent method;
S56:当损失函数不再下降或达到设定训练次数,完成训练。S56: When the loss function no longer decreases or reaches the set training times, the training is completed.
进一步的,在CNN模型训练过程中,所述CNN网络各层级输出的计算公式为:Further, during the CNN model training process, the calculation formula for the output of each level of the CNN network is:
1)卷积层输出:1) Convolutional layer output:
式(3)中,i=1,2,…,I;j=1,2,…,J;I=P-M+1;J=Q-N+1;P、Q为输入矩阵X的尺寸,km,n是M×N的卷积核的第m行,第n列的元素,xi+m-1,j+n-1是输入矩阵X中的第i+m-1行,第j+n-1列的元素。每个卷积核都会产生一个I×J的输出矩阵Y;In formula (3), i=1, 2,...,I; j=1, 2,...,J; I=P-M+1; J=Q-N+1; P and Q are the input matrix X Size, km,n is the m-th row and n-th column element of the M×N convolution kernel, xi+m-1,j+n-1 is the i+m-1-th row in the input matrix X , the element of column j+n-1. Each convolution kernel will produce an I×J output matrix Y;
2)池化层输出(平均池化):2) Pooling layer output (average pooling):
式(4)中,SR为池化区域面积,即池化区域元素个数;Yi,j为卷积层输出Y第i行第j列的元素;In formula (4), SR is the area of the pooling area, that is, the number of elements in the pooling area; Yi,j is the element in the i-th row and j-th column of the convolution layer output Y;
3)全连接层输出:3) Fully connected layer output:
式(5)中,hn为输出神经元的值;vm为输入神经元的值;wmn、bn分别为链接权值与偏置,m为输入神经元的数目;n为输出神经元数目。In formula (5), hn is the value of the output neuron; vm is the value of the input neuron; wmn and bn are the link weight and bias respectively, m is the number of input neurons; n is the output neuron. The number of yuan.
进一步的,在CNN模型训练过程中,所述预测输出值的计算公式及其损失函数的选取分别为:Further, during the CNN model training process, the calculation formula of the predicted output value and the selection of its loss function are:
式中,V和c表示预测输出的权值和偏置;tanh表示双曲正切函数;N为每次送入网络训练的批量大小;K为类型数目;yn,ij表示第i批次及第j类型样本的实际标签,则表示第i批次及第j类型样本的预测输出值。In the formula, V and c represent the weight and bias of the prediction output; tanh represents the hyperbolic tangent function; N is the batch size sent to the network training each time; K is the number of types; yn,ij represents the i-th batch and The actual label of the j-th type sample, Then represents the predicted output value of the i-th batch and j-th type sample.
步骤6:将测试样本输入CNN模型,测试各网络模型的辨识准确率,并通过调节网络参数不断提升辨识准确率以完成网络训练,当辨识准确率达到90%以上时,停止网络参数的调节,完成网络训练获得训练完成的CNN模型。Step 6: Input the test sample into the CNN model, test the recognition accuracy of each network model, and continuously improve the recognition accuracy by adjusting network parameters to complete network training. When the recognition accuracy reaches more than 90%, stop adjusting the network parameters. Complete network training to obtain the trained CNN model.
步骤7:通过滑动时窗获取电力系统振荡信号实测数据。对某段振荡信号采用滑动时窗取样,其过程如式(9)所示:Step 7: Obtain the measured data of the power system oscillation signal through the sliding time window. Sliding time window sampling is used for a certain section of oscillation signal, and the process is as shown in Equation (9):
式(9)中,X为离散的时序序列,长度为N;L为滑窗长度;fs为信号采样频率,由于滑动间隔为0.5s,则滑动步长等于fs/2;L1、L2、…、Ln+1即为通过滑动时窗获取的不同时段振荡信号实测数据。In Equation( 9), ,..., Ln+1 is the measured data of the oscillation signal in different periods obtained through the sliding time window.
依此法获取的待测信号按时序排列,通过对这些待测信号顺序逐一辨识,可以有效分析出待测振荡样本数据各模态的时变特征,以适应电力系统中振荡频繁且复杂的环境。The signals to be measured obtained by this method are arranged in time sequence. By identifying the order of these signals to be measured one by one, the time-varying characteristics of each mode of the oscillation sample data to be measured can be effectively analyzed to adapt to the frequent and complex oscillation environment in the power system. .
步骤8:处理待测振荡信号并输入CNN模型,根据网络输出分析振荡信号的所属类型。根据存在性分类准则,所有网络输出结果均为包含对应模态分段或不包含对应模态分段这两种情况。对于某段待测振荡信号,将其输入各网络模型并综合分析输出结果,即可获得其所包含的振荡类型。Step 8: Process the oscillation signal to be measured and input it into the CNN model, and analyze the type of the oscillation signal according to the network output. According to the existence classification criterion, all network output results contain the corresponding modal segment or do not contain the corresponding modal segment. For a certain oscillation signal to be measured, input it into each network model and comprehensively analyze the output results to obtain the oscillation type it contains.
本发明所述振荡类型辨识方法采用深度学习算法,选取CNN作为振荡样本训练模型,能够从较短时的振荡信号中迅速分析出其所属的类型,实现对电力系统振荡类型的在线快速辨识。The oscillation type identification method of the present invention adopts a deep learning algorithm and selects CNN as the oscillation sample training model, which can quickly analyze the type to which it belongs from short-term oscillation signals, and realize online rapid identification of power system oscillation types.
本发明与现有技术相比,具有如下的优点和有益效果:Compared with the prior art, the present invention has the following advantages and beneficial effects:
(1)本发明引入了人工智能领域的深度学习算法,利用深度神经网络的分类能力实现电力系统振荡类型的辨识,相较于传统的信号模型分析法,辨识更快更准确,更加适用于电力系统中振荡类型的实时辨识。(1) The present invention introduces a deep learning algorithm in the field of artificial intelligence, and uses the classification ability of deep neural networks to realize the identification of power system oscillation types. Compared with the traditional signal model analysis method, the identification is faster and more accurate, and is more suitable for electric power systems. Real-time identification of oscillation types in the system.
(2)本发明充分考虑到真实电网中常见的不稳定与时变振荡信号,采用深度学习算法并且选择定性辨识解决这一问题,使其更加适用于高比例可再生能源并网电力电子化系统中相对复杂的运行工况与振荡环境,相较于传统的信号模型分析法会将其视作稳定与时不变信号处理,得到的辨识结果既不是模态的真实数值也无法体现时变特征更有优势。(2) This invention fully takes into account the common unstable and time-varying oscillation signals in real power grids, uses a deep learning algorithm and selects qualitative identification to solve this problem, making it more suitable for high-proportion renewable energy grid-connected power electronic systems. In the relatively complex operating conditions and oscillation environment, compared with the traditional signal model analysis method, which will be treated as a stable and time-invariant signal, the identification results obtained are neither the real numerical values of the modes nor can they reflect the time-varying characteristics. have more advantages.
附图说明Description of drawings
此处所说明的附图用来提供对本发明实施例的进一步理解,构成本申请的一部分,并不构成对本发明实施例的限定。在附图中:The drawings described here are used to provide a further understanding of the embodiments of the present invention, constitute a part of this application, and do not constitute a limitation to the embodiments of the present invention. In the attached picture:
图1为本发明的工作流程图。Figure 1 is a work flow chart of the present invention.
图2为本发明的CNN模型结构图。Figure 2 is a structural diagram of the CNN model of the present invention.
图3为本发明的CNN训练流程图。Figure 3 is a CNN training flow chart of the present invention.
图4为含噪声情况下理想振荡信号图像。Figure 4 shows the ideal oscillation signal image with noise.
图5为电力系统实测振荡信号图像。Figure 5 shows the measured oscillation signal image of the power system.
具体实施方式Detailed ways
为使本发明的目的、技术方案和优点更加清楚明白,下面结合实施例和附图,对本发明作进一步的详细说明,本发明的示意性实施方式及其说明仅用于解释本发明,并不作为对本发明的限定。In order to make the purpose, technical solutions and advantages of the present invention more clear, the present invention will be further described in detail below in conjunction with the examples and drawings. The schematic embodiments of the present invention and their descriptions are only used to explain the present invention and do not as a limitation of the invention.
本发明中,S51-S56表示步骤5的详细子步骤。In the present invention, S51-S56 represent the detailed sub-steps of step 5.
如图1所示,一种基于CNN的电力系统振荡类型的快速辨识方法,包括如下步骤:步骤1:生成样本数据,即根据电力系统振荡信号数学模型为CNN模型生成训练样本数据与测试样本数据。As shown in Figure 1, a CNN-based rapid identification method of power system oscillation types includes the following steps: Step 1: Generate sample data, that is, generate training sample data and test sample data for the CNN model based on the mathematical model of power system oscillation signals .
所选振荡信号数学模型为EDSs,其公式为:The selected mathematical model of oscillation signal is EDSs, and its formula is:
式(1)中,x(tj)为tj时刻的信号;Ai为幅值;σi为衰减因子;fi为频率;为相位;i=1,2,…,m;j=0,1,2,…,n-1;n为采样点数;m为模型实际阶数;/>为第i阶信号模态加入的时间节点;η(tj)代表噪声信号。In formula (1), x(tj ) is the signal at time tj ; Ai is the amplitude; σi is the attenuation factor; fiis the frequency; is the phase; i=1,2,…,m; j=0,1,2,…,n-1; n is the number of sampling points; m is the actual order of the model;/> is the time node added to the i-th order signal mode; eta (tj ) represents the noise signal.
本实施例中,振荡信号采样时间取1s,采样频率取400Hz,其中采样时间即振荡信号模态辨识所需最短时间。生成信号过程中,振荡信号各模态参数在合理范围内随机取值,最高阶数取到4。此外,为验证本发明抗噪性能,噪声信号η(t)信噪比(signal-noiseratio,SNR)取10dB。In this embodiment, the oscillation signal sampling time is 1 s, and the sampling frequency is 400 Hz, where the sampling time is the shortest time required for mode identification of the oscillation signal. During the signal generation process, each modal parameter of the oscillation signal randomly takes values within a reasonable range, and the highest order is taken to 4. In addition, in order to verify the anti-noise performance of the present invention, the signal-noiseratio (SNR) of the noise signal eta (t) is 10 dB.
步骤2:采用预处理算法处理电力系统振荡样本数据以便后续训练过程中的特征提取。Step 2: Use a preprocessing algorithm to process the power system oscillation sample data for feature extraction in the subsequent training process.
预处理算法采用平铺的方式将序列数据处理为n阶方阵,其公式为:The preprocessing algorithm uses a tiling method to process the sequence data into an n-order square matrix, and its formula is:
式(2)中,X为离散的时序序列,长度为N=400,由于目标矩阵为n阶方阵,因此n=20;In formula (2), X is a discrete time series with a length of N=400. Since the target matrix is an n-order square matrix, n=20;
此外,为了实现对特征矩阵的边缘滤波,对上述20×20矩阵实施零填充操作,即在输入矩阵的边缘使用零值进行填充,使原输入矩阵转换为28阶方阵。In addition, in order to achieve edge filtering of the feature matrix, a zero-filling operation is performed on the above 20×20 matrix, that is, zero values are used to fill the edges of the input matrix, so that the original input matrix is converted into a 28th-order square matrix.
步骤3:确定电力系统振荡类型的分类准则,根据分类准则为样本数据添加类别标签。Step 3: Determine the classification criteria for power system oscillation types, and add category labels to the sample data according to the classification criteria.
分类准则为按存在性分类。一般而言,电力系统振荡信号按不同频率可大致分为四类,分别为:低频振荡(f≤2.5Hz)、次同步振荡(2.5Hz<f≤45Hz)、超同步振荡(50Hz<f<100Hz)、中高频振荡(f≥100Hz);对于某一高阶振荡信号,其可能包含任一类型振荡模态或同时包含多种类型振荡模态,按存在性分类即分别判断该四种振荡模态于某一振荡信号中是否存在,再按存在与不存在分为两类进行辨识,因此本发明需按照四种振荡类型训练四个不同的分类网络。The classification criterion is classification by existence. Generally speaking, power system oscillation signals can be roughly divided into four categories according to different frequencies, namely: low-frequency oscillation (f≤2.5Hz), sub-synchronous oscillation (2.5Hz<f≤45Hz), and super-synchronous oscillation (50Hz<f< 100Hz), medium and high frequency oscillation (f≥100Hz); for a certain high-order oscillation signal, it may contain any type of oscillation mode or multiple types of oscillation modes at the same time, and the four oscillations are judged separately according to the existence classification Whether the mode exists in a certain oscillation signal is divided into two categories for identification based on presence and absence. Therefore, the present invention requires training four different classification networks according to four oscillation types.
该分类方法可以实现不同阶数振荡信号的同时辨识,无需依赖传统辨识方法中难以实现的精确定阶过程,因此辨识结果更加可信。This classification method can realize the simultaneous identification of oscillation signals of different orders without relying on the precise order determination process that is difficult to achieve in traditional identification methods, so the identification results are more credible.
步骤4:根据振荡样本数据特征以及分类需求构建CNN模型。Step 4: Build a CNN model based on the oscillation sample data characteristics and classification requirements.
如图2所示,本发明所搭建的CNN模型具体结构为:输入层-卷积层1-池化层1-卷积层2-池化层2-全局平均池化层-Dropout层1-全连接层1-Dropout层2-全连接层2-Dropout层3-全连接层3-Softmax层-输出层。As shown in Figure 2, the specific structure of the CNN model built by the present invention is: input layer-convolution layer 1-pooling layer 1-convolution layer 2-pooling layer 2-global average pooling layer-Dropout layer 1- Fully connected layer 1 - Dropout layer 2 - Fully connected layer 2 - Dropout layer 3 - Fully connected layer 3 - Softmax layer - Output layer.
其中,输入层用以输入序列数据,根据输入数据格式,将输入层型号设定为28×28×1;卷积层为网络核心部分,负责将原始数据映射到隐层的特征空间,其中卷积层1与卷积层2的卷积核大小均设定为4×4;池化层用以缩小数据尺寸,减少全连接层的参数并加快训练速度,其中池化层1与池化层2均采用最大池化,其大小均设定为2×2;全局平均池化层用以替代全连接层,可以有效减少参数数量并防止过拟合现象发生;全连接层与Softmax层用以完成分布式特征表示到样本标记空间的映射以及结果归一化;Dropout层置于全连接层之后,同样为了防止过拟合现象发生,其随机失活概率为50%;输出层则输出最终的分类结果。Among them, the input layer is used to input sequence data. According to the input data format, the input layer model is set to 28×28×1; the convolution layer is the core part of the network and is responsible for mapping the original data to the feature space of the hidden layer. The convolution layer The convolution kernel sizes of the product layer 1 and the convolution layer 2 are both set to 4×4; the pooling layer is used to reduce the data size, reduce the parameters of the fully connected layer and speed up the training. Among them, the pooling layer 1 and the pooling layer 2 both use maximum pooling, and their sizes are set to 2×2; the global average pooling layer is used to replace the fully connected layer, which can effectively reduce the number of parameters and prevent over-fitting; the fully connected layer and Softmax layer are used to Complete the mapping of distributed feature representation to sample label space and normalization of results; the Dropout layer is placed after the fully connected layer. Also in order to prevent over-fitting, its random deactivation probability is 50%; the output layer outputs the final Classification results.
步骤5:将训练样本输入CNN模型,训练模型的分类能力。Step 5: Input the training samples into the CNN model to train the classification ability of the model.
CNN模型训练过程为一个循环迭代过程,如图3所示,主要通过前向传播算法和反向传播算法更新网络各层级状态与参数,具体包括以下步骤:The CNN model training process is a cyclic iterative process, as shown in Figure 3. It mainly updates the status and parameters of each level of the network through the forward propagation algorithm and the back propagation algorithm, which specifically includes the following steps:
S51:初始化各网络参数,包括各神经元的权重和偏置等;S51: Initialize each network parameter, including the weight and bias of each neuron;
S52:输入训练样本数据;S52: Input training sample data;
S53:计算CNN网络内部各层级输出;S53: Calculate the output of each level within the CNN network;
其中,CNN网络内部各层级输出的计算公式为:Among them, the calculation formula for the output of each level within the CNN network is:
1)卷积层输出:1) Convolutional layer output:
式(3)中,i=1,2,…,I;j=1,2,…,J;I=P-M+1;J=Q-N+1;P、Q为输入矩阵X的尺寸,km,n是M×N的卷积核的第m行,第n列的元素,xi+m-1,j+n-1是输入矩阵X中的第i+m-1行,第j+n-1列的元素。每个卷积核都会产生一个I×J的输出矩阵Y;In formula (3), i=1, 2,...,I; j=1, 2,...,J; I=P-M+1; J=Q-N+1; P and Q are the input matrix X Size, km,n is the m-th row and n-th column element of the M×N convolution kernel, xi+m-1,j+n-1 is the i+m-1-th row in the input matrix X , the element of column j+n-1. Each convolution kernel will produce an I×J output matrix Y;
2)池化层输出(平均池化):2) Pooling layer output (average pooling):
式(4)中,SR为池化区域面积,即池化区域元素个数;Yi,j为卷积层输出Y第i行第j列的元素;In formula (4), SR is the area of the pooling area, that is, the number of elements in the pooling area; Yi,j is the element in the i-th row and j-th column of the convolution layer output Y;
3)全连接层输出:3) Fully connected layer output:
式(5)中,hn为输出神经元的值;vm为输入神经元的值;wmn、bn分别为链接权值与偏置,m为输入神经元的数目;n为输出神经元数目。In formula (5), hn is the value of the output neuron; vm is the value of the input neuron; wmn and bn are the link weight and bias respectively, m is the number of input neurons; n is the output neuron. The number of yuan.
S54:计算当前序列索引预测输出及其损失函数;S54: Calculate the current sequence index prediction output and its loss function;
其中,预测输出值的计算公式及其损失函数的选取分别为:Among them, the calculation formula of the predicted output value and the selection of the loss function are:
式中,V和c表示预测输出的权值和偏置;tanh表示双曲正切函数;N为每次送入网络训练的批量大小;K为类型数目;yn,ij表示第i批次及第j类型样本的实际标签,则表示第i批次及第j类型样本的预测输出值。In the formula, V and c represent the weight and bias of the prediction output; tanh represents the hyperbolic tangent function; N is the batch size sent to the network training each time; K is the number of types; yn,ij represents the i-th batch and The actual label of the j-th type sample, Then represents the predicted output value of the i-th batch and j-th type sample.
S55:计算所有参数基于损失函数的偏导数,通过梯度下降法更新网络所有的参数;S55: Calculate the partial derivatives of all parameters based on the loss function, and update all parameters of the network through the gradient descent method;
S56:当损失函数不再下降或达到设定训练次数,完成训练。S56: When the loss function no longer decreases or reaches the set training times, the training is completed.
步骤6将测试样本输入CNN模型,测试各网络模型的辨识准确率,并通过调节网络参数使得辨识准确率能够达到95%左右,得到训练完成的神经网络模型。Step 6: Input the test sample into the CNN model to test the recognition accuracy of each network model. By adjusting the network parameters, the recognition accuracy can reach about 95%, and the trained neural network model is obtained.
根据上述步骤获得训练完成的神经网络模型,步骤7与步骤8将由两个具体的实施例进行阐释。The trained neural network model is obtained according to the above steps. Steps 7 and 8 will be explained by two specific embodiments.
实施例1:Example 1:
为验证该算法是否能够辨识出系统在振荡过程中叠加新的振荡模态,构造如下理想振荡测试信号:In order to verify whether the algorithm can identify the new oscillation mode superimposed during the oscillation process of the system, the following ideal oscillation test signal is constructed:
式(10)中,ε(t)表示阶跃函数,η(t)表示噪声信号。In Equation (10), ε(t) represents the step function, and η(t) represents the noise signal.
该振荡信号图像如图4所示,信号长度2s,信噪比SNR=10dB。1s前,信号包含两个模态,其中频率f1=0.84Hz,f2=1.21Hz,属于低频振荡。t=1s时,引入一个新的振荡模态,其频率f3=8.8Hz,属于次同步振荡,因此1s后的振荡信号中同时包含低频振荡与次同步振荡两种振荡类型。The oscillation signal image is shown in Figure 4, the signal length is 2s, and the signal-to-noise ratio SNR=10dB. 1s ago, the signal contained two modes, with frequencies f1 =0.84Hz and f2 =1.21Hz, which were low-frequency oscillations. When t=1s, a new oscillation mode is introduced with frequency f3 =8.8Hz, which belongs to subsynchronous oscillation. Therefore, the oscillation signal after 1s contains both low-frequency oscillation and subsynchronous oscillation.
步骤7通过滑动时窗获取待测振荡信号样本。滑窗长度为1s,滑动间隔为0.5s,采样频率为400Hz。为了辨识振荡信号中发生的模态改变,需在1s前后分别使用滑窗取样。本实施例中,选取四段信号作为辨识对象,分别为0~1s、0.5~1.5s与1~2s。其中0.5~1.5s包含了模态发生变化的时间节点,用以验证本发明针对复杂振荡模态辨识的能力。Step 7: Obtain the oscillation signal sample to be measured through the sliding time window. The length of the sliding window is 1s, the sliding interval is 0.5s, and the sampling frequency is 400Hz. In order to identify the modal changes that occur in the oscillation signal, sliding window sampling needs to be used before and after 1s. In this embodiment, four segments of signals are selected as identification objects, namely 0~1s, 0.5~1.5s and 1~2s. Among them, 0.5 to 1.5s includes the time node when the mode changes, which is used to verify the ability of the present invention to identify complex oscillation modes.
步骤8处理待测振荡信号并输入CNN模型,根据网络输出分析振荡类型的辨识结果,如下表1所示:Step 8: Process the oscillation signal to be measured and input it into the CNN model. Analyze the identification results of the oscillation type based on the network output, as shown in Table 1 below:
表1.基于实施例1的振荡类型的辨识结果表Table 1. Table of identification results of oscillation types based on Embodiment 1
参考上表1,可以发现该方法的辨识结果与理想振荡信号的真实类型完全吻合。当所辨识的振荡信号涉及新振荡类型的加入,该方法依旧适用。此外,该方法所需辨识时间仅1s,远快于传统振荡信号辨识方法,因此从各方面考虑,该方法都更具优势。Referring to Table 1 above, it can be found that the identification results of this method are completely consistent with the true type of the ideal oscillation signal. This method is still applicable when the identified oscillation signal involves the addition of a new oscillation type. In addition, the identification time required by this method is only 1 second, which is much faster than the traditional oscillation signal identification method. Therefore, this method has more advantages from all aspects.
实施例2:Example 2:
为了验证本发明的实际辨识效果,从电力系统中获取一段振荡信号实测数据。如图5所示,该段振荡信号由一处小扰动所激发,时间位于4s处。为了辨识该扰动所激发的振荡类型,因此截取扰动结束后的数据作为本实施例中的待测振荡信号。In order to verify the actual identification effect of the present invention, a section of actual measured data of oscillation signals is obtained from the power system. As shown in Figure 5, the oscillation signal in this section is excited by a small disturbance, and the time is located at 4s. In order to identify the type of oscillation excited by the disturbance, the data after the disturbance is intercepted is used as the oscillation signal to be measured in this embodiment.
步骤7:通过滑动时窗获取待测振荡信号样本。滑窗长度为1s,滑动间隔为0.5s,采样频率为400Hz。本实施例中,于扰动结束后滑取两段信号作为辨识对象,分别为4~5s与6~7s,以验证振荡类型是否会发生改变。Step 7: Obtain the oscillation signal sample to be measured through the sliding time window. The length of the sliding window is 1s, the sliding interval is 0.5s, and the sampling frequency is 400Hz. In this embodiment, after the perturbation ends, two sections of signals are slided as identification objects, 4 to 5 s and 6 to 7 s respectively, to verify whether the oscillation type will change.
步骤8:处理待测振荡信号并输入CNN模型,根据网络输出分析振荡类型的辨识结果,如下表2所示:Step 8: Process the oscillation signal to be measured and input it into the CNN model, and analyze the identification results of the oscillation type based on the network output, as shown in Table 2 below:
表2.基于实施例2的振荡类型的辨识结果表Table 2. Identification results of oscillation types based on Embodiment 2
参考上表2,可以判断该系统中由扰动所激发的振荡类型为低频振荡。由于采用滑窗采样,该方法可以伴随振荡信号的产生不断实时返回所产生的振荡类型,判断是否存在新型振荡的加入。此外,1s的识别时间使得该方法能迅速获取辨识结果,因此相较于传统振荡信号辨识方法,本发明更具实用价值。Referring to Table 2 above, it can be judged that the type of oscillation excited by disturbance in this system is low-frequency oscillation. Due to the use of sliding window sampling, this method can continuously return the generated oscillation type in real time as the oscillation signal is generated, and determine whether there is a new type of oscillation added. In addition, the identification time of 1 second enables this method to quickly obtain identification results. Therefore, compared with the traditional oscillation signal identification method, the present invention has more practical value.
以上所述的具体实施方式,对本发明的目的、技术方案和有益效果进行了进一步详细说明,所应理解的是,以上所述仅为本发明的具体实施方式而已,并不用于限定本发明的保护范围,凡在本发明的精神和原则之内,所做的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。The above-described specific embodiments further describe the objectives, technical solutions and beneficial effects of the present invention in detail. It should be understood that the above-mentioned are only specific embodiments of the present invention and are not intended to limit the scope of the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention shall be included in the protection scope of the present invention.
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011050933.3ACN112200038B (en) | 2020-09-29 | 2020-09-29 | A fast identification method for power system oscillation types based on CNN |
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011050933.3ACN112200038B (en) | 2020-09-29 | 2020-09-29 | A fast identification method for power system oscillation types based on CNN |
| Publication Number | Publication Date |
|---|---|
| CN112200038A CN112200038A (en) | 2021-01-08 |
| CN112200038Btrue CN112200038B (en) | 2023-12-05 |
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011050933.3AActiveCN112200038B (en) | 2020-09-29 | 2020-09-29 | A fast identification method for power system oscillation types based on CNN |
| Country | Link |
|---|---|
| CN (1) | CN112200038B (en) |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113486965A (en)* | 2021-07-14 | 2021-10-08 | 西南交通大学 | Training method for abnormity identification model of vehicle network electric coupling data |
| CN115133522A (en)* | 2022-05-24 | 2022-09-30 | 沈阳工业大学 | Fan grid-connected subsynchronous oscillation monitoring network construction and layered reconstruction suppression method |
| CN115187154B (en)* | 2022-09-14 | 2022-12-16 | 华中科技大学 | A Neural Network-Based Risk Prediction Method and System for Regional Power Grid Oscillation Sources |
| CN115377999B (en)* | 2022-09-26 | 2024-07-26 | 华北电力大学 | Subsynchronous oscillation identification method based on broadband measurement data |
| CN115659220B (en)* | 2022-10-25 | 2025-09-30 | 国网河北省电力有限公司电力科学研究院 | SSCI oscillation type identification method, model training method and equipment |
| CN117350170B (en)* | 2023-11-20 | 2024-02-09 | 华北电力大学(保定) | Nonlinear oscillation analysis method based on KOOPMAN deep neural network |
| CN117892130B (en)* | 2024-01-16 | 2025-03-14 | 云南大学 | Industrial process oscillation detection method and system based on lightweight convolutional neural network |
| CN118194239B (en)* | 2024-05-15 | 2024-09-03 | 合肥工业大学 | Multi-modal graph network-based broadband oscillation parameter identification method for power transmission system |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1981453A (en)* | 2004-06-17 | 2007-06-13 | W5网络公司 | Low power wireless communication system and protocol |
| CN104578115A (en)* | 2015-01-26 | 2015-04-29 | 国网四川省电力公司经济技术研究院 | Electric system low frequency oscillation mode identification method based on correlation functions |
| CN105473065A (en)* | 2014-02-11 | 2016-04-06 | 皇家飞利浦有限公司 | Determining return of spontaneous circulation during CPR |
| CN109255394A (en)* | 2018-10-18 | 2019-01-22 | 国网天津市电力公司电力科学研究院 | A kind of forced oscillation recognition methods based on Pattern similarity |
| CN110236536A (en)* | 2019-06-04 | 2019-09-17 | 电子科技大学 | A detection system of EEG high frequency oscillation signal based on convolutional neural network |
| CN110261746A (en)* | 2019-07-08 | 2019-09-20 | 清华大学深圳研究生院 | Electric cable stoppage detection method based on oscillating wave voltage periodic attenuation characteristic |
| CA3105412A1 (en)* | 2018-07-06 | 2020-01-09 | Wobben Properties Gmbh | Wind energy system and method for identifying low-frequency oscillations in an electrical supply network |
| CN111046327A (en)* | 2019-12-18 | 2020-04-21 | 河海大学 | Prony analysis method suitable for low-frequency oscillation and subsynchronous oscillation identification |
| CN111160167A (en)* | 2019-12-18 | 2020-05-15 | 北京信息科技大学 | Spindle fault classification and identification method based on S-transform deep convolutional neural network |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7108192B2 (en)* | 1999-09-17 | 2006-09-19 | Silverbrook Research Pty Ltd | Rotationally symmetric tags |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1981453A (en)* | 2004-06-17 | 2007-06-13 | W5网络公司 | Low power wireless communication system and protocol |
| CN105473065A (en)* | 2014-02-11 | 2016-04-06 | 皇家飞利浦有限公司 | Determining return of spontaneous circulation during CPR |
| CN104578115A (en)* | 2015-01-26 | 2015-04-29 | 国网四川省电力公司经济技术研究院 | Electric system low frequency oscillation mode identification method based on correlation functions |
| CA3105412A1 (en)* | 2018-07-06 | 2020-01-09 | Wobben Properties Gmbh | Wind energy system and method for identifying low-frequency oscillations in an electrical supply network |
| CN109255394A (en)* | 2018-10-18 | 2019-01-22 | 国网天津市电力公司电力科学研究院 | A kind of forced oscillation recognition methods based on Pattern similarity |
| CN110236536A (en)* | 2019-06-04 | 2019-09-17 | 电子科技大学 | A detection system of EEG high frequency oscillation signal based on convolutional neural network |
| CN110261746A (en)* | 2019-07-08 | 2019-09-20 | 清华大学深圳研究生院 | Electric cable stoppage detection method based on oscillating wave voltage periodic attenuation characteristic |
| CN111046327A (en)* | 2019-12-18 | 2020-04-21 | 河海大学 | Prony analysis method suitable for low-frequency oscillation and subsynchronous oscillation identification |
| CN111160167A (en)* | 2019-12-18 | 2020-05-15 | 北京信息科技大学 | Spindle fault classification and identification method based on S-transform deep convolutional neural network |
| Title |
|---|
| A Convolution Neural Network Method for Power System Oscillation Type Identification;Qianyu Li等;《2020 IEEE 4th Conference on Energy Internet and Energy System Integration (EI2)》;499-504* |
| A Novel Neural Network Approach for Power System Low Frequency Oscillation Mode Identification;Zhongting Shen等;《2019 IEEE International Symposium on Circuits and Systems (ISCAS)》;1-5* |
| A Robust Direct Parameter Identification of Exponentially Damped Low-Frequency Oscillation in Power Systems;Zhaobi Chu等;《Journal of Sensors》;1-11* |
| 基于小波变换和神经网络的暂态电能质量扰动自动识别;刘晓芳等;《继电器》(第23期);46-50* |
| 基于改进贝叶斯分类法的电能质量扰动分类方法;张文涛等;《电网技术》(第07期);22-25* |
| 基于神经网络的电力负荷预测方法研究;罗宁等;《自动化与仪器仪表》(第01期);157-160* |
| 多直驱风机经VSC-HVDC并网系统场内/场网次同步振荡特性分析;邵冰冰等;《中国电机工程学报》;第40卷(第12期);3835-3847* |
| Publication number | Publication date |
|---|---|
| CN112200038A (en) | 2021-01-08 |
| Publication | Publication Date | Title |
|---|---|---|
| CN112200038B (en) | A fast identification method for power system oscillation types based on CNN | |
| CN112183368B (en) | LSTM-based rapid identification method for low-frequency oscillation modal characteristics of power system | |
| CN112751345B (en) | Low Frequency Oscillation Mode Identification Method of Power System Based on LSTM and Phase Trajectory | |
| CN110994604B (en) | Power system transient stability assessment method based on LSTM-DNN model | |
| CN111478314B (en) | A Method for Evaluating Transient Stability of Power System | |
| CN108832619A (en) | Power System Transient Stability Evaluation Method Based on Convolutional Neural Network | |
| CN109948833A (en) | A Deterioration Trend Prediction Method of Hydroelectric Units Based on Long Short-Term Memory Network | |
| CN110675712B (en) | A practical training system for power cable oscillating wave partial discharge detection | |
| US11144023B1 (en) | Method for PMU data recovery using an improved cubic spline interpolation and singular value decomposition | |
| CN110176762A (en) | A kind of sub-synchronous oscillation risk online evaluation method and device based on frequency domain impedance | |
| US11170304B1 (en) | Bad data detection algorithm for PMU based on spectral clustering | |
| CN113612237A (en) | Method for positioning resonance-induced subsynchronous oscillation source in offshore wind farm | |
| CN103400009A (en) | Wind electric field dynamic equivalence method based on split level semi-supervised spectral clustering algorithm | |
| CN114049014A (en) | Method, device and system for evaluating operating status of offshore wind turbines | |
| Han et al. | A novel wind farm equivalent model for high voltage ride through analysis based on multi-view incremental transfer clustering | |
| CN113659565A (en) | An online prediction method of frequency situation of new energy power system | |
| CN116505522B (en) | Simulation method, simulation platform and equipment for operation of power system | |
| Chakravorti et al. | Advanced signal processing techniques for multiclass disturbance detection and classification in microgrids | |
| CN111553112A (en) | A method and device for fault identification of power system based on deep belief network | |
| CN112329535B (en) | CNN-based quick identification method for low-frequency oscillation modal characteristics of power system | |
| Yan et al. | YOLOV4-based wind turbine blade crack defect detection | |
| Sudhakar et al. | Faulty diagnostics model for wind power plant application using AI | |
| CN114897006A (en) | Smart meter error classification method and system, equipment and storage medium | |
| CN117390418A (en) | A method, system and equipment for transient stability assessment of wind power grid-connected system | |
| CN117540863A (en) | Photovoltaic power interval prediction method based on photovoltaic output fluctuation rate parting |
| 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 |