discreteRiccatiIterative

Newton's method with exact line search for solving discrete algebraic Riccati equation

Information

This information is part of the Modelica Standard Library maintained by the Modelica Association.

Syntax

           X = Matrices.Utilities.discreteRiccatiIterative(A, B, R, Q, X0);
      (X, r) = Matrices.Utilities.discreteRiccatiIterative(A, B, R, Q, X0, maxSteps, eps);

Description

This function provides a Newton-like method for solving discrete-time algebraic Riccati equations. It uses Exact Line Search to improve the sometimes erratic convergence of Newton's method. Exact line search in this case means, that at each iteration i a Newton step delta_i

  X_i+1 = X_i + delta_i

is taken in the direction to minimize the Frobenius norm of the residual

  r = || A'X_i+1*A - X_i+1 - A'X_i+1*G_i*X_i+1*A + Q ||

with

                       -1
  G_i = B*(R + B'*X_i*B) *B'

Output r is the norm of the residual of the last iteration.

The inputs "maxSteps" and "eps" specify the termination of the iteration. The iteration is terminated if either maxSteps iteration steps have been performed or the relative change delta_i/X_i became smaller than eps.

With an appropriate initial value X0 a sufficiently accurate solution might be reach with a few iteration steps. Although a Lyapunov equation of order n (n is the order of the Riccati equation) is to be solved at each iteration step, the algorithm might be faster than a direct method like Matrices.discreteRiccati, since direct methods have to solve the 2*n-order Hamiltonian system equation. The algorithm is taken from [1] and [2].

References

  [1] Benner, P., Byers, R.
      An Exact Line Search Method for Solving Generalized Continuous-Time Algebraic Riccati Equations
      IEEE Transactions On Automatic Control, Vol. 43, No. 1, pp. 101-107, 1998.
  [2] Datta, B.N.
      Numerical Methods for Linear Control Systems
      Elsevier Academic Press, 2004.

Example

     A  = [0.9970,    0.0000,    0.0000,    0.0000;
           1.0000,    0.0000,    0.0000,    0.0000;
           0.0000,    1.0000,    0.0000,    0.0000;
           0.0000,    0.0000,    1.0000,    0.0000];

     B  = [0.0150;
           0.0000;
           0.0000;
           0.0000];

     R = [0.2500];

     Q = [0, 0, 0, 0;
          0, 0, 0, 0;
          0, 0, 0, 0;
          0, 0, 0, 1];

    X0=identity(4);

    (X,r) = Matrices.Utilities.discreteRiccatiIterative(A, B, R, Q, X0);

  //  X = [30.625, 0.0, 0.0, 0.0;
            0.0,   1.0, 0.0, 0.0;
            0.0,   0.0, 1.0, 0.0;
            0.0,   0.0, 0.0, 1.0];

  // r =   3.10862446895044E-015

See also

Matrices.Utilities.continuousRiccatiIterative
Matrices.discreteRiccati

Syntax

(X, r) = discreteRiccatiIterative(A, B, R, Q, X0, maxSteps, eps)

Inputs (7)

A

Type: Real[:,size(A, 1)]

Description: Matrix A of discrete Riccati equation

B

Type: Real[size(A, 1),:]

Description: Matrix B of discrete Riccati equation

R

Default Value: identity(size(B, 2))

Type: Real[size(B, 2),size(B, 2)]

Description: Matrix R of discrete Riccati equation

Q

Default Value: identity(size(A, 1))

Type: Real[size(A, 1),size(A, 2)]

Description: Matrix Q of discrete Riccati equation

X0

Default Value: identity(size(A, 1))

Type: Real[size(A, 1),size(A, 2)]

Description: Initial approximate solution discrete Riccati equation

maxSteps

Default Value: 10

Type: Integer

Description: Maximal number of iteration steps

eps

Default Value: Matrices.frobeniusNorm(A) * 1e-9

Type: Real

Description: Tolerance for stop criterion

Outputs (2)

X

Type: Real[size(X0, 1),size(X0, 2)]

r

Type: Real