JP4470855B2 - Model identification device - Google Patents

Description

  The present invention relates to a model identification apparatus for identifying a mathematical model from a frequency response function indicating a frequency characteristic of a dynamic system.

The problem of deriving the transfer function G (z) indicating the characteristic from the frequency response function H (z) obtained as a priori information by the measurement is as follows: By defining an unknown variable of the polynomial as an unknown variable and defining a simultaneous equation with the frequency characteristic indicated by the frequency response function as a constant, it can be reduced to a problem of solving the linear simultaneous equation.

If there are more equations than unknown variables, the solution is not uniquely determined, and the square of the difference at each frequency point between the frequency characteristic of the above transfer function and the frequency characteristic indicated by the above frequency response function. Is a least-squares problem (Equation 2) for finding a solution of a linear simultaneous equation that minimizes.

Further, in order to eliminate the solution bias (linearization error) due to linearization, which is a problem in Equation 2, for example, weights using the solution are added as shown in Equation 3, and the iterative calculation is performed while updating the weights. Find the convergence solution. Let this converged solution be the coefficient of the transfer function.

In the conventional least squares problem in Patent Document 1 described above, the coefficient term matrix A constituting the simultaneous equations expressed by Equation 1 has no error and is composed of frequency characteristic data at each frequency point obtained by measurement. The least square solution is obtained on the assumption that an error exists only in the column vector d.
However, as shown in Equations 1 and 2, the matrix A includes the frequency characteristics obtained by the measurement due to the linearization, and is obtained on the assumption that the matrix A has no error. The obtained least squares solution has a drawback that it includes a bias error due to a measurement error.

In the iterative calculation based on Equation 3, since the weight is set based on the solution including the bias error, an error other than the measurement error is included in the matrix A, and the iterative calculation does not converge, Alternatively, there is a disadvantage that the linearization error cannot be eliminated even if a convergent solution is obtained.
In order to compensate for these drawbacks, a total least square method (Total) is a method that can obtain a least squares solution of simultaneous equations in consideration of the fact that errors are included not only in the column vector d but also in the coefficient term matrix A in recent years. A method called “Last Square” has been proposed (see Non-Patent Document 2).

Satoshi Noda and Naoto Nagaoka, "Development of ARMA Transient Phenomenon Calculation Model by Linearized Least Squares", IEEJ Transactions B, 114, 1994, P-396-401 Sabin Van Huffel, Joos Vandwalle, “The Total Least Squares Programmable Aspects and Analysis”, Society for Industrial 91 and Applied.

  Applying the comprehensive least square method in the above-mentioned patent document to Equations 2 and 3 can be easily performed by performing singular value decomposition on a matrix composed of coefficient term matrix A and column vector b. However, the overall least squares method is more susceptible to excessive errors contained in the data than the conventional least squares method, so it is sensitive to deterioration of the accuracy of the solution and changes in weight when performing iterative calculations. There is a problem that the convergence solution cannot be obtained.

  The present invention has been made in view of such a problem, and even when a transfer function is identified based on a frequency response function including an error, it is possible to obtain an accurate least-squares solution in consideration of the error, It is another object of the present invention to provide a model identification apparatus that can avoid deterioration of the accuracy of the least square solution without being affected by the linearization error in the process of iterative calculation, and that can obtain a convergent solution stably. Yes.

To achieve the above object, a model identification apparatus according to claim 1 of the present invention identifies a transfer function in the form of a rational function of the s operator or z operator from a frequency response function indicating the frequency characteristics of a dynamic system. In the identification device, in order to minimize the square of the difference at each frequency point between the frequency characteristic of the transfer function and the frequency characteristic indicated by the frequency response function, the unknown coefficient of the denominator polynomial and the numerator polynomial constituting the transfer function Is an unknown variable, and a first simultaneous equation having a frequency characteristic indicated by the frequency response function as a constant is created, and a first solution x obtained by using the general least square method is obtained for the first simultaneous equation. A second solution obtained by using the synthetic least-squares method for the second simultaneous equation that is generated by substituting the first solution x into the first simultaneous equation and generating the residual as a constant term. Derived e , Characterized by comprising a transfer function identification means for the coefficients of these two solutions x, the transfer function of the sum y of e.
According to this configuration, an identification process in which the sum y of the two solutions x and e is a transfer function coefficient, that is, an identification process in which the sum of the two solutions x and e is a true total least squares solution y is performed. Thus, even when the transfer function G (z) is identified based on the frequency response function H (z) including an error, an accurate least-squares solution considering the error can be obtained.

