CN105960614A - 系统辨识装置 - Google Patents

Abstract

系统辨识装置针对属于被指定的系统维数的搜索范围的各个维数，通过递归方法来辨识线性离散时间系统，计算对与各维数对应的线性离散时间系统应用了实际的辨识用输入数据的情况下的系统输出而作为系统特性，将该系统输出和动态系统的实际的辨识用输出数据在时域中的误差平方和的范数分布(41)在阈值(42)以下的维数中的最小的维数决定为系统维数n，并且根据动态系统的输入输出矢量、以及使用所决定的系统维数而生成的状态矢量，辨识线性离散时间系统的系统矩阵。

Description

系统辨识装置

技术领域

本发明涉及根据对作为对象的动态系统施加了伪随机输入的情况下的系统输入输出来构筑该系统的数学模型的系统辨识装置。

背景技术

作为以往的基于伪随机输入的系统辨识装置，有基于例如非专利文献1记载的N4SID法的装置。在该N4SID法中，根据在对记述于线性离散时间系统(A_d，B_d，C_d，D_d)的动态系统施加了伪随机输入的情况下的系统输入输出，生成与系统输入有关的块汉克尔(Hankel)矩阵(U_p，U_f)以及与系统输出有关的块汉克尔矩阵(Y_p，Y_f)，进而根据块汉克尔矩阵(U_f，Y_f)，生成输入输出矢量(^～U_K|K，^～Y_K|K)。另外，关于“^～”的标记，本来应该是在“U”的字符的上部划横线(上划线(overbar))，但无法实现该标记。因此，在本说明书中，除了用图片插入的公式部分以外，用“^～”来代替横线(上划线)。

接下来，对结合上述块汉克尔矩阵而得到的数据矩阵进行LQ分解，根据通过LQ分解而得到的子矩阵和块汉克尔矩阵U_p、Y_p，生成平行投影Θ。对该平行投影Θ进行奇异值分解，将具有有效值的奇异值的个数决定为系统维数，根据奇异值分解的结果和所决定的系统维数，计算动态系统的状态矢量(^～X_K，^～X_K+1)。最后，通过对输入输出矢量(^～U_K|K，^～Y_K|K)以及状态矢量(^～X_K，^～X_K+1)应用最小二乘法，辨识记述动态系统的线性离散时间系统(A_d，B_d，C_d，D_d)。

另外，作为以往的基于伪随机输入的系统辨识装置的其他例子，有例如专利文献1记载的曝光装置和除振装置、系统辨识装置及其方法。

在该曝光装置和除振装置、系统辨识装置及其方法中，根据在对作为对象的动态系统施加了伪随机输入的情况下的系统输入输出，通过以N4SID法为代表的子空间法，辨识动态系统的状态方程式。此时，通过使要辨识的状态方程式的系统维数与根据动态系统的运动方程式决定的系统维数一致，利用基于运动方程式的特性方程式和基于辨识出的状态方程式的特性方程式的比较，辨识运动方程式所包含的未知的物理参数。

专利文献1：日本特开2000-82662号公报

非专利文献1：システム同定－部分空間法からのアプローチ(系统辨识-根据子空间法的方法)，朝仓书店pp.117-120

发明内容

在这样的基于伪随机输入的系统辨识装置中，根据具有有效值的奇异值的个数、或者根据从动态系统的运动方程式决定的系统维数，决定作为对象的动态系统的系统维数。

但是，根据现实的系统输入输出计算出的平行投影Θ的奇异值为平滑的单调减少的情况较多，在该情况下，具有有效值的奇异值和成为可忽略的微小的值的奇异值的边界变得不明确。因此，在非专利文献1记载的以往的系统辨识装置中，存在如下问题：系统维数的决定依赖于作业者的判断，未必始终决定出出最佳的系统维数、或者关于系统维数的决定需要反复试验。

另外，在通过动态系统的建模得到的运动方程式中，难以记载动态系统的全部实际的动态特性，一般，“根据运动方程式决定的系统维数<动态系统的实际的系统维数”。因此，在专利文献1记载的以往的系统辨识装置中，原本存在无法决定记述动态系统的最佳的系统维数这样的问题。

另外，在以往的基于伪随机输入的系统辨识装置中，作为辨识结果得到的线性离散时间系统(A_d，B_d，C_d，D_d)的稳定性完全未被考虑，所以还存在尽管现实的动态系统是稳定的，有时还是被辨识为不稳定系统这样的问题。

本发明是鉴于上述情况而完成的，其目的在于得到一种系统辨识装置，即使在根据现实的系统输入输出计算出的平行投影Θ的奇异值为平滑的单调减少，从而具有有效值的奇异值和成为可忽略的微小的值的奇异值的边界变得不明确的情况下，也能够从系统维数的决定来排除反复试验，决定最佳的系统维数。

另外，本发明的目的在于得到一种系统辨识装置，在现实的动态系统是稳定的情形明确的情况下，能够限定于稳定系统而进行辨识。

