CN103760768A - 基于正常运行数据的发电机组多变量系统辨识方法 - Google Patents

Abstract

本发明提供一种基于正常运行数据的发电机组多变量系统辨识方法，对正常运行时的发电机组多变量系统进行阶跃响应测试，获得发电机组多变量系统的各个回路的输入输出信号；进行信号分解和频谱分析，获取发电机组多变量系统的所有回路中的最大临界频率；根据最大临界频率计算发电机组多变量系统在多个频率点上的频率特性；根据频率特性构建阵发电机组多变量系统的传递函数矩阵模型。本发明的技术方案，通过将一个强耦合的发电机组多变量系统的辨识过程分解成多个单输入单输出系统的辨识过程，通过系统正常运行时的阶跃响应测试，在发电机组多变量系统正常运行状态下确定系统在重要频率段的频率特性，辨识过程简单、抗干扰能力强、辨识精度高。

Description

基于正常运行数据的发电机组多变量系统辨识方法

技术领域

本发明涉及多变量系统应用技术领域，特别是涉及一种基于正常运行数据的发电机组多变量系统辨识方法。

背景技术

现代火电厂的大型单元机组具有大时滞、时变、非线性、强耦合、机炉动态特性差异大的特点，是一个典型的多输入多输出过程对象。单元机组协调控制系统是电厂自动控制系统中最复杂的系统，它的任务是使机组负荷紧密跟踪外界负荷需求，并保持机前压力的稳定。机组采用协调控制系统运行时，将锅炉、汽机和发电机作为一个统一整体进行控制。

目前的大多数控制方法是在传统的PID控制基础上综合了锅炉和汽轮机的特点，提出各种控制方案，如前馈控制、解耦控制、内模控制、预测控制等。由于机组协调被控对象具有非线性，因而在不同负荷段，系统的动态特性是不同的。即使在标称负荷点得到了满意的控制器参数,随着负荷点的变化，系统性能也会变得很差。协调控制器的参数整定和负荷适应性问题是一个具有普遍性和现实意义的问题。

学者李少远等结合自整定技术、自适应技术、增益调度方法、给出了智能解耦的控制策略，学者谭文等从经典解耦控制的角度对协调控制器进行了设计，并给出了相应的参数整定方法。

虽然这些先进控制算法具有良好的跟踪性能和稳定性，但由于针对现有系统闭环辨识技术，在辨识单元机组的模型传递函数时，需要引入外加测试信号来辨识出多变量过程模型，对过程装置正常运行的带来了扰动，导致了多产品质量的影响，降低了辨识精度，对于模型变化以及非线性的适应性还是有局限性。

发明内容

基于此，有必要针对现有系统闭环辨识技术辨识精度低的问题，提供一种基于正常运行数据的发电机组多变量系统辨识方法。

一种基于正常运行数据的发电机组多变量系统辨识方法，包括如下步骤：

对正常运行时的发电机组多变量系统进行阶跃响应测试，获得所述发电机组多变量系统的各个回路的输入输出信号；

对所述输入输出信号进行信号分解和频谱分析，获取发电机组多变量系统的所有回路中的最大临界频率；

根据所述最大临界频率计算发电机组多变量系统在多个频率点上的频率特性；

根据所述频率特性构建阵发电机组多变量系统的传递函数矩阵模型。

上述基于正常运行数据的发电机组多变量系统辨识方法，通过将一个强耦合的发电机组多变量系统的辨识过程分解成多个单输入单输出系统的辨识过程，通过系统正常运行时的阶跃响应测试，对系统的输入输出信号进行信号分解和频谱分析，得到全部回路最大临界频率并获取到多个频率点的频率特性，进而获得该发电机组多变量系统的传递函数矩阵模型。该技术方案无需引入外加测试信号，在发电机组多变量系统正常运行状态下确定系统在重要频率段的频率特性，辨识过程简单、抗干扰能力强、辨识精度高，为发电机组多变量系统的协调控制器的参数自适应、克服机组非线性特性奠定了基础。

附图说明

图1为一个实施例的基于正常运行数据的发电机组多变量系统辨识方法的流程图；

