Purpose
To compute the matrices of the positive feedback controller | Ak | Bk | K = |----|----| | Ck | Dk | for the shaped plant | A | B | G = |---|---| | C | D | in the McFarlane/Glover Loop Shaping Design Procedure.Specification
SUBROUTINE SB10ID( N, M, NP, A, LDA, B, LDB, C, LDC, D, LDD, $ FACTOR, NK, AK, LDAK, BK, LDBK, CK, LDCK, $ DK, LDDK, RCOND, IWORK, DWORK, LDWORK, BWORK, $ INFO ) C .. Scalar Arguments .. INTEGER INFO, LDA, LDAK, LDB, LDBK, LDC, LDCK, LDD, $ LDDK, LDWORK, M, N, NK, NP DOUBLE PRECISION FACTOR C .. Array Arguments .. INTEGER IWORK( * ) LOGICAL BWORK( * ) DOUBLE PRECISION A( LDA, * ), AK( LDAK, * ), B( LDB, * ), $ BK( LDBK, * ), C( LDC, * ), CK( LDCK, * ), $ D( LDD, * ), DK( LDDK, * ), DWORK( * ), $ RCOND( 2 )Arguments
Input/Output Parameters
N (input) INTEGER The order of the plant. N >= 0. M (input) INTEGER The column size of the matrix B. M >= 0. NP (input) INTEGER The row size of the matrix C. NP >= 0. A (input) DOUBLE PRECISION array, dimension (LDA,N) The leading N-by-N part of this array must contain the system state matrix A of the shaped plant. LDA INTEGER The leading dimension of the array A. LDA >= max(1,N). B (input) DOUBLE PRECISION array, dimension (LDB,M) The leading N-by-M part of this array must contain the system input matrix B of the shaped plant. LDB INTEGER The leading dimension of the array B. LDB >= max(1,N). C (input) DOUBLE PRECISION array, dimension (LDC,N) The leading NP-by-N part of this array must contain the system output matrix C of the shaped plant. LDC INTEGER The leading dimension of the array C. LDC >= max(1,NP). D (input) DOUBLE PRECISION array, dimension (LDD,M) The leading NP-by-M part of this array must contain the system matrix D of the shaped plant. LDD INTEGER The leading dimension of the array D. LDD >= max(1,NP). FACTOR (input) DOUBLE PRECISION = 1 implies that an optimal controller is required; > 1 implies that a suboptimal controller is required, achieving a performance FACTOR less than optimal. FACTOR >= 1. NK (output) INTEGER The order of the positive feedback controller. NK <= N. AK (output) DOUBLE PRECISION array, dimension (LDAK,N) The leading NK-by-NK part of this array contains the controller state matrix Ak. LDAK INTEGER The leading dimension of the array AK. LDAK >= max(1,N). BK (output) DOUBLE PRECISION array, dimension (LDBK,NP) The leading NK-by-NP part of this array contains the controller input matrix Bk. LDBK INTEGER The leading dimension of the array BK. LDBK >= max(1,N). CK (output) DOUBLE PRECISION array, dimension (LDCK,N) The leading M-by-NK part of this array contains the controller output matrix Ck. LDCK INTEGER The leading dimension of the array CK. LDCK >= max(1,M). DK (output) DOUBLE PRECISION array, dimension (LDDK,NP) The leading M-by-NP part of this array contains the controller matrix Dk. LDDK INTEGER The leading dimension of the array DK. LDDK >= max(1,M). RCOND (output) DOUBLE PRECISION array, dimension (2) RCOND(1) contains an estimate of the reciprocal condition number of the X-Riccati equation; RCOND(2) contains an estimate of the reciprocal condition number of the Z-Riccati equation.Workspace
IWORK INTEGER array, dimension (max(2*N,N*N,M,NP)) DWORK DOUBLE PRECISION array, dimension (LDWORK) On exit, if INFO = 0, DWORK(1) contains the optimal value of LDWORK. LDWORK INTEGER The dimension of the array DWORK. LDWORK >= 4*N*N + M*M + NP*NP + 2*M*N + N*NP + 4*N + max( 6*N*N + 5 + max(1,4*N*N+8*N), N*NP + 2*N ). For good performance, LDWORK must generally be larger. An upper bound of LDWORK in the above formula is LDWORK >= 10*N*N + M*M + NP*NP + 2*M*N + 2*N*NP + 4*N + 5 + max(1,4*N*N+8*N). BWORK LOGICAL array, dimension (2*N)Error Indicator
INFO INTEGER = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; = 1: the X-Riccati equation is not solved successfully; = 2: the Z-Riccati equation is not solved successfully; = 3: the iteration to compute eigenvalues or singular values failed to converge; = 4: the matrix Ip - D*Dk is singular; = 5: the matrix Im - Dk*D is singular; = 6: the closed-loop system is unstable.Method
The routine implements the formulas given in [1].References
[1] McFarlane, D. and Glover, K. A loop shaping design procedure using H_infinity synthesis. IEEE Trans. Automat. Control, vol. AC-37, no. 6, pp. 759-769, 1992.Numerical Aspects
The accuracy of the results depends on the conditioning of the two Riccati equations solved in the controller design (see the output parameter RCOND).Further Comments
NoneExample
Program Text
* SB10ID EXAMPLE PROGRAM TEXT * Copyright (c) 2002-2017 NICONET e.V. * * .. Parameters .. INTEGER NIN, NOUT PARAMETER ( NIN = 5, NOUT = 6 ) INTEGER NMAX, MMAX, PMAX PARAMETER ( NMAX = 10, MMAX = 10, PMAX = 10 ) INTEGER LDA, LDAK, LDB, LDBK, LDC, LDCK, LDD, LDDK PARAMETER ( LDA = NMAX, LDAK = NMAX, LDB = NMAX, $ LDBK = NMAX, LDC = PMAX, LDCK = MMAX, $ LDD = PMAX, LDDK = MMAX ) INTEGER LIWORK PARAMETER ( LIWORK = MAX( 2*NMAX, NMAX*NMAX, MMAX, PMAX ) ) INTEGER LDWORK PARAMETER ( LDWORK = 4*NMAX*NMAX + MMAX*MMAX + PMAX*PMAX + $ 2*MMAX*NMAX + NMAX*PMAX + 4*NMAX + $ MAX( 10*NMAX*NMAX + 8*NMAX + 5, $ NMAX*PMAX + 2*NMAX ) ) * .. Local Scalars .. DOUBLE PRECISION FACTOR INTEGER I, INFO, J, M, N, NK, NP * .. Local Arrays .. LOGICAL BWORK(2*NMAX) INTEGER IWORK(LIWORK) DOUBLE PRECISION A(LDA,NMAX), AK(LDA,NMAX), B(LDB,MMAX), $ BK(LDBK,PMAX), C(LDC,NMAX), CK(LDCK,NMAX), $ D(LDD,MMAX), DK(LDDK,PMAX), DWORK(LDWORK), $ RCOND( 2 ) * .. External Subroutines .. EXTERNAL SB10ID * .. 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, NP IF ( N.LT.0 .OR. N.GT.NMAX ) THEN WRITE ( NOUT, FMT = 99990 ) N ELSE IF ( M.LT.0 .OR. M.GT.MMAX ) THEN WRITE ( NOUT, FMT = 99989 ) M ELSE IF ( NP.LT.0 .OR. NP.GT.PMAX ) THEN WRITE ( NOUT, FMT = 99988 ) NP 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 = * ) ( ( C(I,J), J = 1,N ), I = 1,NP ) READ ( NIN, FMT = * ) ( ( D(I,J), J = 1,M ), I = 1,NP ) READ ( NIN, FMT = * ) FACTOR CALL SB10ID( N, M, NP, A, LDA, B, LDB, C, LDC, D, LDD, $ FACTOR, NK, AK, LDAK, BK, LDBK, CK, LDCK, $ DK, LDDK, RCOND, IWORK, DWORK, LDWORK, $ BWORK, INFO ) IF ( INFO.EQ.0 ) THEN WRITE ( NOUT, FMT = 99997 ) DO 10 I = 1, NK WRITE ( NOUT, FMT = 99992 ) ( AK(I,J), J = 1,NK ) 10 CONTINUE WRITE ( NOUT, FMT = 99996 ) DO 20 I = 1, NK WRITE ( NOUT, FMT = 99992 ) ( BK(I,J), J = 1,NP ) 20 CONTINUE WRITE ( NOUT, FMT = 99995 ) DO 30 I = 1, M WRITE ( NOUT, FMT = 99992 ) ( CK(I,J), J = 1,NK ) 30 CONTINUE WRITE ( NOUT, FMT = 99994 ) DO 40 I = 1, M WRITE ( NOUT, FMT = 99992 ) ( DK(I,J), J = 1,NP ) 40 CONTINUE WRITE( NOUT, FMT = 99993 ) WRITE( NOUT, FMT = 99991 ) ( RCOND(I), I = 1, 2 ) ELSE WRITE( NOUT, FMT = 99998 ) INFO END IF END IF STOP * 99999 FORMAT (' SB10ID EXAMPLE PROGRAM RESULTS',/1X) 99998 FORMAT (/' INFO on exit from SB10ID =',I2) 99997 FORMAT (/' The controller state matrix AK is'/) 99996 FORMAT (/' The controller input matrix BK is'/) 99995 FORMAT (/' The controller output matrix CK is'/) 99994 FORMAT (/' The controller matrix DK is'/) 99993 FORMAT (/' The estimated condition numbers are'/) 99992 FORMAT (10(1X,F9.4)) 99991 FORMAT ( 2(1X,D12.5)) 99990 FORMAT (/' N is out of range.',/' N = ',I5) 99989 FORMAT (/' M is out of range.',/' M = ',I5) 99988 FORMAT (/' NP is out of range.',/' NP = ',I5) ENDProgram Data
SB10ID EXAMPLE PROGRAM DATA 6 2 3 -1.0 0.0 4.0 5.0 -3.0 -2.0 -2.0 4.0 -7.0 -2.0 0.0 3.0 -6.0 9.0 -5.0 0.0 2.0 -1.0 -8.0 4.0 7.0 -1.0 -3.0 0.0 2.0 5.0 8.0 -9.0 1.0 -4.0 3.0 -5.0 8.0 0.0 2.0 -6.0 -3.0 -4.0 2.0 0.0 -5.0 -7.0 4.0 -6.0 -3.0 9.0 1.0 -2.0 1.0 -1.0 2.0 -4.0 0.0 -3.0 -3.0 0.0 5.0 -1.0 1.0 1.0 -7.0 5.0 0.0 -8.0 2.0 -2.0 1.0 -2.0 0.0 4.0 5.0 -3.0 1.0Program Results
SB10ID EXAMPLE PROGRAM RESULTS The controller state matrix AK is -39.0671 9.9293 22.2322 -27.4113 43.8655 -6.6117 3.0006 11.0878 -11.4130 15.4269 33.6805 -6.6934 -23.9953 14.1438 -33.4358 -32.3191 9.7316 25.4033 -24.0473 42.0517 -44.1655 18.7767 34.8873 -42.4369 50.8437 The controller input matrix BK is -10.2905 -16.5382 -10.9782 -4.3598 -8.7525 -5.1447 6.5962 1.8975 6.2316 -9.8770 -14.7041 -11.8778 -9.6726 -22.7309 -18.2692 The controller output matrix CK is -0.6647 -0.0599 -1.0376 0.5619 1.7297 -8.4202 3.9573 7.3094 -7.6283 10.6768 The controller matrix DK is 0.8466 0.4979 -0.6993 -1.2226 -4.8689 -4.5056 The estimated condition numbers are 0.13861D-01 0.90541D-02