为了解决上述课题并实现目的，本发明的系统辨识装置将对作为辨识对象的动态系统施加了伪随机输入的情况下的系统输入输出以及被指定的系统维数的搜索范围作为输入，该系统辨识装置的特征在于，具备：系统输入输出提取部，从所述动态系统的系统输入输出提取应用于辨识的辨识用输入输出数据；块汉克尔矩阵生成部，根据所述辨识用输入输出数据，生成块汉克尔矩阵；输入输出矢量生成部，根据所述块汉克尔矩阵，生成所述动态系统的输入矢量以及输出矢量；LQ分解部，结合所述块汉克尔矩阵而生成数据矩阵，输出对该数据矩阵进行LQ分解而得到的子矩阵；平行投影生成部，根据所述子矩阵和所述块汉克尔矩阵，生成平行投影；奇异值分解部，通过所述平行投影的奇异值分解，输出将所述平行投影的奇异矢量作为列矢量的第1正交矩阵、将该平行投影的右奇异矢量作为列矢量的第2正交矩阵以及该平行投影的奇异值；系统维数决定部，根据所述第2正交矩阵及所述奇异值、所述动态系统的输入矢量及输出矢量、和所述搜索范围，针对属于该搜索范围的各个维数，辨识记述动态系统的线性离散时间系统的系统矩阵，进而根据基于该系统矩阵计算出的线性离散时间系统的系统特性和动态系统的实际的系统特性的比较，决定系统维数；状态矢量生成部，根据所述第2正交矩阵及奇异值、和所决定的所述系统维数，生成所述动态系统的状态矢量；和系统矩阵辨识部，根据所述动态系统的输入矢量及输出矢量以及所述动态系统的状态矢量，辨识记述该动态系统的线性离散时间系统的系统矩阵，所述系统辨识装置将辨识出的所述系统矩阵作为记述所述动态系统的线性离散时间系统输出。

根据本发明，在作为辨识对象的动态系统中，即使在根据现实的系统输入输出计算出的平行投影的奇异值为平滑的单调减少，从而具有有效值的奇异值和成为可忽略的微小的值的奇异值的边界变得不明确的情况下，也能够从系统维数的决定来排除反复试验，始终进行最佳的系统维数的决定、和记述动态系统的线性离散时间系统的辨识。

附图说明

图1是示出实施方式1以及实施方式2的系统辨识装置的整体结构的框图。

图2是示出实施方式1的系统辨识装置中的系统输入输出的时间波形的概略图。

图3是示出实施方式1以及实施方式2的系统辨识装置中的平行投影的奇异值与维数的关系的概略图。

图4是示出实施方式1的系统辨识装置中的系统维数决定部的内部结构的框图。

图5是示出在实施方式1以及实施方式2的系统辨识装置中，辨识出的线性离散时间系统在时域或者频域中的误差平方和的范数(norm)与维数的关系的概略图。

图6是示出对实施方式2的系统辨识装置中的动态系统进行了M序列加振的情况下的系统输入输出的时间波形的概略图。

图7是示出实施方式2的系统辨识装置中的系统维数决定部的内部结构的框图。

图8是示出实施方式3的整体结构的框图。

(符号说明)

1：系统输入输出提取部；2：块汉克尔矩阵生成部；3：输入输出矢量生成部；4：LQ分解部；5：平行投影生成部；6：奇异值分解部；7：系统维数决定部；8：状态矢量生成部；9：系统矩阵辨识部；10：系统辨识装置；11：系统输入；12：系统输出；13：系统输入阈值；21：(理想系统输入输出中的平行投影的)奇异值分布；22：(现实的系统输入输出中的平行投影的)奇异值分布；31：递归系统矩阵推测部；32：系统特性推测部；33：系统维数推测部；34：系统稳定性评价部；41：(时域或者频域中的)误差平方和的范数分布；42：(时域或者频域中的)误差平方和范数阈值；51：DC伺服马达。

具体实施方式

以下，参照附图，说明本发明的实施方式的系统辨识装置。另外，本发明不限于以下所示的实施方式。

实施方式1.

图1是示出实施方式1的系统辨识装置的整体结构的框图，图2是示出实施方式1的系统辨识装置中的系统输入输出的时间波形的概略图。

在实施方式1的系统辨识装置10中，如图1以及图2所示，将对作为辨识对象的动态系统施加了伪随机输入的情况下的系统输入11(u(jT_S)(j＝0，1，2，…))以及系统输出12(y(jT_S)(j＝0，1，2，…))作为输入。

系统输入输出提取部1针对用预先设定的比例阈值和系统输入11的最大值相乘而得到的值决定的系统输入阈值13，将系统输入11的绝对值在系统输入阈值13以上的时刻的最小值作为伪随机输入施加时刻(在图2中为3T_S)，将伪随机输入施加时刻以后的系统输入11以及系统输出12分别作为辨识用输入数据(u_id(jT_S)(j＝0，1，2，…))以及辨识用输出数据(y_id(jT_S)(j＝0，1，2，…))提取并输出。

块汉克尔矩阵生成部2根据从系统输入输出提取部1输出的辨识用输入数据u_id(jT_S)(j＝0，1，2，…)以及辨识用输出数据y_id(jT_S)(j＝0，1，2，…)，生成块汉克尔矩阵U_p、U_f以及Y_p、Y_f。

输入输出矢量生成部3根据块汉克尔矩阵U_p、U_f、Y_p、Y_f，生成动态系统的输入矢量^～U_K|K以及输出矢量^～Y_K|K。

LQ分解部4生成结合块汉克尔矩阵U_p、U_f、Y_p、Y_f而得到的数据矩阵，生成并输出对该数据矩阵进行了LQ分解的子矩阵L₂₂、L₃₂。

平行投影生成部5根据从LQ分解部4输出的子矩阵L₂₂、L₃₂、和从块汉克尔矩阵生成部2输出的块汉克尔矩阵U_p、Y_p，生成动态系统的平行投影Θ。

奇异值分解部6对从平行投影生成部5输出的平行投影Θ进行奇异值分解，输出将平行投影Θ的左奇异矢量作为列矢量的第1正交矩阵U、将平行投影Θ的右奇异矢量作为列矢量的第2正交矩阵V以及平行投影Θ的奇异值σ_i(i＝1，2，3…)。