图2为机组协调控制系统框图。

具体实施方式

下面结合附图对本发明的基于正常运行数据的发电机组多变量系统辨识方法的具体实施方式作详细描述。

参考图1所示，图1为一个实施例的基于正常运行数据的发电机组多变量系统辨识方法的流程图，包括如下步骤：

步骤S1：对正常运行时的发电机组多变量系统进行阶跃响应测试，获得所述发电机组多变量系统的各个回路的输入输出信号。

在一个实施例中，步骤S1的阶跃响应测试过程包括如下：

在发电机组多变量系统正常运行时，依次对各个回路的参考输入施加设定值的阶跃测试信号，并保持其它参考输入不变，获取每次阶跃输入测试时各个回路对象的输入输出信号。

具体的，阶跃响应测试过程可以进一步如下：

（1）对第一个回路参考输入r1施加阶跃测试信号，记录此时多输入多输出对象的输入输出信号

（2）对第二个回路参考输入r2施加阶跃测试信号，记录此时多输入多输出对象的输入输出信号

（3）以此类推，对余下各项参考输入依次进行阶跃信号测试，直到记录第M个回路参考输入信号阶跃变化后的对象的输入输出信号

步骤S2：对所述输入输出信号进行信号分解和频谱分析，获取发电机组多变量系统的所有回路中的最大临界频率。

在一个实施例中，步骤S2的获取最大临界频率过程包括如下：

首先，对所述输入输出信号进行信号进行分解成瞬态部分和稳态部分；获取发电机组多变量系统的系统时域响应等式。然后根据所述系统时域响应等式进行频谱分析，获得各个回路频率点的系统频率响应等式。并根据所述频率响应等式计算各个回路频率点对应的相位，并根据所述相位和迭代公式计算获得发电机组多变量系统的所有回路中的最大临界频率的最大值。

具体的，获取最大临界频率过程进一步可以包括如下：

将阶跃测试获取的信号都分解为瞬态部分和稳态部分，假定先改变给定值r₁，而保持其它给定值不变，Y¹(s)和U¹(s)可表示为：

$Y^1\left(s\right)=\left[\begin{array}{c}{\∫}_0^{{t_{y1}}^1}\mathrm{\Δ}{y_1}^1\left(t\right)e^{-\mathrm{st}}\mathrm{dt}+\frac{{y_1}^1(\∞)}{s}\\ {\∫}_0^{{t_{y2}}^1}\mathrm{\Δ}{y_2}^1\left(t\right)e^{-\mathrm{st}}\mathrm{dt}+\frac{{y_2}^1(\∞)}{s}\end{array}\right]$

$U^1\left(s\right)=\left[\begin{array}{c}{\∫}_0^{{t_{u1}}^1}\mathrm{\Δ}{u_1}^1\left(t\right)e^{-\mathrm{st}}\mathrm{dt}+\frac{{u_1}^1(\∞)}{s}\\ {\∫}_0^{{t_{u2}}^1}\mathrm{\Δ}{u_2}^1\left(t\right)e^{-\mathrm{st}}\mathrm{dt}+\frac{{u_2}^1(\∞)}{s}\end{array}\right]$

它们满足如下关系：

Y¹(s)=G(s)U¹(s)

而后保持其它给定值不变，然后改变给定值r₂，保持其他给定值不变，以此类推，最后获得下面等式：

[Y¹(s)…Y^m(s)]=G(s)[U¹(s)…U^m(s)]

式中，G(s)为系统时域响应函数，s表示拉普拉斯变量。

将s=jω代入上式，可获得发电机组多变量系统的各个回路频率点的频率响应公式：

由上式可得到各个频率的相位，再通过迭代公式：

${\mathrm{\ω}}_{n+1}={\mathrm{\ω}}_n-(\mathrm{\π}+{\mathrm{\φ}}_n)\frac{{\mathrm{\ω}}_n-{\mathrm{\ω}}_{n-1}}{{\mathrm{\φ}}_n-{\mathrm{\φ}}_{n-1}}$

即可确定发电机组多变量系统所有回路临界频率中的最大值ω_max；式中，ω_n为频率点n的频率，φ_n为频率点n对应的相位。