Further, in the model identification apparatus according to claim 2 of the present invention, the transfer function identification means applies an iterative calculation of a weighted least squares method to the first simultaneous equations, and the weight is set to the first equation. A first solution x obtained by using the general least squares method for the simultaneous equations, and a second simultaneous system in which the residual obtained by substituting the first solution x into the second simultaneous equations is used as a constant term The frequency of the transfer function obtained by setting the equation y based on the sum y with the solution e obtained using the general least squares method and substituting the sum y of the two solutions x and e. When the relative deviation between the response characteristic and the frequency characteristic of the frequency response function is less than the set threshold value, or the relative error between the previous solution obtained and the solution obtained this time is less than the set threshold value. It is characterized by ending in case.
According to this configuration, the weight w is updated and the transfer function coefficient (true total least squares solution) y is obtained again, and is not affected by the linearization error in the process of iterative calculation, thereby avoiding deterioration of the accuracy of the solution. Furthermore, by updating the weight w only with the true total least squares solution y, a convergent solution can be obtained stably.

  As described above, according to the present invention, even when a transfer function is identified on the basis of a frequency response function including an error, an accurate least square solution considering the error can be obtained, and in the process of iterative calculation. Thus, there is an effect that the accuracy of the least squares solution can be avoided without being affected by the linearization error, and a convergent solution can be obtained more stably.

Hereinafter, embodiments of the present invention will be described with reference to the drawings. However, parts corresponding to each other in all the drawings in this specification are denoted by the same reference numerals, and description of the overlapping parts will be omitted as appropriate.
(First embodiment)
FIG. 1 is a diagram showing a configuration of a model identification apparatus according to the first embodiment of the present invention.
A model identification apparatus 10 shown in FIG. 1 includes a function reading unit 12 and a transfer function identification unit 14, and the transfer function identification unit 14 includes a first simultaneous equation creation unit 16 and a first matrix creation unit 18. A first singular value decomposition unit 20, a first total least squares solution calculation unit 22, a residual calculation unit 24, a second simultaneous equation creation unit 26, a second matrix creation unit 28, 2 singular value decomposition unit 30, second total least square solution calculation unit 32, and total least square solution addition unit 34.

The function reading unit 12 reads a frequency response function H (z) indicating the frequency characteristic of the target dynamic system.
Based on the read frequency response function H (z), the first simultaneous equation creating unit 16 defines a transfer function G (z) assuming that the order is p, and constructs this transfer function G (z). A simultaneous equation (Formula 1) is created with the unknown coefficients of the denominator polynomial and the numerator polynomial as unknown variables and the frequency characteristics indicated by the frequency response function H (z) as constants.
The first matrix creation unit 18 creates an m-row 2p + 2-column matrix B (Formula 4) using the coefficient term matrix A and the constant term column vector d shown in Formula 1.

The first singular value decomposition unit 20 performs singular value decomposition of the matrix B shown in Expression 4 (Expression 5).

Here, the matrix U is an m-by-m orthogonal matrix composed of left singular column vectors, the matrix V is a 2p + 2 × 2p + 2 column orthogonal matrix composed of right singular vectors, and the matrix S has singular values in diagonal elements. This is a diagonal matrix with m rows and 2p + 2 columns.
The first total least squares solution calculation unit 22 calculates the first total least squares solution x of Formula 1 using the 2p + 2 column vector of the matrix V (Formula 6).

Here, v _{i, j} represents the element in the i-th row and j-th column of the matrix V.
The residual calculation unit 24 substitutes the solution x obtained by the first total least squares solution calculation unit 22 into Equation 1 to calculate the residual r of Equation 1 (Equation 7).

The second simultaneous equation creating unit 26 creates a new simultaneous equation with the variable e using the residual r obtained by Equation 7 and the coefficient term matrix A of Equation 1 (Equation 8).

Like the first matrix creation unit 18, the second matrix creation unit 28 creates a matrix C (formula 9) of m rows 2p + 2 columns using the coefficient term matrix A and the constant term column vector r. .

Like the first singular value decomposition unit 20, the second singular value decomposition unit 30 performs singular value decomposition of the matrix C shown in Expression 9 (Expression 10).

Here, the matrix Uc is an m-row m-column orthogonal matrix composed of left singular column vectors, the matrix Vc is a 2p + 2 × 2p + 2-column orthogonal matrix composed of right singular vectors, and the matrix Sc has diagonal elements with singular values. This is a diagonal matrix with m rows and 2p + 2 columns.
The second total least squares solution calculation unit 32 calculates the second total least squares solution e of Formula 8 using the 2p + 2 column column vector of the matrix Vc (Formula 11).