系统维数决定部7根据从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)、从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K、和作业者指定的系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)，针对属于该搜索范围的各个维数n_i(i＝1，2，…，a)，辨识记述动态系统的线性离散时间系统的系统矩阵。进而，根据该系统矩阵，计算在对与属于该搜索范围的各个维数n_i(i＝1，2，…，a)对应的线性离散时间系统应用了实际的辨识用输入数据u_id(jT_S)(j＝0，1，2，…)的情况下的系统输出，根据与动态系统的实际的辨识用输出数据y_id(jT_S)(j＝0，1，2，…)(在图1中记载为动态系统的系统特性)的比较，决定系统维数n。

在状态矢量生成部8中，根据从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)、和从系统维数决定部7输出的系统维数n，生成动态系统的状态矢量^～X_K+1、^～X_K。

系统矩阵辨识部9根据从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K、和从状态矢量生成部8输出的动态系统的状态矢量^～X_K+1、^～X_K，辨识并输出记述动态系统的线性离散时间系统的系统矩阵A_d、B_d、C_d、D_d。

图3是示出实施方式1的系统辨识装置10中的平行投影Θ的奇异值σ_i与维数(i＝1，2，3…)的关系的概略图，图4是示出实施方式1的系统辨识装置10中的系统维数决定部7的内部结构的框图，图5是示出实施方式1的系统辨识装置10中的辨识出的线性离散时间系统的系统输出和动态系统的实际的系统输出在时域中的误差平方和的范数||e_n||与维数n_i(i＝1，2，…，a)的关系的概略图。

如图3所示，根据动态系统的系统输入输出计算出的平行投影Θ的奇异值σ_i(i＝1，2，3…)理想地相对维数(i＝1，2，3…)成为例如奇异值分布21所示的关系。在该情况下，能够明确地规定具有有效值的奇异值的个数，该个数与动态系统的系统维数n对应(在图3的情况下系统维数n＝4)。

另一方面，根据受到观测噪声等的影响的现实的系统输入输出计算出的奇异值σ_i相对维数(i＝1，2，3…)成为例如奇异值分布22所示的关系，所以具有有效值的奇异值和成为可忽略的微小的值的奇异值的边界变得不明确，未必始终决定出最佳的系统维数n。因此，产生关于系统维数n的决定需要反复试验这样的问题。

因此，在实施方式1的系统辨识装置10中，通过系统维数决定部7执行图4所示的处理。具体而言，如以下所述。

在系统维数决定部7中设置递归系统矩阵推测部31、系统特性推测部32以及系统维数推测部33。

递归系统矩阵推测部31关于与属于作业者预先指定的系统维数的搜索范围n_i(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的第1维数n_i对应的系统矩阵的辨识，使用与比第1维数n_i低1个等级的第2维数n_i-1对应的系统矩阵A_d，n_i-1、B_d，n_i-1、C_d，n_i-1、D_d，n_i-1的辨识结果、从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)中分别比第2维数n_i-1大且在第1维数n_i以下的右奇异矢量v_j及奇异值σ_j(j＝n_i-1+1，n_i-1+2，…，n_i)、和从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K，通过递归方法，辨识与第1维数n_i对应的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i。

接下来，系统特性推测部32针对属于系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的各个维数，根据从递归系统矩阵推测部31输出的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i，计算对辨识出的线性离散时间系统应用了实际的辨识用输入数据u_id(jT_S)(j＝0，1，2，…)的情况下的系统输出。

另外，使i递增并直至i＝a为止执行递归系统矩阵推测部31以及系统特性推测部32的处理。

系统维数推测部33构成为：计算从系统特性推测部32输出的线性离散时间系统的系统输出、和动态系统的实际的辨识用输出数据y_id(jT_S)(j＝0，1，2，…)(在图4中记载为动态系统的系统特性)在时域中的误差平方和e_ni(i＝1，2，…，a)，如图5所示，将误差平方和的范数||e_ni||的分布41在预先设定的误差平方和范数阈值42以下的维数中的最小的维数决定为系统维数n并输出(在图5的情况下系统维数n＝n₆)。

接下来，说明实施方式1的系统辨识装置的动作。

能够如下式那样，用1输入P输出的n维线性离散时间系统来记述作为辨识对象的动态系统。

[式1]

x((j+1)T_s)＝A_dx(jT_s)+B_du(jT_s)

y(jT_s)＝C_dx(jT_s)+D_du(jT_s)

其中，状态矢量：x∈Rⁿ

系统输入：u∈R

系统输出：y∈R^P

系统矩阵：A_d∈R^n×n，B_d∈Rⁿ，C_d∈R^P×n，D_d∈R^P

如果用伪随机输入构成向上述动态系统的系统输入u(jT_S)，则该系统输入u(jT_S)以及与上述[式1]式对应的系统输出y(jT_S)为例如图2所示的系统输入11以及系统输出12那样的时间波形。

此处，如在图1以及图2的说明中叙述的那样，将预先设定的比例阈值和系统输入11(u(jT_S))的最大值相乘而得到的下式用作系统输入阈值13。

[式2]

系统输入比例阈值·max(u(jT_s))

系统输入输出提取部1将系统输入11的绝对值在系统输入阈值13以上的时刻的最小值确定为伪随机输入施加时刻j_minT_S(在图2的情况下j_minT_S＝3T_S)。

另外，系统输入输出提取部1使用下式来提取伪随机输入施加时刻j_minT_S以后的系统输入11以及系统输出12。

[式3]

u_id(jT_s)＝u((j_min+j)T_s)(j＝0，1，2，…)

y_id(jT_s)＝y((j_min+j)T_s)(j＝0，1，2，…)

进而，系统输入输出提取部1将使用上述[式3]提取出的各个值作为辨识用输入数据u_id(jT_S)以及辨识用输出数据y_id(jT_S)，从而从作为对象的动态系统的系统输入输出去除伪随机输入施加之前的系统静止时域数据。

