SUBROUTINE MVZERO (L,M,N,NMAX,LMAX,MAX,A,B,C,D,G,H,NZER,NRANK, + ZR,ZI,BETA,SUM,DUMMY) C C ***** PURPOSE: C THIS SUBROUTINE COMPUTES THE FINITE MULTIVARIABLE ZEROS OF THE C LINEAR SYSTEM C C X'(T) = A*X(T) + B*U(T) C Y(T) = C*X(T) + D*U(T) C C WHERE A IS N X N, B IS N X M, C IS L X N, AND D IS L X M. THE C ALGORITHM EMPLOYED IS DESCRIBED BELOW UNDER ALGORITHM NOTES. C C ***** CALLING SEQUENCE PARAMETERS: INTEGER L,M,N,NMAX,LMAX,MAX,NZER,NRANK DOUBLE PRECISION A(NMAX,N),B(NMAX,M),C(LMAX,N),D(LMAX,M),G(MAX,1), + H(MAX,1),ZR(1),ZI(1),BETA(1),SUM(1),DUMMY(1) C C ***** LOCAL VARIABLES: INTEGER I,IERR,J DOUBLE PRECISION EPS,EPSP1,EPS1,GNORM,T LOGICAL MATZ C C ***** FORTRAN FUNCTIONS CALLED: DOUBLE PRECISION DABS C C ***** SUBROUTINES CALLED: C ZEROS,QZHES(EISPACK),QZIT(EISPACK),QZVAL(EISPACK) C C ------------------------------------------------------------------ C C ***** PARAMETER DESCRIPTION: C ON INPUT: C C L,M,N INTEGERS C NUMBERS OF OUTPUTS, INPUTS, AND STATES, RESPEC- C TIVELY, OF THE LINEAR SYSTEM WHOSE ZEROS ARE TO BE C DETERMINED. C C NMAX,LMAX INTEGERS C ROW DIMENSIONS OF THE ARRAYS CONTAINING A,B AND C,D, C RESPECTIVELY, AS DECLARED IN THE MAIN CALLING C PROGRAM DIMENSION STATEMENT. C C MAX INTEGER C ROW DIMENSION OF THE (N+L) X (N+M) ARRAYS G,H AS C DECLARED IN THE MAIN CALLING PROGRAM DIMENSION C STATEMENT. C *** NOTE THAT MAX MUST BE .GE. (N+L) C C A,B,C,D DOUBLE PRECISION (NMAX,N),(NMAX,M),(LMAX,N),(LMAX,M) C ARRAYS CONTAINING THE SYSTEM MATRICES; NOTE THAT C THESE ARRAYS ARE NOT MODIFIED BY THE SUBROUTINE. C C G,H DOUBLE PRECISION (MAX,1),(MAX,1) C SCRATCH ARRAYS OF SIZE (N+L) X (N+M). C C BETA,SUM, DOUBLE PRECISION (1),(1),(1) C DUMMY SCRATCH VECTORS OF LENGTHS N, MAX(L,M), MAX(L,M,N), C RESPECTIVELY. NOTE: LENGTH N WILL SUFFICE FOR ALL C THREE VECTORS ASSUMING THERE ARE AT LEAST AS MANY C STATES AS INPUTS OR OUTPUTS. C C ON OUTPUT: C C NZER INTEGER C = K ON NORMAL RETURN WITH K (FROM 0 TO N) C FINITE MULTIVARIABLE ZEROS. C = -K ON ABNORMAL RETURN IF THE QZ ALGORITHM C APPLIED TO THE FINITE PENCIL FAILED TO C CONVERGE AT ROW K. C C NRANK INTEGER C THE NORMAL RANK OF THE SYSTEM TRANSFER FUNCTION; C GENERICALLY, THIS NUMBER WILL BE MIN(L,M); IF C NRANK IS RETURNED AS 0 THE SYSTEM IS SAID TO BE C DEGENERATE BUT MAY STILL HAVE MEANINGFUL ZEROS (CF. C REFERENCES (1) AND (2)). C C ZR,ZI DOUBLE PRECISION (1),(1) C VECTORS (OF LENGTH AT MOST N) CONTAINING THE REAL C AND IMAGINARY PARTS, RESPECTIVELY, OF THE NZER C FINITE MULTIVARIABLE ZEROS. C C ***** ALGORITHM NOTES: C 1. THE ZEROS ARE DETERMINED BY MEANS OF THE QZ ALGORITHM C (SEE REF. (3)) APPLIED TO A FINITE MATRIX PENCIL WHICH C IS DETERMINED AS A RESULT OF AN ORTHOGONAL STRUCTURE C ALGORITHM APPLIED TO THE ORIGINAL SYSTEM PENCIL; REFS. C (1),(2), AND (3) SHOULD ALL BE CONSULTED FOR FURTHER C ALGORITHMIC DETAILS. C 2. TWO IMPORTANT ARRAYS ARE PRODUCED BY THE QZ ALGORITHM: C THE ALPHA'S WHICH ARE COMPLEX AND THE BETA'S WHICH ARE C REAL. THE ZEROS ARE DETERMINED AS RATIOS (ALPHA/BETA) C OF THESE NUMBERS. THE BETA'S ARE AVAILABLE IF DESIRED C IN THE SCRATCH VECTOR BETA AND THE ALPHA'S MAY BE C RECOVERED BY MULTIPLYING COMPONENTS OF ZR AND ZI BY THE C RESPECTIVE COMPONENTS OF BETA. C C ***** REFERENCES: C (1) EMAMI-NAEINI, A., AND P. VANDOOREN, COMPUTATION OF ZEROS C OF LINEAR MULTIVARIABLE SYSTEMS, REPORT NO. NA-80-03, C DEP'T. OF COMP. SCI., STANFORD UNIVERSITY, JULY, 1980. C C (2) LAUB, A.J., AND B.C. MOORE, CALCULATION OF TRANSMISSION C ZEROS USING QZ TECHNIQUES, AUTOMATICA, 14(1978),557-566. C C (3) MOLER, C.B., AND G.W. STEWART, AN ALGORITHM FOR GENER- C ALIZED MATRIX EIGENVALUE PROBLEMS, SIAM J. NUMER. ANAL., C 10(1973),241-256. C C ***** SOFTWARE HISTORY: C WRITTEN BY ALAN J. LAUB C DEP'T. OF ELEC. ENGRG. - SYSTEMS C UNIV. OF SOUTHERN CALIFORNIA C LOS ANGELES, CA 90089 C PH.: (213) 743-5535 C JANUARY, 1981. C MOST RECENT VERSION: AUG. 3, 1983. C C ------------------------------------------------------------------ C MATZ = .FALSE. EPS1 = -1.0D0 C C DETERMINE MACHINE PRECISION C EPS = 1.0D0 10 CONTINUE EPS = EPS/2.0D0 EPSP1 = EPS+1.0D0 IF (EPSP1 .GT. 1.0D0) GO TO 10 EPS = 2.0D0*EPS C C COMPUTE SQRT(EPS) * 1-NORM OF (A B) C (C D) C GNORM = 0.0D0 DO 40 J = 1,N T = 0.0D0 DO 20 I = 1,N T = T+DABS(A(I,J)) 20 CONTINUE DO 30 I = 1,L T = T+DABS(C(I,J)) 30 CONTINUE IF (T .GT. GNORM) GNORM = T 40 CONTINUE DO 70 J = 1,M T = 0.0D0 DO 50 I = 1,N T = T+DABS(B(I,J)) 50 CONTINUE DO 60 I = 1,L T = T+DABS(D(I,J)) 60 CONTINUE IF (T .GT. GNORM) GNORM = T 70 CONTINUE EPS = DSQRT(EPS)*GNORM C C GET REDUCED ORDER PENCIL WHICH WILL CONTAIN THE FINITE ZEROS C CALL ZEROS (A,B,C,D,M,NMAX,N,LMAX,L,MAX,EPS,H,G, + NZER,NRANK,SUM,DUMMY) C C APPLY QZ ALGORITHM TO THE FINITE PENCIL C IF (NZER .EQ. 0) RETURN IF (NZER .EQ. 1) GO TO 100 CALL QZHES (MAX,NZER,G,H,MATZ,G) CALL QZIT (MAX,NZER,G,H,EPS1,MATZ,G,IERR) IF (IERR .EQ. 0) GO TO 80 NZER = -IERR RETURN 80 CONTINUE CALL QZVAL (MAX,NZER,G,H,ZR,ZI,BETA,MATZ,G) C C COMPUTE THE MULTIVARIABLE ZEROS C DO 90 I = 1,NZER ZR(I) = ZR(I)/BETA(I) ZI(I) = ZI(I)/BETA(I) 90 CONTINUE RETURN 100 CONTINUE ZR(1) = G(1,1)/H(1,1) ZI(1) = 0.0D0 RETURN END