**Purpose**

To compute a reduced order model by using singular perturbation approximation formulas.

SUBROUTINE AB09DD( DICO, N, M, P, NR, A, LDA, B, LDB, C, LDC, $ D, LDD, RCOND, IWORK, DWORK, INFO ) C .. Scalar Arguments .. CHARACTER DICO INTEGER INFO, LDA, LDB, LDC, LDD, M, N, NR, P DOUBLE PRECISION RCOND C .. Array Arguments .. DOUBLE PRECISION A(LDA,*), B(LDB,*), C(LDC,*), D(LDD,*), DWORK(*) INTEGER IWORK(*)

**Mode Parameters**

DICO CHARACTER*1 Specifies the type of the original system as follows: = 'C': continuous-time system; = 'D': discrete-time system.

N (input) INTEGER The dimension of the state vector, i.e. the order of the matrix A; also the number of rows of matrix B and the number of columns of the matrix C. N >= 0. M (input) INTEGER The dimension of input vector, i.e. the number of columns of matrices B and D. M >= 0. P (input) INTEGER The dimension of output vector, i.e. the number of rows of matrices C and D. P >= 0. NR (input) INTEGER The order of the reduced order system. N >= NR >= 0. A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the leading N-by-N part of this array must contain the state dynamics matrix of the original system. On exit, the leading NR-by-NR part of this array contains the state dynamics matrix Ar of the reduced order system. 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 input/state matrix of the original system. On exit, the leading NR-by-M part of this array contains the input/state matrix Br of the reduced order system. LDB INTEGER The leading dimension of array B. LDB >= MAX(1,N). C (input/output) DOUBLE PRECISION array, dimension (LDC,N) On entry, the leading P-by-N part of this array must contain the state/output matrix of the original system. On exit, the leading P-by-NR part of this array contains the state/output matrix Cr of the reduced order system. LDC INTEGER The leading dimension of array C. LDC >= MAX(1,P). D (input/output) DOUBLE PRECISION array, dimension (LDD,M) On entry, the leading P-by-M part of this array must contain the input/output matrix of the original system. On exit, the leading P-by-M part of this array contains the input/output matrix Dr of the reduced order system. If NR = 0 and the given system is stable, then D contains the steady state gain of the system. LDD INTEGER The leading dimension of array D. LDD >= MAX(1,P). RCOND (output) DOUBLE PRECISION The reciprocal condition number of the matrix A22-g*I (see METHOD).

IWORK INTEGER array, dimension (2*(N-NR)) DWORK DOUBLE PRECISION array, dimension (4*(N-NR))

INFO INTEGER = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; = 1: if the matrix A22-g*I (see METHOD) is numerically singular.

Given the system (A,B,C,D), partition the system matrices as ( A11 A12 ) ( B1 ) A = ( ) , B = ( ) , C = ( C1 C2 ), ( A21 A22 ) ( B2 ) where A11 is NR-by-NR, B1 is NR-by-M, C1 is P-by-NR, and the other submatrices have appropriate dimensions. The matrices of the reduced order system (Ar,Br,Cr,Dr) are computed according to the following residualization formulas: -1 -1 Ar = A11 + A12*(g*I-A22) *A21 , Br = B1 + A12*(g*I-A22) *B2 -1 -1 Cr = C1 + C2*(g*I-A22) *A21 , Dr = D + C2*(g*I-A22) *B2 where g = 0 if DICO = 'C' and g = 1 if DICO = 'D'.

None

**Program Text**