Here, v _{ci, j} indicates the element in the i-th row and j-th column of the matrix Vc.
The total least square solution adding unit 34 calculates the sum y of the first and second total least square solutions x and e, and sets this y as the true total least square solution of Equation 1.
The operation of the transfer function G (z) identification process by the model identification apparatus 10 having such a configuration will be described with reference to the flowchart shown in FIG.

  In step S1, the function reading unit 12 reads a frequency response function H (z) indicating the frequency characteristic of the target dynamic system. In step S2, a transfer function G (z) assuming the order as p is defined based on the frequency response function H (z) read in step S1, and a denominator polynomial and a numerator polynomial constituting the transfer function are defined. The unknown coefficient is set as an unknown variable, and simultaneous equations Ax = d of Formula 1 are created with the frequency characteristic indicated by the frequency response function H (z) as a constant.

Subsequently, in step S3, the coefficient term matrix A and the constant term column vector d shown in Formula 1 are used in the first matrix creation unit 18, and the matrix B = [A d] of m rows 2p + 2 columns in Formula 4 is obtained. Created. In step S4, the singular value decomposition of the matrix B shown in Equation 4 is performed as B = USV ^{T in} Equation 5.
In step S5, the first total least squares solution calculation unit 22 uses the column vector in the 2p + 2 column of the matrix V described above in Formula 5 to obtain the first total least squares solution x in Formula 1. It is calculated as Equation 6 below.

Next, in step S6, the solution x obtained in step S5 is substituted into Equation 1, and the residual r in Equation 1 is calculated as r = d−Ax in Equation 7.
In step S 7, a new simultaneous equation with variable e is created as Ae = r in Equation 8 using the residual r obtained in Equation 7 and the coefficient term matrix A in Equation 1.
Next, in step S8, the coefficient term matrix A and the constant term column vector r are used in the second matrix creation unit 28 as in step S3, and the matrix of C = [A r] in Expression 9 is obtained. Created.

In step S9, as in step S4, the singular value decomposition of the matrix C shown in equation 9 is performed as C = UcScVc ^{T in} equation 10.
In step S10, the second total least squares solution e of Expression 8 is calculated as Expression 11 using the column vector of the 2p + 2 column of the matrix Vc described above in Expression 9.
In step S11, the total least squares solution adding unit 34 calculates a true total least squares solution y that is the sum of the two solutions x and e obtained in steps S5 and S10, and the process ends.

  As described above, according to the model identification apparatus 10 of the first embodiment, the frequency characteristic indicated by the frequency response function H (z) is a constant, and the denominator polynomial and the numerator polynomial constituting the transfer function G (z) to be identified are A first total least squares solution x is obtained for simultaneous equations with unknown coefficients as unknown variables, and a residue is calculated based on the solution x, and a second total least squares solution e for the residue is obtained from the residue. When the transfer function G (z) is identified on the basis of the frequency response function H (z) including an error by performing a process of calculating and adding the sum of the two solutions x and e to a true total least squares solution y In the case of, it is possible to obtain an accurate least-squares solution considering an error.

(Second Embodiment)
FIG. 5 is a diagram showing a configuration of a model identification apparatus according to the second embodiment of the present invention.
The model identification device 40 shown in FIG. 5 includes a function reading unit 12, a weight initial setting unit 42, and a transfer function identification unit 44. The transfer function identification unit 44 includes the first simultaneous equation creation unit 16 and The first matrix creation unit 18, the first singular value decomposition unit 20, the first total least squares solution calculation unit 22, the residual calculation unit 24, the second simultaneous equation creation unit 26, and the second Matrix creation unit 28, second singular value decomposition unit 30, second total least square solution calculation unit 32, total least square solution addition unit 34, calculation convergence determination unit 46, weight update unit 48, It is configured with.

The weight initial setting unit 42 performs initial setting of weights (weight coefficients). However, the unit matrix is used for the initial setting of weights.
The calculation convergence determination unit 46 performs the convergence determination of the calculation for obtaining the true total least squares solution y as will be described later. When it is determined that the convergence has not been achieved, the unconvergence determination result is used as the weight update unit. Output to 48.
The convergence determination of the iterative calculation is performed by determining the convergence of the identification or the convergence of the solution and performing an abort determination. For example, when the relative deviation between the frequency response characteristic of the transfer function obtained by substituting the true total least squares solution y and the frequency characteristic of the frequency response function is below a set threshold (Formula 12), or an iterative process When the relative error between the solution obtained in the previous calculation and the solution obtained this time is below the set threshold (Formula 13), the iterative calculation is terminated.

  When the unconvergence determination result is input from the calculation convergence determination unit 46, the weight update unit 48 updates the weight (weight coefficient) and outputs it to the first simultaneous equation creation unit 16. As shown in Equation 3, the weight is a true total least squares solution y obtained in the previous calculation in the process of iterative calculation performed in the order of the components 16 to 34, 46, and 48 of the transfer function identification unit 44. Set based on. After updating this weight, the process returns to the first simultaneous equation creating unit 16 again, and the true total least squares solution y is sequentially calculated by the same process.

