## SB10MD

### D-step in the D-K iteration for continuous-time case

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

Purpose

```  To perform the D-step in the D-K iteration. It handles
continuous-time case.

```
Specification
```      SUBROUTINE SB10MD( NC, MP, LENDAT, F, ORD, MNB, NBLOCK, ITYPE,
\$                   QUTOL, A, LDA, B, LDB, C, LDC, D, LDD, OMEGA,
\$                   MJU, IWORK, LIWORK, DWORK, LDWORK, ZWORK,
\$                   LZWORK, INFO )
C     .. Scalar Arguments ..
INTEGER           F, INFO, LDA, LDAD, LDB, LDBD, LDC, LDCD, LDD,
\$                  LDDD, LDWORK, LENDAT, LIWORK, LZWORK, MNB, MP,
\$                  NC, ORD, TOTORD
DOUBLE PRECISION  QUTOL
C     .. Array Arguments ..
INTEGER           ITYPE(*), IWORK(*), NBLOCK(*)
\$                  C(LDC, *), CD(LDCD, *), D(LDD, *), DD(LDDD, *),
\$                  DWORK(*), MJU(*), OMEGA(*)
COMPLEX*16        ZWORK(*)

```
Arguments

Input/Output Parameters

```  NC      (input) INTEGER
The order of the matrix A.  NC >= 0.

MP      (input) INTEGER
The order of the matrix D.  MP >= 0.

LENDAT  (input) INTEGER
The length of the vector OMEGA.  LENDAT >= 2.

F       (input) INTEGER
The number of the measurements and controls, i.e.,
the size of the block I_f in the D-scaling system.
F >= 0.

ORD     (input/output) INTEGER
The MAX order of EACH block in the fitting procedure.
ORD <= LENDAT-1.
On exit, if ORD < 1 then ORD = 1.

MNB     (input) INTEGER
The number of diagonal blocks in the block structure of
the uncertainty, and the length of the vectors NBLOCK
and ITYPE.  1 <= MNB <= MP.

NBLOCK  (input) INTEGER array, dimension (MNB)
The vector of length MNB containing the block structure
of the uncertainty. NBLOCK(I), I = 1:MNB, is the size of
each block.

ITYPE   (input) INTEGER array, dimension (MNB)
The vector of length MNB indicating the type of each
block.
For I = 1 : MNB,
ITYPE(I) = 1 indicates that the corresponding block is a
real block. IN THIS CASE ONLY MJU(JW) WILL BE ESTIMATED
CORRECTLY, BUT NOT D(S)!
ITYPE(I) = 2 indicates that the corresponding block is a
complex block. THIS IS THE ONLY ALLOWED VALUE NOW!
NBLOCK(I) must be equal to 1 if ITYPE(I) is equal to 1.

QUTOL   (input) DOUBLE PRECISION
The acceptable mean relative error between the D(jw) and
the frequency responce of the estimated block
[ADi,BDi;CDi,DDi]. When it is reached, the result is
taken as good enough.
A good value is QUTOL = 2.0.
If QUTOL < 0 then only mju(jw) is being estimated,
not D(s).

A       (input/output) DOUBLE PRECISION array, dimension (LDA,NC)
On entry, the leading NC-by-NC part of this array must
contain the A matrix of the closed-loop system.
On exit, if MP > 0, the leading NC-by-NC part of this
array contains an upper Hessenberg matrix similar to A.

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

B       (input/output) DOUBLE PRECISION array, dimension (LDB,MP)
On entry, the leading NC-by-MP part of this array must
contain the B matrix of the closed-loop system.
On exit, the leading NC-by-MP part of this array contains
the transformed B matrix corresponding to the Hessenberg
form of A.

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

C       (input/output) DOUBLE PRECISION array, dimension (LDC,NC)
On entry, the leading MP-by-NC part of this array must
contain the C matrix of the closed-loop system.
On exit, the leading MP-by-NC part of this array contains
the transformed C matrix corresponding to the Hessenberg
form of A.

LDC     INTEGER
The leading dimension of the array C.  LDC >= MAX(1,MP).

D       (input) DOUBLE PRECISION array, dimension (LDD,MP)
The leading MP-by-MP part of this array must contain the
D matrix of the closed-loop system.

LDD     INTEGER
The leading dimension of the array D.  LDD >= MAX(1,MP).

OMEGA   (input) DOUBLE PRECISION array, dimension (LENDAT)
The vector with the frequencies.

TOTORD  (output) INTEGER
The TOTAL order of the D-scaling system.
TOTORD is set to zero, if QUTOL < 0.

The leading TOTORD-by-TOTORD part of this array contains
the A matrix of the D-scaling system.
Not referenced if QUTOL < 0.

LDAD >= MAX(1,MP*ORD), if QUTOL >= 0;
LDAD >= 1,             if QUTOL <  0.

BD      (output) DOUBLE PRECISION array, dimension (LDBD,MP+F)
The leading TOTORD-by-(MP+F) part of this array contains
the B matrix of the D-scaling system.
Not referenced if QUTOL < 0.

LDBD    INTEGER
The leading dimension of the array BD.
LDBD >= MAX(1,MP*ORD), if QUTOL >= 0;
LDBD >= 1,             if QUTOL <  0.

CD      (output) DOUBLE PRECISION array, dimension (LDCD,MP*ORD)
The leading (MP+F)-by-TOTORD part of this array contains
the C matrix of the D-scaling system.
Not referenced if QUTOL < 0.

LDCD    INTEGER
The leading dimension of the array CD.
LDCD >= MAX(1,MP+F), if QUTOL >= 0;
LDCD >= 1,           if QUTOL <  0.

DD      (output) DOUBLE PRECISION array, dimension (LDDD,MP+F)
The leading (MP+F)-by-(MP+F) part of this array contains
the D matrix of the D-scaling system.
Not referenced if QUTOL < 0.

LDDD    INTEGER
The leading dimension of the array DD.
LDDD >= MAX(1,MP+F), if QUTOL >= 0;
LDDD >= 1,           if QUTOL <  0.

MJU     (output) DOUBLE PRECISION array, dimension (LENDAT)
The vector with the upper bound of the structured
singular value (mju) for each frequency in OMEGA.

```
Workspace
```  IWORK   INTEGER array, dimension (LIWORK)

LIWORK  INTEGER
The length of the array IWORK.
LIWORK >= MAX( NC, 4*MNB-2, MP, 2*ORD+1 ), if QUTOL >= 0;
LIWORK >= MAX( NC, 4*MNB-2, MP ),          if QUTOL <  0.

DWORK   DOUBLE PRECISION array, dimension (LDWORK)
On exit, if INFO = 0, DWORK(1) returns the optimal value
of LDWORK, DWORK(2) returns the optimal value of LZWORK,
and DWORK(3) returns an estimate of the minimum reciprocal
of the condition numbers (with respect to inversion) of
the generated Hessenberg matrices.

LDWORK  INTEGER
The length of the array DWORK.
LDWORK >= MAX( 3, LWM, LWD ), where
LWM = LWA + MAX( NC + MAX( NC, MP-1 ),
2*MP*MP*MNB - MP*MP + 9*MNB*MNB +
MP*MNB + 11*MP + 33*MNB - 11 );
LWD = LWB + MAX( 2, LW1, LW2, LW3, LW4, 2*ORD ),
if QUTOL >= 0;
LWD = 0,         if QUTOL <  0;
LWA = MP*LENDAT + 2*MNB + MP - 1;
LWB = LENDAT*(MP + 2) + ORD*(ORD + 2) + 1;
LW1 = 2*LENDAT + 4*HNPTS;  HNPTS = 2048;
LW2 =   LENDAT + 6*HNPTS;  MN  = MIN( 2*LENDAT, 2*ORD+1 );
LW3 = 2*LENDAT*(2*ORD + 1) + MAX( 2*LENDAT, 2*ORD + 1 ) +
MAX( MN + 6*ORD + 4, 2*MN + 1 );
LW4 = MAX( ORD*ORD + 5*ORD, 6*ORD + 1 + MIN( 1, ORD ) ).

ZWORK   COMPLEX*16 array, dimension (LZWORK)

LZWORK  INTEGER
The length of the array ZWORK.
LZWORK >= MAX( LZM, LZD ), where
LZM = MAX( MP*MP + NC*MP + NC*NC + 2*NC,
6*MP*MP*MNB + 13*MP*MP + 6*MNB + 6*MP - 3 );
LZD = MAX( LENDAT*(2*ORD + 3), ORD*ORD + 3*ORD + 1 ),
if QUTOL >= 0;
LZD = 0,         if QUTOL <  0.

```
Error Indicator
```  INFO    (output) INTEGER
=  0:  successful exit;
<  0:  if INFO = -i, the i-th argument had an illegal
value;
=  1:  if one or more values w in OMEGA are (close to
some) poles of the closed-loop system, i.e., the
matrix jw*I - A is (numerically) singular;
=  2:  the block sizes must be positive integers;
=  3:  the sum of block sizes must be equal to MP;
=  4:  the size of a real block must be equal to 1;
=  5:  the block type must be either 1 or 2;
=  6:  errors in solving linear equations or in matrix
inversion;
=  7:  errors in computing eigenvalues or singular values.
= 1i:  INFO on exit from SB10YD is i. (1i means 10 + i.)

```
Method
```  I.   First, W(jw) for the given closed-loop system is being
estimated.
II.  Now, AB13MD SLICOT subroutine can obtain the D(jw) scaling
system with respect to NBLOCK and ITYPE, and colaterally,
mju(jw).
If QUTOL < 0 then the estimations stop and the routine exits.
III. Now that we have D(jw), SB10YD subroutine can do block-by-
block fit. For each block it tries with an increasing order
of the fit, starting with 1 until the
between the Dii(jw) and the estimated frequency responce
of the block becomes less than or equal to the routine
argument QUTOL, or the order becomes equal to ORD.
IV.  Arrange the obtained blocks in the AD, BD, CD and DD
matrices and estimate the total order of D(s), TOTORD.
V.   Add the system I_f to the system obtained in IV.

```
References
```  [1] Balas, G., Doyle, J., Glover, K., Packard, A. and Smith, R.
Mu-analysis and Synthesis toolbox - User's Guide,
The Mathworks Inc., Natick, MA, USA, 1998.

```
```  None
```
Example

Program Text

```  None
```
Program Data
```  None
```
Program Results
```  None
```