块汉克尔矩阵生成部2根据从系统输入输出提取部1输出的辨识用输入数据u_id(jT_S)(j＝0，1，2，…)以及辨识用输出数据y_id(jT_S)(j＝0，1，2，…)，生成用下式给出的块汉克尔矩阵U_p、U_f、Y_p、Y_f。

[式4]

输入输出矢量生成部3根据该块汉克尔矩阵U_p、U_f、Y_p、Y_f，生成用下式给出的动态系统的输入矢量^～U_K|K以及输出矢量^～Y_K|K。

[式5]

${\stackrel{\&OverBar;}{U}}_{K|K}=\&lsqb;\begin{array}{cccc}u\left({\mathrm{KT}}_s\right)& u\left(\right(K+1\left)T_s\right)& ...& u\left(\right(K+N-2)T_s)\end{array}\&rsqb;=U_f(1,1:N-1)\&Element;R^{1\&times;(N-1)}$

${\stackrel{\&OverBar;}{Y}}_{K|K}=\&lsqb;\begin{array}{cccc}y\left({\mathrm{KT}}_s\right)& y\left(\right(K+1\left)T_s\right)& ...& y\left(\right(K+N-2)T_s)\end{array}\&rsqb;=Y_f(1:P,1:N-1)\&Element;R^{P\&times;(N-1)}$

LQ分解部4生成用结合了块汉克尔矩阵U_p、U_f、Y_p、Y_f的下式给出的数据矩阵。

[式6]

$\left[\begin{array}{c}{U}_f\\ U_p\\ Y_p\\ Y_f\end{array}\right]$

另外，LQ分解部4对上述数据矩阵如下式那样进行LQ分解，并且从进行了LQ分解的矩阵的元素输出子矩阵L₂₂、L₃₂。

[式7]

$\left[\begin{array}{c}{U}_f\\ U_p\\ Y_p\\ Y_f\end{array}\right]=\left[\begin{array}{ccc}{L}_{11}& 0& 0\\ L_{21}& L_{22}& 0\\ L_{31}& L_{32}& L_{33}\end{array}\right]\left[\begin{array}{c}{Q}_1^T\\ Q_2^T\\ Q_3^T\end{array}\right]$

其中，正交矩阵：Q₁∈R^N×K，Q₂∈R^N×K(1+P)，Q₃∈R^N×KP

块下三角矩阵：L₁₁∈R^K×K，L₂₂∈R^{K(1+P)×K(1+P)}，L₃₃∈R^KP×KP

L₂₁∈R^K(1+P)×K，L₃₁∈R^KP×K，L₃₂∈R^KP×K(1+P)

平行投影生成部5根据从LQ分解部4输出的子矩阵L₂₂、L₃₂、和从块汉克尔矩阵生成部2输出的块汉克尔矩阵U_p、Y_p，生成用下式定义的动态系统的平行投影Θ。

[式8]

奇异值分解部6通过对用上式表示的平行投影Θ进行奇异值分解，输出用下式给出的将平行投影Θ的左奇异矢量u_j作为列矢量的第1正交矩阵U、将平行投影Θ的右奇异矢量v_j作为列矢量的第2正交矩阵V、以及平行投影Θ的奇异值σ_i(i＝1，2，3...)。

[式9]

Θ＝U∑V^r

其中，第1正交矩阵：U＝[u₁ u₂ … u_KP]∈R^KP×KP

第2正交矩阵：V＝[v₁ v₂ … v_N]∈R^N×N

平行投影Θ的奇异值σ₁≥σ₂≥…≥σ_n≥σ_n+1≥σ_n-2≥…

作为对象的动态系统的系统维数n能够根据以下的关系来决定，以下的关系是：在平行投影Θ的奇异值中，具有有效值的奇异值是n个，相比于这n个奇异值，在第n+1个以后充分小。

[式10]

如图3所示，根据动态系统的系统输入输出计算出的平行投影Θ的奇异值σ_i理想地相对维数(i＝1，2，3…)成为例如奇异值分布21所示的关系。在该情况下，能够明确地规定具有有效值的奇异值的个数，能够根据该个数决定动态系统的系统维数n(在图3的情况下系统维数n＝4)。另一方面，根据受到观测噪声等的影响的现实的系统输入输出计算出的奇异值σ_i相对维数(i＝1，2，3…)成为例如奇异值分布22所示的关系，所以具有有效值的奇异值和成为可忽略的微小的值的奇异值的边界σ_n>>σ_n+1变得不明确。因此，在以往的方法中，存在如下课题：未必始终决定出最佳的系统维数n，为了决定出最佳的系统维数n而需要反复试验。

因此，在实施方式1的系统辨识装置10中，当“在时域中最适合于实际的系统输入输出”这样的前提下，在系统维数决定部7中决定最佳的系统维数n。在系统维数决定部7中，如图1所示，根据从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)、从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K、和作业者指定的系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)，针对属于该搜索范围的各个维数n_i(i＝1，2，…，a)，辨识记述动态系统的线性离散时间系统的系统矩阵。进而，根据该系统矩阵，计算对与属于该搜索范围的各个维数n_i(i＝1，2，…，a)对应的线性离散时间系统应用了实际的辨识用输入数据u_id(jT_S)(j＝0，1，2，…)的情况下的系统输出，根据与动态系统的实际的辨识用输出数据y_id(jT_S)(j＝0，1，2，…)(在图1中记载为动态系统的系统特性)的比较，决定系统维数n。

