Purpose
To compute for a controllable matrix pair ( A, B ) a matrix G such that the matrix A - B*G has the desired eigenstructure, specified by desired eigenvalues and free eigenvector elements. The pair ( A, B ) should be given in orthogonal canonical form as returned by the SLICOT Library routine AB01ND.Specification
SUBROUTINE SB01DD( N, M, INDCON, A, LDA, B, LDB, NBLK, WR, WI, $ Z, LDZ, Y, COUNT, G, LDG, TOL, IWORK, DWORK, $ LDWORK, INFO ) C .. Scalar Arguments .. INTEGER COUNT, INDCON, INFO, LDA, LDB, LDG, LDWORK, $ LDZ, M, N DOUBLE PRECISION TOL C .. Array Arguments .. INTEGER IWORK( * ), NBLK( * ) DOUBLE PRECISION A( LDA, * ), B( LDB, * ), DWORK( * ), $ G( LDG, * ), WI( * ), WR( * ), Y( * ), $ Z( LDZ, * )Arguments
Input/Output Parameters
N (input) INTEGER The order of the matrix A and the number of rows of the matrix B. N >= 0. M (input) INTEGER The number of columns of the matrix B. M >= 0. INDCON (input) INTEGER The controllability index of the pair ( A, B ). 0 <= INDCON <= N. A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the leading N-by-N part of this array must contain the N-by-N matrix A in orthogonal canonical form, as returned by SLICOT Library routine AB01ND. On exit, the leading N-by-N part of this array contains the real Schur form of the matrix A - B*G. The elements below the real Schur form of A are set to zero. LDA INTEGER The leading dimension of the 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 N-by-M matrix B in orthogonal canonical form, as returned by SLICOT Library routine AB01ND. On exit, the leading N-by-M part of this array contains the transformed matrix B. LDB INTEGER The leading dimension of the array B. LDB >= max(1,N). NBLK (input) INTEGER array, dimension (N) The leading INDCON elements of this array must contain the orders of the diagonal blocks in the orthogonal canonical form of A, as returned by SLICOT Library routine AB01ND. The values of these elements must satisfy the following conditions: NBLK(1) >= NBLK(2) >= ... >= NBLK(INDCON), NBLK(1) + NBLK(2) + ... + NBLK(INDCON) = N. WR (input) DOUBLE PRECISION array, dimension (N) WI (input) DOUBLE PRECISION array, dimension (N) These arrays must contain the real and imaginary parts, respectively, of the desired poles of the closed-loop system, i.e., the eigenvalues of A - B*G. The poles can be unordered, except that complex conjugate pairs of poles must appear consecutively. The elements of WI for complex eigenvalues are modified internally, but restored on exit. Z (input/output) DOUBLE PRECISION array, dimension (LDZ,N) On entry, the leading N-by-N part of this array must contain the orthogonal matrix Z generated by SLICOT Library routine AB01ND in the reduction of ( A, B ) to orthogonal canonical form. On exit, the leading N-by-N part of this array contains the orthogonal transformation matrix which reduces A - B*G to real Schur form. LDZ INTEGER The leading dimension of the array Z. LDZ >= max(1,N). Y (input) DOUBLE PRECISION array, dimension (M*N) Y contains elements which are used as free parameters in the eigenstructure design. The values of these parameters are often set by an external optimization procedure. COUNT (output) INTEGER The actual number of elements in Y used as free eigenvector and feedback matrix elements in the eigenstructure design. G (output) DOUBLE PRECISION array, dimension (LDG,N) The leading M-by-N part of this array contains the feedback matrix which assigns the desired eigenstructure of A - B*G. LDG INTEGER The leading dimension of the array G. LDG >= max(1,M).Tolerances
TOL DOUBLE PRECISION The tolerance to be used in rank determination when transforming (A, B). If the user sets TOL > 0, then the given value of TOL is used as a lower bound for the reciprocal condition number (see the description of the argument RCOND in the SLICOT routine MB03OD); a (sub)matrix whose estimated condition number is less than 1/TOL is considered to be of full rank. If the user sets TOL <= 0, then an implicitly computed, default tolerance, defined by TOLDEF = N*N*EPS, is used instead, where EPS is the machine precision (see LAPACK Library routine DLAMCH).Workspace
IWORK INTEGER array, dimension (M) DWORK DOUBLE PRECISION array, dimension (LDWORK) On exit, if INFO = 0, DWORK(1) returns the optimal value of LDWORK. LDWORK INTEGER The length of the array DWORK. LDWORK >= MAX(M*N,M*M+2*N+4*M+1). For optimum performance LDWORK should be larger.Error Indicator
INFO INTEGER = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; = 1: if the pair ( A, B ) is not controllable or the free parameters are not set appropriately.Method
The routine implements the method proposed in [1], [2].References
[1] Petkov, P.Hr., Konstantinov, M.M., Gu, D.W. and Postlethwaite, I. Optimal pole assignment design of linear multi-input systems. Report 96-11, Department of Engineering, Leicester University, 1996. [2] Petkov, P.Hr., Christov, N.D. and Konstantinov, M.M. A computational algorithm for pole assignment of linear multi input systems. IEEE Trans. Automatic Control, vol. AC-31, pp. 1044-1047, 1986.Numerical Aspects
The method implemented is backward stable.Further Comments
The eigenvalues of the real Schur form matrix As, returned in the array A, are very close to the desired eigenvalues WR+WI*i. However, the eigenvalues of the closed-loop matrix A - B*G, computed by the QR algorithm using the matrices A and B, given on entry, may be far from WR+WI*i, although the relative error norm( Z'*(A - B*G)*Z - As )/norm( As ) is close to machine accuracy. This may happen when the eigenvalue problem for the matrix A - B*G is ill-conditioned.Example
Program Text
* SB01DD 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, LDG, LDZ PARAMETER ( LDA = NMAX, LDB = NMAX, LDG = MMAX, $ LDZ = NMAX ) INTEGER LDWORK PARAMETER ( LDWORK = MAX( 3*NMAX, MMAX*NMAX, $ MMAX*MMAX + 2*NMAX + 4*MMAX + 1 ) $ ) * .. Local Scalars .. DOUBLE PRECISION TOL INTEGER COUNT, I, INDCON, INFO1, INFO2, J, M, N, NCONT CHARACTER*1 JOBZ * .. Local Arrays .. INTEGER IWORK(MMAX), NBLK(NMAX) DOUBLE PRECISION A(LDA,NMAX), B(NMAX,MMAX), DWORK(LDWORK), $ G(LDG,NMAX), WI(NMAX), WR(NMAX), Y(MMAX*NMAX), $ Z(LDZ,NMAX) * .. External Subroutines .. EXTERNAL AB01ND, SB01DD * .. 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, JOBZ IF ( N.LT.0 .OR. N.GT.NMAX ) THEN WRITE ( NOUT, FMT = 99994 ) N ELSE IF ( M.LT.0 .OR. M.GT.MMAX ) THEN WRITE ( NOUT, FMT = 99993 ) M ELSE READ ( NIN, FMT = * ) ( ( A(I,J), J = 1,N ), I = 1,N ) READ ( NIN, FMT = * ) ( ( B(I,J), J = 1,M ), I = 1,N ) READ ( NIN, FMT = * ) ( WR(I), I = 1,N ) READ ( NIN, FMT = * ) ( WI(I), I = 1,N ) READ ( NIN, FMT = * ) ( Y(I), I = 1,M*N ) * First reduce the given system to canonical form. CALL AB01ND( JOBZ, N, M, A, LDA, B, LDB, NCONT, INDCON, $ NBLK, Z, LDZ, DWORK, TOL, IWORK, DWORK(N+1), $ LDWORK-N, INFO1 ) * IF ( INFO1.EQ.0 ) THEN * Find the state feedback matrix G. CALL SB01DD( N, M, INDCON, A, LDA, B, LDB, NBLK, WR, WI, $ Z, LDZ, Y, COUNT, G, LDG, TOL, IWORK, DWORK, $ LDWORK, INFO2 ) * IF ( INFO2.NE.0 ) THEN WRITE ( NOUT, FMT = 99997 ) INFO2 ELSE WRITE ( NOUT, FMT = 99996 ) DO 10 I = 1, M WRITE ( NOUT, FMT = 99995 ) ( G(I,J), J = 1,N ) 10 CONTINUE END IF ELSE WRITE ( NOUT, FMT = 99998 ) INFO1 END IF END IF END IF STOP * 99999 FORMAT (' SB01DD EXAMPLE PROGRAM RESULTS',/1X) 99998 FORMAT (' INFO on exit from AB01ND =',I2) 99997 FORMAT (' INFO on exit from SB01DD =',I2) 99996 FORMAT (' The state feedback matrix G 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) ENDProgram Data
SB01DD EXAMPLE PROGRAM DATA 4 2 0.0 I -1.0 0.0 2.0 -3.0 1.0 -4.0 3.0 -1.0 0.0 2.0 4.0 -5.0 0.0 0.0 -1.0 -2.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 -1.0 -1.0 -1.0 -1.0 0.0 0.0 0.0 0.0 1.0 2.0 2.0 1.0 -1.0 -2.0 3.0 1.0Program Results
SB01DD EXAMPLE PROGRAM RESULTS The state feedback matrix G is -5.2339 3.1725 -15.7885 21.7043 -1.6022 0.8504 -5.1914 6.2339