SB06ND

Minimum norm feedback matrix for "deadbeat control" of a state-space representation

[Specification] [Arguments] [Method] [References] [Comments] [Example]

Purpose

  To construct the minimum norm feedback matrix F to perform
  "deadbeat control" on a (A,B)-pair of a state-space model (which
  must be preliminarily reduced to upper "staircase" form using
  SLICOT Library routine AB01OD) such that the matrix R = A + BFU'
  is nilpotent.
  (The transformation matrix U reduces R to upper Schur form with
  zero blocks on its diagonal (of dimension KSTAIR(i)) and
  therefore contains bases for the i-th controllable subspaces,
  where i = 1,...,KMAX).

Specification
      SUBROUTINE SB06ND( N, M, KMAX, A, LDA, B, LDB, KSTAIR, U, LDU, F,
     $                   LDF, DWORK, INFO )
C     .. Scalar Arguments ..
      INTEGER           INFO, KMAX, LDA, LDB, LDF, LDU, M, N
C     .. Array Arguments ..
      INTEGER           KSTAIR(*)
      DOUBLE PRECISION  A(LDA,*), B(LDB,*), DWORK(*), F(LDF,*), U(LDU,*)

Arguments

Input/Output Parameters

  N       (input) INTEGER
          The actual state dimension, i.e. the order of the
          matrix A.  N >= 0.

  M       (input) INTEGER
          The actual input dimension.  M >= 0.

  KMAX    (input) INTEGER
          The number of "stairs" in the staircase form as produced
          by SLICOT Library routine AB01OD.  0 <= KMAX <= N.

  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the leading N-by-N part of this array must
          contain the transformed state-space matrix of the
          (A,B)-pair with triangular stairs, as produced by SLICOT
          Library routine AB01OD (with option STAGES = 'A').
          On exit, the leading N-by-N part of this array contains
          the matrix U'AU + U'BF.

  LDA     INTEGER
          The leading dimension of array A.  LDA >= MAX(1,N).

  B       (input/output) DOUBLE PRECISION array, dimension (LDB,M)
          On entry, the leading N-by-M part of this array must
          contain the transformed triangular input matrix of the
          (A,B)-pair as produced by SLICOT Library routine AB01OD
          (with option STAGES = 'A').
          On exit, the leading N-by-M part of this array contains
          the matrix U'B.

  LDB     INTEGER
          The leading dimension of array B.  LDB >= MAX(1,N).

  KSTAIR  (input) INTEGER array, dimension (KMAX)
          The leading KMAX elements of this array must contain the
          dimensions of each "stair" as produced by SLICOT Library
          routine AB01OD.

  U       (input/output) DOUBLE PRECISION array, dimension (LDU,N)
          On entry, the leading N-by-N part of this array must
          contain either a transformation matrix (e.g. from a
          previous call to other SLICOT routine) or be initialised
          as the identity matrix.
          On exit, the leading N-by-N part of this array contains
          the product of the input matrix U and the state-space
          transformation matrix which reduces A + BFU' to real
          Schur form.

  LDU     INTEGER
          The leading dimension of array U.  LDU >= MAX(1,N).

  F       (output) DOUBLE PRECISION array, dimension (LDF,N)
          The leading M-by-N part of this array contains the
          deadbeat feedback matrix F.

  LDF     INTEGER
          The leading dimension of array F.  LDF >= MAX(1,M).

Workspace
  DWORK   DOUBLE PRECISION array, dimension (2*N)

Error Indicator
  INFO    INTEGER
          = 0:  successful exit;
          < 0:  if INFO = -i, the i-th argument had an illegal
                value.

Method
  Starting from the (A,B)-pair in "staircase form" with "triangular"
  stairs, dimensions KSTAIR(i+1) x KSTAIR(i), (described by the
  vector KSTAIR):

                 | B | A      *  . . .  *  |
                 |  1|  11       .      .  |
                 |   | A     A     .    .  |
                 |   |  21    22     .  .  |
                 |   |    .      .     .   |
   [ B | A ]  =  |   |      .      .    *  |
                 |   |        .      .     |
                 | 0 |   0                 |
                 |   |          A      A   |
                 |   |           r,r-1  rr |

  where the i-th diagonal block of A has dimension KSTAIR(i), for
  i = 1,2,...,r, the feedback matrix F is constructed recursively in
  r steps (where the number of "stairs" r is given by KMAX). In each
  step a unitary state-space transformation U and a part of F are
  updated in order to achieve the final form:

                    | 0   A      *   . . .  *   |
                    |      12      .        .   |
                    |                .      .   |
                    |     0    A       .    .   |
                    |           23       .  .   |
                    |         .      .          |
  [ U'AU + U'BF ] = |           .      .    *   | .
                    |             .      .      |
                    |                           |
                    |                     A     |
                    |                      r-1,r|
                    |                           |
                    |                       0   |

References
  [1] Van Dooren, P.
      Deadbeat control: a special inverse eigenvalue problem.
      BIT, 24, pp. 681-699, 1984.

Numerical Aspects
  The algorithm requires O((N + M) * N**2) operations and is mixed
  numerical stable (see [1]).

Further Comments
  None
Example

Program Text