具体而言，如图4所示，递归系统矩阵推测部31关于与属于作业者指定的系统维数的搜索范围n_i(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的第1维数n_i对应的系统矩阵的辨识，使用与比第1维数n_i低1个等级的第2维数n_i-1对应的系统矩阵A_d，n_i-1、B_d，n_i-1、C_d，n_i-1、D_d，n_i-1的辨识结果、从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)中分别比第2维数n_i-1大且在第1维数n_i以下的右奇异矢量v_j及奇异值σ_j(j＝n_i-1+1，n_i-1+2，…，n_i)、和从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K，通过下式所示的递归方法，辨识与第1维数n_i对应的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i。

[式11]

与第1维数n_i对应的状态矢量：

与第1维数n_i对应的系统矩阵

其中，

系统特性推测部32针对属于系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的各个维数n_i，根据从递归系统矩阵推测部31输出的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i，计算对辨识出的线性离散时间系统应用了实际的辨识用输入数据u_id(jT_S)(j＝0，1，2，…)(参照[式3])的情况下的系统输出^y_id，n_i(jT_S)(j＝0，1，2，…)。另外，“^y”的标记表示在“y”这样的字符的上部附加了“^”的记号的代替标记。

另外，系统维数推测部33使用下式来计算从系统特性推测部32输出的线性离散时间系统的系统输出^y_id，n_i(jT_S)(j＝0，1，2，…)、和动态系统的实际的辨识用输出数据y_id(jT_S)(j＝0，1，2，…)(在图4中记载为动态系统的系统特性)在时域中的误差平方和en_i(i＝1，2，…，a)。

[式12]

使由上式表示的误差平方和的范数||en_i||成为最小的维数n_i是“在时域中最适合于实际的系统输入输出”的系统维数n。另一方面，如果观测噪声是白噪声，则无论其噪声等级如何，伴随维数n_i的增加，实际的范数||en_i||都单调减少，如图5所示，在某个维数以上成为大致恒定值。因此，此处，为了避免系统维数n的推测值超出必要地成为高维数，规定用下式给出的误差平方和范数阈值42。

[式13]

误差平方和容许值·

系统维数推测部33将误差平方和的范数||en_i||的分布41在上述误差平方和范数阈值42以下的维数中的最小的维数决定为系统维数n而输出(在图5的例子中，系统维数n＝n₆)。

状态矢量生成部8根据从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)、和从系统维数决定部7输出的系统维数n，依照下式生成动态系统的状态矢量^～X_K、^～X_K+1。

[式14]

$\begin{array}{c}{X}_f=\&lsqb;\begin{array}{cccc}x\left({\mathrm{KT}}_s\right)& x\left(\right(K+1\left)T_s\right)& ...& x\left(\right(K+N-1\left)T_s\right)\end{array}\&rsqb;\\ \&ap;{\&Sigma;}_n^{1/2}{V}_n^T=\&Sigma;{(1:n,1:n)}^{1/2}V{(:,1:n)}^T\&Element;R^{n\&times;N}\end{array}$

${\stackrel{\&OverBar;}{X}}_{K+1}=\&lsqb;\begin{array}{cccc}x\left(\right(K+1\left)T_s\right)& x\left(\right(K+2\left)T_s\right)& ...& x\left(\right(K+N-1\left)T_s\right)\end{array}\&rsqb;=X_f(:,2:N)\&Element;R^{n\&times;(N-1)}$

${\stackrel{\&OverBar;}{X}}_K=\&lsqb;\begin{array}{cccc}x\left({\mathrm{KT}}_s\right)& x\left(\right(K+1\left)T_s\right)& ...& x\left(\right(K+N-2\left)T_s\right)\end{array}\&rsqb;=X_f(:,1:N-1)\&Element;R^{n\&times;(N-1)}$

最后，系统矩阵辨识部9根据从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K、和从状态矢量生成部8输出的动态系统的状态矢量^～X_K、^～X_K+1，使用下式来辨识记述动态系统的线性离散时间系统的系统矩阵A_d、B_d、C_d、D_d并输出。

[式15]

$\left[\begin{array}{cc}{A}_d& B_d\\ C_d& D_d\end{array}\right]=\left(\left[\begin{array}{c}{\stackrel{\&OverBar;}{X}}_{K+1}\\ {\stackrel{\&OverBar;}{Y}}_{K|K}\end{array}\right]{\left[\begin{array}{c}{\stackrel{\&OverBar;}{X}}_K\\ {\stackrel{\&OverBar;}{U}}_{K|K}\end{array}\right]}^T\right){\left(\left[\begin{array}{c}{\stackrel{\&OverBar;}{X}}_K\\ {\stackrel{\&OverBar;}{U}}_{K|K}\end{array}\right]{\left[\begin{array}{c}{\stackrel{\&OverBar;}{X}}_K\\ {\stackrel{\&OverBar;}{U}}_{K|K}\end{array}\right]}^T\right)}^{-1}\&Element;R^{(P+n)\&times;(1-n)}$

其中A_d∈R^n×n，B_d∈Rⁿ，C_d∈R^P×n，D_d∈R^P

这样，根据实施方式1的系统辨识装置10，即使在根据现实的系统输入输出计算出的平行投影Θ的奇异值σ_i(i＝1，2，3…)为平滑的单调减少，从而具有有效值的奇异值和成为在辨识中可忽略的微小的值的奇异值的边界变得不明确的情况下，也能够从系统维数n的决定来排除反复试验，决定相对现实的动态系统在时域中一致度高的系统维数n，能够进行记述动态系统的线性离散时间系统的辨识。

另外，还能够通过从动态系统的现实的系统输入输出去除伪随机输入施加之前的系统静止时域数据来提高辨识精度。

进而，由于递归系统矩阵推测部31的存在，能够降低用于决定相对现实的动态系统一致度高的系统维数n的运算量。