具体的，将ω_n-1和φ_n-1的初始值设为零，而ω_n取为一个尽量小的数，如10^-3。迭代运算公式具有二次收敛速度，在几个迭代运算后，ω_max就可获得至少99%的准确度。

步骤S3：根据所述最大临界频率计算发电机组多变量系统在多个频率点上的频率特性。

在一个实施例中，步骤S3计算频率特性的方法包括如下：

首先根据需要辨识的频率响应点的数目M，将频率范围(0,ω_max)划分为M个离散频率点，得到频率点0,Δω,2Δω,…,(M-1)Δω，其中Δω=ω_max/M；然后根据所述系统频率响应等式计算各个离散频率点的频率特性G(jω_l)，包括幅值和相角，其中，l=1,2,…M。

步骤S4：根据所述频率特性构建阵发电机组多变量系统的传递函数矩阵模型。

在一个实施例中，步骤S4的构建阵发电机组多变量系统的传递函数矩阵模型的方法包括如下：

根据所述频率特性在各个频率响应匹配获得二价加纯滞后模型的参数，其中，传递函数矩阵模型的元素为二价加纯滞后模型；利用最小二乘法并根据所述二价加纯滞后模型的参数获取发电机组多变量系统的加纯滞后的传递函数矩阵模型。

具体的，即通过前述步骤得到的多变量过程对象的频率响应G(jω_l),l=1,2,…M，可获取传递函数矩阵。进一步可以包括如下：

传递函数矩阵的每个元素采用二阶加纯滞后模型：

$g^{\′}\left(s\right)=\frac{1}{{\mathrm{as}}^2+\mathrm{bs}+c}{e}^{-\mathrm{Ls}},$

上述模型可以表示单调的、振荡的和非最小相位过程对象，二阶加纯滞后模型的参数可通过在ω_l,l=1,2,...,M的频率响应点匹配g'(jω)和g(jω)获取，即

$g\left({\mathrm{j\ω}}_l\right)=\frac{l_0}{{\left({\mathrm{j\ω}}_l\right)}^{2}a+j{\mathrm{\ω}}_{l}b+c}{e}^{-j{\mathrm{\ω}}_{l}L},l=\mathrm{1,2},...,M.$

上式中的参数a,b,c和频率点数L可通过频率特性G(jω_l)的幅值条件和相位条件确定。

$\begin{array}{c}{{\mathrm{\ω}}_l}^4{\left|g\left({\mathrm{j\ω}}_l\right)\right|}^2{a}^2+{{\mathrm{\ω}}_l}^2{\left|g\left({\mathrm{j\ω}}_l\right)\right|}^2(b^2-2\mathrm{ac})+{\left|g\left({\mathrm{j\ω}}_l\right)\right|}^2{c}^2=1\\ -\arg \left[g\left({\mathrm{j\ω}}_l\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_l}{c-{{\mathrm{a\ω}}_l}^2}\right)={\mathrm{\ω}}_{l}L\end{array}$

当l=1,2,...,M时，上述等式可写成矩阵形式：

Φθ＝Γ

其中