*     SB06ND EXAMPLE PROGRAM TEXT
*     Copyright (c) 2002-2017 NICONET e.V.
*
*     .. Parameters ..
      INTEGER          NIN, NOUT
      PARAMETER        ( NIN = 5, NOUT = 6 )
      INTEGER          NMAX, MMAX
      PARAMETER        ( NMAX = 20, MMAX = 20 )
      INTEGER          LDA, LDB, LDU, LDV, LDF
      PARAMETER        ( LDA = NMAX, LDB = NMAX, LDU = NMAX,
     $                   LDV = MMAX, LDF = MMAX )
      INTEGER          LIWORK
      PARAMETER        ( LIWORK = MMAX )
      INTEGER          LDWORK
      PARAMETER        ( LDWORK = NMAX + MAX(NMAX,3*MMAX) )
*     PARAMETER        ( LDWORK = 4*NMAX)
      DOUBLE PRECISION ZERO, ONE
      PARAMETER        ( ZERO = 0.0D0, ONE = 1.0D0 )
*     .. Local Scalars ..
      DOUBLE PRECISION TOL
      INTEGER          I, INFO, J, KMAX, M, N, NCONT
      CHARACTER*1      JOBU, JOBV
*     .. Local Arrays ..
      DOUBLE PRECISION A(LDA,NMAX), B(LDB,MMAX), DWORK(LDWORK),
     $                 F(LDF,NMAX), U(LDU,NMAX), V(LDV,MMAX)
      INTEGER          IWORK(LIWORK), KSTAIR(NMAX)
C     .. External Functions ..
      LOGICAL          LSAME
      EXTERNAL         LSAME
*     .. External Subroutines ..
      EXTERNAL         AB01OD, DLASET, SB06ND
*     .. Intrinsic Functions ..
      INTRINSIC        MAX
*     .. Executable Statements ..
*
      WRITE ( NOUT, FMT = 99999 )
*     Skip the heading in the data file and read the data.
      READ ( NIN, FMT = '()' )
      READ ( NIN, FMT = * ) N, M, TOL, JOBU, JOBV
      IF ( N.LT.0 .OR. N.GT.NMAX ) THEN
         WRITE ( NOUT, FMT = 99994 ) N
      ELSE
         READ ( NIN, FMT = * ) ( ( A(I,J), I = 1,N ), J = 1,N )
         IF ( M.LT.0 .OR. M.GT.MMAX ) THEN
            WRITE ( NOUT, FMT = 99993 ) M
         ELSE
            READ ( NIN, FMT = * ) ( ( B(I,J), J = 1,M ), I = 1,N )
*           First put (A,B) into staircase form with triangular pivots
*           and determine the stairsizes.
            CALL AB01OD( 'A', JOBU, JOBV, N, M, A, LDA, B, LDB, U,
     $                   LDU, V, LDV, NCONT, KMAX, KSTAIR, TOL, IWORK,
     $                   DWORK, LDWORK, INFO )
*
            IF ( INFO.EQ.0 ) THEN
               IF( LSAME( JOBU, 'N' ) ) THEN
*                 Initialize U as the identity matrix.
                  CALL DLASET( 'Full', N, N, ZERO, ONE, U, LDU )
               END IF
*              Perform "deadbeat control" to give F.
               CALL SB06ND( N, M, KMAX, A, LDA, B, LDB, KSTAIR, U, LDU,
     $                      F, LDF, DWORK, INFO )
*
               IF ( INFO.NE.0 ) THEN
                  WRITE ( NOUT, FMT = 99997 ) INFO
               ELSE
                  WRITE ( NOUT, FMT = 99996 )
                  DO 60 I = 1, M
                     WRITE ( NOUT, FMT = 99995 ) ( F(I,J), J = 1,N )
   60             CONTINUE
               END IF
            ELSE
               WRITE ( NOUT, FMT = 99998 ) INFO
            END IF
         END IF
      END IF
      STOP
*
99999 FORMAT (' SB06ND EXAMPLE PROGRAM RESULTS',/1X)
99998 FORMAT (' INFO on exit from AB01OD = ',I2)
99997 FORMAT (' INFO on exit from SB06ND = ',I2)
99996 FORMAT (' The deadbeat feedback matrix F is ')
99995 FORMAT (20(1X,F8.4))
99994 FORMAT (/' N is out of range.',/' N = ',I5)
99993 FORMAT (/' M is out of range.',/' M = ',I5)
      END
Program Data
 SB06ND EXAMPLE PROGRAM DATA
   5     2     0.0     N     N
  -17.0   24.0   41.0   68.0   15.0
   23.0  -35.0   27.0   14.0   16.0
   34.0   26.0  -13.0   20.0   22.0
   10.0   12.0   19.0  -21.0   63.0
   11.0   18.0   25.0   52.0  -29.0
  -31.0   14.0
   74.0  -69.0
  -59.0   16.0
   16.0  -25.0
  -25.0   36.0
Program Results
 SB06ND EXAMPLE PROGRAM RESULTS

 The deadbeat feedback matrix F is 
  -0.4819  -0.5782  -2.7595  -3.1093   0.0000
   0.2121  -0.4462   0.7698  -1.5421  -0.5773

Return to index