另外，在实施方式1的系统辨识装置10中，计算对线性离散时间系统应用了实际的辨识用输入数据的情况下的系统输出而作为系统特性，将该系统输出和动态系统的实际的辨识用输出数据在时域中的误差平方和的范数分布41在阈值42以下的维数中的最小的维数决定为系统维数n。但是，本发明不限于此，也可以计算线性离散时间系统的系统特性而作为频率响应，根据该频率响应和从动态系统的辨识用输入输出数据得到的实际的频率响应在频域中的误差平方和来决定系统维数n。在该情况下，也可以进而根据动态系统的实际的频率响应来决定权重函数，根据对线性离散时间系统的频率响应、和动态系统的实际的频率响应在频域中的误差平方值乘以该权重函数而得到的值的相加值，来决定系统维数n。

实施方式2.

接下来，说明实施方式2的系统辨识装置。另外，示出实施方式2的系统辨识装置的整体结构的框图、示出平行投影Θ的奇异值σ_i与维数(i＝1，2，3…)的关系的概略图、以及示出辨识出的线性离散时间系统的频率响应和动态系统的实际的频率响应在频域中的误差平方和的范数||en_i||与维数n_i(i＝1，2，…，a)的关系的概略图分别与在实施方式1的说明中使用的图1、图3以及图5相同。

图6是示出对实施方式2的系统辨识装置中的动态系统进行了M序列加振的情况下的系统输入输出的时间波形的概略图。

如图6所示，在实施方式2的系统辨识装置10中，根据对作为辨识对象的动态系统施加了M序列信号的情况下的系统输入11(u(jT_S)(j＝0，1，2，…))以及系统输出12(y(jT_S)(j＝0，1，2，…))，辨识记述该系统的线性离散时间系统。

图7是示出实施方式2的系统辨识装置中的系统维数决定部7的内部结构的框图。在图7中，附加与图4相同的符号的构成要素是与实施方式1相同或者等同的构成要素，并追加系统稳定性评价部34。

在实施方式2的系统维数决定部7中，如图7所示，针对属于作业者指定的系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的各个维数n_i，根据由递归系统矩阵推测部31辨识出的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i，由系统稳定性评价部34评价线性离散时间系统的稳定性。

系统特性推测部32针对由系统稳定性评价部34判断为稳定系统的维数，根据从递归系统矩阵推测部31输出的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i，计算与辨识出的线性离散时间系统有关的频率响应。

系统维数推测部33根据从动态系统的系统输入输出得到的实际的频率响应(在图7中记载为动态系统的系统特性)来决定权重函数，计算对从系统特性推测部32输出的线性离散时间系统的频率响应、和动态系统的实际的频率响应在频域中的误差平方值乘以该权重函数而得到的值的相加值en_i(n_i：成为稳定系统的维数)，将如图5所示该相加值的范数||en_i||的分布41在预先设定的误差平方和范数阈值42以下的维数中的最小的维数决定为系统维数n而输出(在图5的情况下系统维数n＝n₆)。

接下来，说明实施方式2的系统辨识装置的动作。

设为作为辨识对象的动态系统能够通过[式1]记述为1输入P输出的n维线性离散时间系统。如果用M序列信号构成向该动态系统的系统输入u(jT_S)，则该系统输入u(jT_S)以及与[式1]对应的系统输出y(jT_S)成为例如图6所示的系统输入11以及系统输出12那样的时间波形。

在实施方式2的系统辨识装置10中，如图1以及图6所示，系统输入输出提取部1将预先设定的比例阈值和系统输入11(u(jT_S))的最大值相乘而得到的[式2]作为系统输入阈值13，将系统输入11的绝对值在系统输入阈值13以上的时刻的最小值确定为M序列信号施加时刻j_minT_s(图6的例子是j_minT_s＝2T_s)。

另外，系统输入输出提取部1使用[式3]提取M序列信号施加时刻j_minT_s以后的系统输入11以及系统输出12，将提取出的各个作为辨识用输入数据u_id(jT_S)以及辨识用输出数据y_id(jT_S)，从而从作为对象的动态系统的系统输入输出去除M序列信号施加之前的系统静止时域数据。

接下来，与实施方式1同样地，块汉克尔矩阵生成部2生成用[式4]给出的块汉克尔矩阵U_p、U_f、Y_p、Y_f，输入输出矢量生成部3生成用[式5]给出的动态系统的输入矢量^～U_K|K以及输出矢量^～Y_K|K，LQ分解部4对结合块汉克尔矩阵U_p、U_f、Y_p、Y_f而得到的数据矩阵([式6])进行LQ分解([式7])，输出子矩阵L₂₂、L₃₂。

平行投影生成部5生成使用[式8]定义的动态系统的平行投影Θ，奇异值分解部6对所生成的平行投影Θ进行奇异值分解，从而输出用[式9]给出的第1正交矩阵U、第2正交矩阵V以及奇异值σ_i(i＝1，2，3…)。

在系统维数决定部7中，执行图7所示的处理。首先，递归系统矩阵推测部31关于属于作业者指定的系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的各个维数n_i，通过[式11]所示的递归方法，辨识对应的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i。

接下来，系统稳定性评价部34针对属于作业者指定的系统维数的搜索范围n_i＝(n₁，n₂，…，n_a)(其中n₁<n₂<…<n_a)的各个维数n_i，根据递归系统矩阵推测部31辨识出的系统矩阵A_d，n_i，通过以下的内容来评价线性离散时间系统的稳定性。

[式16]

维数n_i的线性离散时间系统稳定

系统矩阵A_d，n_i的所有固有值的绝对值小于1

系统矩阵A_d，n_i的所有固有值存在于单位圆内