* AB09DD 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 = 20, MMAX = 20, PMAX = 20 ) INTEGER LDA, LDB, LDC, LDD PARAMETER ( LDA = NMAX, LDB = NMAX, LDC = PMAX, $ LDD = PMAX ) INTEGER LIWORK PARAMETER ( LIWORK = 2*NMAX ) INTEGER LDWORK PARAMETER ( LDWORK = 4*NMAX ) * .. Local Scalars .. DOUBLE PRECISION RCOND INTEGER I, INFO, J, M, N, NR, P CHARACTER*1 DICO * .. Local Arrays .. DOUBLE PRECISION A(LDA,NMAX), B(LDB,MMAX), C(LDC,NMAX), $ D(LDD,MMAX), DWORK(LDWORK) INTEGER IWORK(LIWORK) * .. External Subroutines .. EXTERNAL AB09DD * .. 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, P, NR, DICO IF ( N.LT.0 .OR. N.GT.NMAX ) THEN WRITE ( NOUT, FMT = 99990 ) N ELSE READ ( NIN, FMT = * ) ( ( A(I,J), J = 1,N ), I = 1,N ) IF ( M.LT.0 .OR. M.GT.MMAX ) THEN WRITE ( NOUT, FMT = 99989 ) M ELSE READ ( NIN, FMT = * ) ( ( B(I,J), J = 1,M ), I = 1, N ) IF ( P.LT.0 .OR. P.GT.PMAX ) THEN WRITE ( NOUT, FMT = 99988 ) P ELSE READ ( NIN, FMT = * ) ( ( C(I,J), J = 1,N ), I = 1,P ) READ ( NIN, FMT = * ) ( ( D(I,J), J = 1,M ), I = 1,P ) * Find a reduced ssr for (A,B,C). CALL AB09DD( DICO, N, M, P, NR, A, LDA, B, LDB, C, LDC, $ D, LDD, RCOND, IWORK, DWORK, INFO ) * IF ( INFO.NE.0 ) THEN WRITE ( NOUT, FMT = 99998 ) INFO ELSE WRITE ( NOUT, FMT = 99997 ) RCOND IF( NR.GT.0 ) WRITE ( NOUT, FMT = 99996 ) DO 20 I = 1, NR WRITE ( NOUT, FMT = 99995 ) ( A(I,J), J = 1,NR ) 20 CONTINUE IF( NR.GT.0 ) WRITE ( NOUT, FMT = 99993 ) DO 40 I = 1, NR WRITE ( NOUT, FMT = 99995 ) ( B(I,J), J = 1,M ) 40 CONTINUE IF( NR.GT.0 ) WRITE ( NOUT, FMT = 99992 ) DO 60 I = 1, P WRITE ( NOUT, FMT = 99995 ) ( C(I,J), J = 1,NR ) 60 CONTINUE WRITE ( NOUT, FMT = 99991 ) DO 70 I = 1, P WRITE ( NOUT, FMT = 99995 ) ( D(I,J), J = 1,M ) 70 CONTINUE END IF END IF END IF END IF STOP * 99999 FORMAT (' AB09DD EXAMPLE PROGRAM RESULTS',/1X) 99998 FORMAT (' INFO on exit from AB09DD = ',I2) 99997 FORMAT (' The computed reciprocal condition number = ',1PD12.5) 99996 FORMAT (/' The reduced state dynamics matrix Ar is ') 99995 FORMAT (20(1X,F8.4)) 99993 FORMAT (/' The reduced input/state matrix Br is ') 99992 FORMAT (/' The reduced state/output matrix Cr is ') 99991 FORMAT (/' The reduced input/output matrix Dr is ') 99990 FORMAT (/' N is out of range.',/' N = ',I5) 99989 FORMAT (/' M is out of range.',/' M = ',I5) 99988 FORMAT (/' P is out of range.',/' P = ',I5) END

AB09DD EXAMPLE PROGRAM DATA (Continuous system) 7 2 3 5 C -0.04165 4.9200 -4.9200 0 0 0 0 0 -3.3300 0 0 0 3.3300 0 0.5450 0 0 -0.5450 0 0 0 0 0 4.9200 -0.04165 4.9200 0 0 0 0 0 0 -3.3300 0 3.3300 -5.2100 0 0 0 0 -12.5000 0 0 0 0 -5.2100 0 0 -12.5000 0 0 0 0 0 0 0 0 0 0 12.5000 0 0 12.5000 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0

AB09DD EXAMPLE PROGRAM RESULTS The computed reciprocal condition number = 1.00000D+00 The reduced state dynamics matrix Ar is -0.0416 4.9200 -4.9200 0.0000 0.0000 -1.3879 -3.3300 0.0000 0.0000 0.0000 0.5450 0.0000 0.0000 -0.5450 0.0000 0.0000 0.0000 4.9200 -0.0416 4.9200 0.0000 0.0000 0.0000 -1.3879 -3.3300 The reduced input/state matrix Br is 0.0000 0.0000 3.3300 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 3.3300 The reduced state/output matrix Cr is 1.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 1.0000 0.0000 0.0000 0.0000 0.0000 0.0000 1.0000 0.0000 The reduced input/output matrix Dr is 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000