The operation of the transfer function G (z) identification process by the model identification device 40 having such a configuration will be described with reference to the flowchart shown in FIG.
In step S21, the function reading unit 12 reads the frequency response function H (z), and in step S22, the weight initial setting unit 42 performs initial setting with the weight w = 1. Thereafter, in steps S23 to S32, processing is performed in the same manner as in steps S2 to S11 of the first embodiment. When the true total least squares solution y is obtained in this process, the convergence determination is performed by the calculation convergence determination unit 46 in step S33. As a result, if it is determined that it has not converged, the true total least squares solution obtained in the previous calculation of the iterative process as shown in Formula 3 in the weight update unit 48 in step S34. Updated based on y. After this update, the process returns to step S23, and after the updated weight w is set in the first simultaneous equation creating unit 16, the true total least squares solution y is again obtained in the same manner as described above. Then, in step S33, calculation convergence is determined. As a result, if any of the threshold values is below, the iterative calculation is terminated and a predetermined transfer function G (z) is obtained.

  As described above, according to the model identification device 40 of the second embodiment, the error is taken into account when the transfer function is identified based on the frequency response function including the error as in the first embodiment. An accurate least squares solution can be obtained. In addition, it is possible to avoid deterioration of the accuracy of the solution without being affected by the linearization error in the process of iterative calculation to obtain the true total least squares solution y again by updating the weight w, and further, the true total By updating the weight w with only the least squares solution y, a convergent solution can be obtained stably.

It is a figure which shows the structure of the model identification apparatus which concerns on the 1st Embodiment of this invention. It is a flowchart for demonstrating operation | movement of the identification process of the transfer function by the model identification apparatus which concerns on 1st Embodiment. It is a figure which shows the structure of the model identification apparatus which concerns on the 2nd Embodiment of this invention. It is a flowchart for demonstrating operation | movement of the identification process of the transfer function by the model identification apparatus which concerns on 2nd Embodiment.

Explanation of symbols

DESCRIPTION OF SYMBOLS 10,40 Model identification apparatus 12 Function reading part 14,44 Transfer function identification part 16 1st simultaneous equation creation part 18 1st matrix creation part 20 1st singular value decomposition part 22 1st total least squares solution calculation part 24 residual calculation unit 26 second simultaneous equation creation unit 28 second matrix creation unit 30 second singular value decomposition unit 32 second total least square solution calculation unit 34 total least square solution addition unit 42 weight initial setting unit 46 Calculation convergence determination unit 48 Weight update unit

Claims (2)

In a model identification apparatus for identifying a transfer function in a rational function form of an s operator or a z operator from a frequency response function indicating a frequency characteristic of a dynamic system,
In order to minimize the square of the difference at each frequency point between the frequency characteristic of the transfer function and the frequency characteristic indicated by the frequency response function, unknown coefficients of the denominator polynomial and the numerator polynomial constituting the transfer function are used as unknown variables. A first simultaneous equation having a frequency characteristic indicated by the frequency response function as a constant is created, and a first solution x obtained by using the general least square method is derived for the first simultaneous equation, Deriving a second solution e obtained by using the general least squares method for a second simultaneous equation created by substituting the first solution x into the first simultaneous equation as a constant term. A model identification apparatus comprising: a transfer function identification unit that uses a sum y of these two solutions x and e as a coefficient of a transfer function. The transfer function identifying means applies an iterative calculation of a weighted least square method to the first simultaneous equations,
A first solution x obtained by using an overall least square method for the first simultaneous equations and a weight obtained by substituting the first solution x into the second simultaneous equations are constants. Based on the sum y with the solution e obtained using the general least squares method for the second simultaneous equation created as a term,
When the relative deviation between the frequency response characteristic of the transfer function obtained by substituting the sum y of the two solutions x and e and the frequency characteristic of the frequency response function is less than a set threshold value, or The model identification apparatus according to claim 1, wherein the model identification apparatus ends when a relative error between a solution obtained in the previous iteration and a solution obtained this time is below a set threshold value.

Images