系统特性推测部32针对由系统稳定性评价部34判断为稳定系统的维数，根据递归系统矩阵推测部31生成的系统矩阵A_d，n_i、B_d，n_i、C_d，n_i、D_d，n_i，计算辨识出的线性离散时间系统的频率响应^Hn_i(kΔf)(k＝0，1，2，…，N/2-1)。

在实施方式2的系统辨识装置10中，当“在频域中最适合于实际的频率响应”这样的前提下，通过系统维数推测部33决定最佳的系统维数n。具体而言，如以下所述。

首先，根据用下式给出的辨识用输入输出数据u_id(jT_S)、y_id(jT_S)的有限离散化傅立叶变换U_id(kΔf)、Y_id(kΔf)(k＝0，1，2，…，N/2-1)，计算使用以下各式得到的动态系统的实际的频率响应H(kΔf)(k＝0，1，2，…，N/2-1)(在图7中记载为动态系统的系统特性)。

[式17]

$U_{id}\left(k\mathrm{\&Delta;}f\right)=\frac{T}{N}\underset{j=0}{\overset{N-1}{\&Sigma;}}{u}_{id}\left({\mathrm{jT}}_s\right)\exp \&lsqb;-i\&CenterDot;2\mathrm{\&pi;}\frac{jk}{N}\&rsqb;,(k=0,1,2,...,\frac{N}{2}-1)$

$Y_{id}\left(k\mathrm{\&Delta;}f\right)=\frac{T}{N}\underset{j=0}{\overset{N-1}{\&Sigma;}}{y}_{id}\left({\mathrm{jT}}_s\right)\exp \&lsqb;-i\&CenterDot;2\mathrm{\&pi;}\frac{jk}{N}\&rsqb;,\left(k=0,1,2,...,\frac{N}{2}-1\right)$

其中，采样周期：

采样频率：

频率分辨率：

时刻：

频率：

[式18]

$H\left(k\mathrm{\&Delta;}f\right)=\frac{Y_{id}\left(k\mathrm{\&Delta;}f\right)U_{id}{\left(k\mathrm{\&Delta;}f\right)}^{*}}{U_{id}\left(k\mathrm{\&Delta;}f\right)U_{id}{\left(k\mathrm{\&Delta;}f\right)}^{*}},(k=0,1,2,...,\frac{N}{2}-1)$

接下来，根据频率响应H(kΔf)(k＝0，1，2，…，N/2-1)，决定例如由对高增益且低频区域附加权重而得到的下式表示那样的权重函数W(kΔf)(k＝0，1，2，…，N/2-1)。

[式19]

$W\left(k\mathrm{\&Delta;}f\right)=\frac{\left|H\left(k\mathrm{\&Delta;}f\right)\right|}{k\mathrm{\&Delta;}f},(k=0,1,2,...,\frac{N}{2}-1)$

然后，通过下式，计算对从系统特性推测部32输出的线性离散时间系统的频率响应^Hn_i(kΔf)和动态系统的实际的频率响应H(kΔf)在频域中的误差平方值乘以权重函数W(kΔf)而得到的值的相加值en_i(n_i：成为稳定系统的维数)。

[式20]

使该加权误差平方和的范数||en_i||成为最小的维数n_i为“根据权重函数，在频域中最适合于实际的频率响应的”稳定的系统维数n。此处，将如图5所示加权误差平方和的范数||en_i||的分布41在用[式13]给出的误差平方和范数阈值42以下的维数中的最小的维数决定为系统维数n而输出(在图5的例子中系统维数n＝n₆)。

状态矢量生成部8根据从奇异值分解部6输出的第2正交矩阵V及奇异值σ_i(i＝1，2，3…)、和从系统维数决定部7输出的系统维数n，使用[式14]来生成动态系统的状态矢量^～X_K、^～X_K+1。

最后，系统矩阵辨识部9根据从输入输出矢量生成部3输出的动态系统的输入矢量^～U_K|K及输出矢量^～Y_K|K、和从状态矢量生成部8输出的动态系统的状态矢量^～X_K、^～X_K+1，使用[式15]来辨识并输出记述动态系统的线性离散时间系统的系统矩阵A_d、B_d、C_d、D_d。

这样，根据实施方式2的系统辨识装置10，即使在根据现实的系统输入输出计算出的平行投影Θ的奇异值σ_i(i＝1，2，3…)为平滑的单调减少，从而具有有效值的奇异值和成为在辨识中可忽略的微小的值的奇异值的边界变得不明确的情况下，也能够从系统维数n的决定来排除反复试验，根据频域中的权重函数决定相对现实的动态系统一致度高的系统维数n，能够进行记述动态系统的线性离散时间系统的辨识。

另外，还能够通过从动态系统的现实的系统输入输出去除M序列信号施加之前的系统静止时域数据，提高辨识精度。

进而，由于递归系统矩阵推测部31的存在，能够降低用于决定相对现实的动态系统一致度高的系统维数n的运算量。

另外，由于系统稳定性评价部34的存在，在现实的动态系统是稳定系统的情形明确的情况下，能够进行限定于稳定系统的线性离散时间系统的辨识。

另外，在实施方式2的系统辨识装置10中，计算线性离散时间系统的系统特性而作为频率响应，将该频率响应和从动态系统的辨识用输入输出数据得到的实际的频率响应在频域中的误差平方和的范数分布41在预先设定的阈值42以下的维数中的最小的维数决定为系统维数n。但是，本发明不限于此，也可以计算对线性离散时间系统应用了实际的辨识用输入数据的情况下的系统输出而作为系统特性，根据该系统输出和动态系统的实际的辨识用输出数据在时域中的误差平方和，决定系统维数n。

实施方式3.