$\mathrm{\Φ}=\left[\begin{array}{ccc}{{\mathrm{\ω}}_1}^4{\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2& {{\mathrm{\ω}}_1}^2{\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2& {\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2\\ {{\mathrm{\ω}}_2}^4{\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2& {{\mathrm{\ω}}_2}^2{\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2& {\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2\\ \·& \·& \·\\ \·& \·& \·\\ \·& \·& \·\\ {{\mathrm{\ω}}_M}^4{\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2& {{\mathrm{\ω}}_M}^2{\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2& {\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2\end{array}\right]$

$\mathrm{\Γ}=\left[\begin{array}{c}{\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2\\ {\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2\\ \·\\ \·\\ \·\\ {\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2\end{array}\right]$

$\mathrm{\θ}=\left[\begin{array}{c}{\mathrm{\θ}}_1\\ {\mathrm{\θ}}_2\\ {\mathrm{\θ}}_3\end{array}\right]=\left[\begin{array}{c}{a}^2\\ b^2-2\mathrm{ac}\\ c^2\end{array}\right].$

式中，θ可通过最小二乘法得到

θ＝(Φ^TΦ)^-1Φ^TΓ

则模型参数可从θ中通过下式获得：

$\left[\begin{array}{c}a\\ b\\ c\end{array}\right]=\left[\begin{array}{c}\sqrt{{\mathrm{\θ}}_1}\\ \sqrt{{\mathrm{\θ}}_2+2\sqrt{{\mathrm{\θ}}_1{\mathrm{\θ}}_3}}\\ \sqrt{{\mathrm{\θ}}_3}\end{array}\right]$

参数a,b,c确定后，纯滞后L最终可通过最小二乘法获得：

$\left[\begin{array}{c}{\mathrm{\ω}}_1\\ {\mathrm{\ω}}_2\\ \·\\ \·\\ \·\\ {\mathrm{\ω}}_M\end{array}\right]L=\left[\begin{array}{c}-\arg \left[g\left({\mathrm{j\ω}}_1\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_1}{c-{{\mathrm{a\ω}}_1}^2}\right)\\ -\arg \left[g\left({\mathrm{j\ω}}_2\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_2}{c-{{\mathrm{a\ω}}_2}^2}\right)\\ \·\\ \·\\ \·\\ -\arg \left[g\left({\mathrm{j\ω}}_M\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_M}{c-{{\mathrm{a\ω}}_M}^2}\right)\end{array}\right]$

通过上述获取的传递函数矩阵模型，可以用于发电机组多变量系统的自适应调整中，该模型去除了任何外加测试信号，仅需要闭环系统正常运行时的多变量过程输入输出数据，消除了对过程装置正常运行的扰动，具有高的辨识精度，从而使得发电机组多变量系统实现了精确控制。

综合上述各个实施例，本发明的技术方案，将一个强耦合的发电机组多变量系统的辨识过程分解成多个单输入单输出系统的辨识过程，通过系统正常运行时的阶跃响应测试，对系统的输入输出信号进行信号分解和频谱分析，得到全部回路最大临界频率并获取到多个频率点的频率特性，获得该多变量系统的传递函数矩阵模型。

与现有的发电机组多变量系统的辨识技术方案相比，本发明的技术方案实现简单、抗干扰能力强、辨识精度高，不需要中断发电机组多变量系统正常的运行状态进行试验，能得到发电机组多变量系统在重要频率段的频率特性，然后采用最小二乘法在幅值和相位两方面拟合出单元机组在当前负荷点上的动态模型，得到多参数的二阶加纯滞后的辨识模型。为发电机组多变量系统的协调控制器的参数自适应、克服机组非线性特性奠定了基础。

为了更加清晰本发明的技术方案，下面结合附图举出基于本发明的方案实现的一个应用实例。

参考图2所示，图2为机组协调控制系统框图，其中，K(s)为机组协调控制器，G(s)表示单元机组过程对象。P_SP、N_SP为系统输入控制量，B为燃料量、μ为调节阀开度，P_T为主蒸汽压力、N为单元机组的输出功率。

考虑一个火电单元机组协调控制系统，如图2所示，单元机组的模型传递函数为：

$G\left(s\right)=\left[\begin{array}{cc}\frac{{12.8e}^{-s}}{1+16.7s}& \frac{{-18.9e}^{-3s}}{1+21s}\\ \frac{{6.6e}^{-7s}}{1+10.9s}& \frac{{-19.4e}^{-3s}}{1+14.4s}\end{array}\right]$

具体辨识的步骤如下：

第一步，对第一个回路参考输入r1施加阶跃测试信号，记录系统在NSR=20%噪声环境下的输入输出信号

对第二个回路参考输入r2施加阶跃测试信号，记录系统在NSR=20%噪声环境下的输入输出信号

第二步，将第一步数据带入频率响应公式G(jω_n)和迭代公式求出系统所有回路最大极限频率为ω_max=0.56rad/s.

第三步，从0到ω_max划分10个频率点，依次为0、0.002、0.0044、0.0088、0.0175、0.035、0.07、0.14、0.28、0.56，根据频率响应公式G(jω_n)，依次求出该频率在各个系统回路中的幅值和相角。

第四步，将各个回路频率点的特性构建回路特性矩阵，最终得到系统的二阶加纯滞后模型矩阵：

$G\left(s\right)=\left[\begin{array}{cc}\frac{e^{-1.0082s}}{1.2631+0.078s}& \frac{{-e}^{-3.0339s}}{1.1069s+0.0526}\\ \frac{e^{-6.7246s}}{0.8587s^2+1.7491s+0.1487}& \frac{{-e}^{-2.4787s}}{0.7411s+0.051}\end{array}\right]$

通过上述应用实例，可以看出，本发明的基于正常运行数据的发电机组多变量系统辨识方法，不需要任何过程对象的先验知识，基于系统正常运行的数据在线辨识出重要频率范围内的频率响应特性，进而获得二阶加纯滞后模型的传递函数矩阵，上述实施实例表明该辨识方法对火电单元机组过程对象是有效的，且在噪声环境下辨识结果具有较高的精度。

以上所述实施例仅表达了本发明的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求为准。

Claims (7)

1.一种基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，包括如下步骤：

对正常运行时的发电机组多变量系统进行阶跃响应测试，获得所述发电机组多变量系统的各个回路的输入输出信号；

对所述输入输出信号进行信号分解和频谱分析，获取发电机组多变量系统的所有回路中的最大临界频率；

根据所述最大临界频率计算发电机组多变量系统在多个频率点上的频率特性；

根据所述频率特性构建阵发电机组多变量系统的传递函数矩阵模型。

2.根据权利要求1所述的基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，所述对正常运行时的发电机组多变量系统进行阶跃响应测试，获得所述发电机组多变量系统的各个回路的输入输出信号的步骤包括：

在发电机组多变量系统正常运行时，依次对各个回路的参考输入施加设定值的阶跃测试信号，并保持其它参考输入不变，获取每次阶跃输入测试时各个回路对象的输入输出信号。

3.根据权利要求1所述的基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，所述对所述输入输出信号进行信号分解和频谱分析，获得发电机组多变量系统的所有回路中的最大临界频率的步骤包括：

对所述输入输出信号进行信号进行分解成瞬态部分和稳态部分；获取发电机组多变量系统的系统时域响应等式；

根据所述系统时域响应等式进行频谱分析，获得各个回路频率点的系统频率响应等式；

根据所述频率响应等式计算各个回路频率点对应的相位，并根据所述相位和迭代公式计算获得发电机组多变量系统的所有回路中的最大临界频率的最大值。

4.根据权利要求3所述的基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，所述系统时域响应等式如下：

[Y¹(s)…Y^m(s)]=G(s)[U¹(s)…U^m(s)]

式中，G(s)为系统时域响应函数；

所述系统频率响应等式如下：

式中，G(jω)为系统频率响应函数，s=jω；

所述迭代公式如下：

${\mathrm{\ω}}_{n+1}={\mathrm{\ω}}_n-(\mathrm{\π}+{\mathrm{\φ}}_n)\frac{{\mathrm{\ω}}_n-{\mathrm{\ω}}_{n-1}}{{\mathrm{\φ}}_n-{\mathrm{\φ}}_{n-1}},$

φ_n=min{arg[G(jω_n)]}

式中，ω_n为频率点n的频率，φ_n为频率点n对应的相位。

5.根据权利要求3或4所述的基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，所述根据所述最大临界频率计算发电机组多变量系统在多个频率点上的频率特性的步骤包括：

根据需要辨识的频率响应点的数目M，将频率范围(0,ω_max)划分为M个离散频率点，得到频率点0,Δω,2Δω,…,(M-1)Δω，其中Δω=ω_max/M；

根据所述系统频率响应等式计算各个离散频率点的频率特性G(jω_l)，包括幅值和相角，其中，l=1,2,…M。

6.根据权利要求1所述的基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，所述根据所述频率特性构建阵发电机组多变量系统的传递函数矩阵模型的步骤包括：

根据所述频率特性在各个频率响应匹配获得二价加纯滞后模型的参数，其中，传递函数矩阵模型的元素为二价加纯滞后模型；

利用最小二乘法并根据所述二价加纯滞后模型的参数获取发电机组多变量系统的加纯滞后的传递函数矩阵模型。

7.根据权利要求6所述的基于正常运行数据的发电机组多变量系统辨识方法，其特征在于，所述获取加纯滞后的传递函数矩阵模型的方法包括如下：

$g\left({\mathrm{j\ω}}_l\right)=\frac{l_0}{{\left({\mathrm{j\ω}}_l\right)}^{2}a+j{\mathrm{\ω}}_{l}b+c}{e}^{-j{\mathrm{\ω}}_{l}L},l=\mathrm{1,2},...,M.$

其中

$\begin{array}{c}{{\mathrm{\ω}}_l}^4{\left|g\left({\mathrm{j\ω}}_l\right)\right|}^2{a}^2+{{\mathrm{\ω}}_l}^2{\left|g\left({\mathrm{j\ω}}_l\right)\right|}^2(b^2-2\mathrm{ac})+{\left|g\left({\mathrm{j\ω}}_l\right)\right|}^2{c}^2=1\\ -\arg \left[g\left({\mathrm{j\ω}}_l\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_l}{c-{{\mathrm{a\ω}}_l}^2}\right)={\mathrm{\ω}}_{l}L\end{array}$

当l=1,2,...,M时，传递函数矩阵模型的矩阵形式如下：

Φθ＝Γ

其中

$\mathrm{\Φ}=\left[\begin{array}{ccc}{{\mathrm{\ω}}_1}^4{\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2& {{\mathrm{\ω}}_1}^2{\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2& {\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2\\ {{\mathrm{\ω}}_2}^4{\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2& {{\mathrm{\ω}}_2}^2{\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2& {\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2\\ \·& \·& \·\\ \·& \·& \·\\ \·& \·& \·\\ {{\mathrm{\ω}}_M}^4{\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2& {{\mathrm{\ω}}_M}^2{\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2& {\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2\end{array}\right]$

$\mathrm{\Γ}=\left[\begin{array}{c}{\left|g\left({\mathrm{j\ω}}_1\right)\right|}^2\\ {\left|g\left({\mathrm{j\ω}}_2\right)\right|}^2\\ \·\\ \·\\ \·\\ {\left|g\left({\mathrm{j\ω}}_M\right)\right|}^2\end{array}\right]$

$\mathrm{\θ}=\left[\begin{array}{c}{\mathrm{\θ}}_1\\ {\mathrm{\θ}}_2\\ {\mathrm{\θ}}_3\end{array}\right]=\left[\begin{array}{c}{a}^2\\ b^2-2\mathrm{ac}\\ c^2\end{array}\right]$

式中的θ计算如下：

θ=(Φ^TΦ)^-1Φ^TΓ

$\left[\begin{array}{c}a\\ b\\ c\end{array}\right]=\left[\begin{array}{c}\sqrt{{\mathrm{\θ}}_1}\\ \sqrt{{\mathrm{\θ}}_2+2\sqrt{{\mathrm{\θ}}_1{\mathrm{\θ}}_3}}\\ \sqrt{{\mathrm{\θ}}_3}\end{array}\right]$

$\left[\begin{array}{c}{\mathrm{\ω}}_1\\ {\mathrm{\ω}}_2\\ \·\\ \·\\ \·\\ {\mathrm{\ω}}_M\end{array}\right]L=\left[\begin{array}{c}-\arg \left[g\left({\mathrm{j\ω}}_1\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_1}{c-{{\mathrm{a\ω}}_1}^2}\right)\\ -\arg \left[g\left({\mathrm{j\ω}}_2\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_2}{c-{{\mathrm{a\ω}}_2}^2}\right)\\ \·\\ \·\\ \·\\ -\arg \left[g\left({\mathrm{j\ω}}_M\right)\right]-{\tan }^{-1}\left(\frac{{\mathrm{b\ω}}_M}{c-{{\mathrm{a\ω}}_M}^2}\right)\end{array}\right]$

式中，θ为二价加纯滞后模型的参数。

Images