在实施方式3中，说明作为辨识对象的动态系统是DC伺服马达的情况。图8是示出实施方式3的整体结构的框图。在本实施方式中，图8所示的系统辨识装置10为与图1所示的实施方式1的系统辨识装置10相同或者等同的结构。在本实施方式中，作为DC伺服马达51的输入电流[A_rms]，输入例如M序列信号等伪随机信号，将该伪随机信号作为动态系统的系统输入11(u(jT_S)(j＝0，1，2，…))。另外，取得角速度[rad/s]作为动态系统的系统输出12(y(jT_S)(j＝0，1，2，…))。系统辨识装置10将这些系统输入输出以及系统维数的搜索范围作为输入，辨识记述DC伺服马达51的线性离散时间系统。此时，关于系统维数的搜索范围，如n_i＝(1，2，…，50)那样，相对预测的系统维数取充分的宽度来设定即可。在系统辨识装置10中，能够进行相对现实的动态系统一致度高的系统维数的决定、和记述动态系统的线性离散时间系统的辨识，所以能够将该线性离散时间系统用于伺服马达的控制系中的参数设计、滤波器的参数设计等。

Claims (7)

1.一种系统辨识装置，将对作为辨识对象的动态系统施加了伪随机输入的情况下的系统输入输出以及被指定的系统维数的搜索范围作为输入，该系统辨识装置的特征在于，具备：

系统输入输出提取部，从所述动态系统的系统输入输出提取应用于辨识的辨识用输入输出数据；

块汉克尔矩阵生成部，根据所述辨识用输入输出数据，生成块汉克尔矩阵；

输入输出矢量生成部，根据所述块汉克尔矩阵，生成所述动态系统的输入矢量以及输出矢量；

LQ分解部，结合所述块汉克尔矩阵而生成数据矩阵，输出对该数据矩阵进行LQ分解而得到的子矩阵；

平行投影生成部，根据所述子矩阵和所述块汉克尔矩阵，生成平行投影；

奇异值分解部，通过所述平行投影的奇异值分解，输出将所述平行投影的左奇异矢量作为列矢量的第1正交矩阵、将该平行投影的右奇异矢量作为列矢量的第2正交矩阵以及该平行投影的奇异值；

系统维数决定部，根据所述第2正交矩阵及所述奇异值、所述动态系统的输入矢量及输出矢量、和所述搜索范围，针对属于该搜索范围的各个维数，辨识记述动态系统的线性离散时间系统的系统矩阵，进而根据基于该系统矩阵计算出的线性离散时间系统的系统特性和动态系统的实际的系统特性的比较，决定系统维数；

状态矢量生成部，根据所述第2正交矩阵及奇异值、和所决定的所述系统维数，生成所述动态系统的状态矢量；和

系统矩阵辨识部，根据所述动态系统的输入矢量及输出矢量以及所述动态系统的状态矢量，辨识记述该动态系统的线性离散时间系统的系统矩阵，

所述系统辨识装置将辨识出的所述系统矩阵作为记述所述动态系统的线性离散时间系统输出。

2.根据权利要求1所述的系统辨识装置，其特征在于，

所述系统维数决定部具备：

系统特性推测部，针对属于所述搜索范围的各个维数，计算关于辨识出的线性离散时间系统应用了实际的辨识用输入数据的情况下的系统输出，将该系统输出作为线性离散时间系统的系统特性输出；和

系统维数推测部，将所述线性离散时间系统的系统输出和所述动态系统的实际的辨识用输出数据在时域中的误差平方和的范数在设定阈值以下的维数中的最小的维数决定为系统维数而输出。

3.根据权利要求1所述的系统辨识装置，其特征在于，

所述系统维数决定部具备：

系统特性推测部，针对属于所述搜索范围的各个维数，计算辨识出的线性离散时间系统的频率响应，将该频率响应作为线性离散时间系统的系统特性输出；和

系统维数推测部，将所述线性离散时间系统的频率响应和从所述动态系统的系统输入输出得到的实际的频率响应在频域中的误差平方和的范数在设定阈值以下的维数中的最小的维数决定为系统维数而输出。

4.根据权利要求3所述的系统辨识装置，其特征在于，

所述系统维数推测部根据从所述动态系统的系统输入输出得到的实际的频率响应，决定权重函数，

所述系统维数推测部计算对所述线性离散时间系统的频率响应和所述动态系统的实际的频率响应在频域中的误差平方值乘以该权重函数而得到的值的相加值，

所述系统维数推测部将该相加值的范数在设定阈值以下的维数中的最小的维数决定为系统维数而输出。

5.根据权利要求1所述的系统辨识装置，其特征在于，

所述系统维数决定部具备递归系统矩阵推测部，该递归系统矩阵推测部关于与属于所述搜索范围的第1维数对应的系统矩阵的辨识，使用与在该搜索范围中比第1维数低1个等级的第2维数对应的系统矩阵的辨识结果、所述第2正交矩阵及奇异值中分别比所述第2维数大且在所述第1维数以下的右奇异矢量及奇异值、以及所述动态系统的输入矢量及输出矢量，通过递归方法来辨识与所述第1维数对应的系统矩阵。

6.根据权利要求1所述的系统辨识装置，其特征在于，

所述系统输入输出提取部将所设定的比例阈值和系统输入的最大值相乘而得到的值作为系统输入阈值，将系统输入的绝对值在系统输入阈值以上的时刻的最小值作为伪随机输入施加时刻，将伪随机输入施加时刻以后的系统输入以及系统输出分别作为辨识用输入数据以及辨识用输出数据提取。

7.根据权利要求1所述的系统辨识装置，其特征在于，

所述系统维数决定部具备系统稳定性评价部，该系统稳定性评价部针对属于所述搜索范围的各个维数来评价线性离散时间系统的稳定性，

根据与成为稳定系统的维数对应的线性离散时间系统的系统特性，决定系统维数。