Purpose
To preprocess the input-output data for estimating the matrices of a linear time-invariant dynamical system and to find an estimate of the system order. The input-output data can, optionally, be processed sequentially.Specification
SUBROUTINE IB01AD( METH, ALG, JOBD, BATCH, CONCT, CTRL, NOBR, M, $ L, NSMP, U, LDU, Y, LDY, N, R, LDR, SV, RCOND, $ TOL, IWORK, DWORK, LDWORK, IWARN, INFO ) C .. Scalar Arguments .. DOUBLE PRECISION RCOND, TOL INTEGER INFO, IWARN, L, LDR, LDU, LDWORK, LDY, M, N, $ NOBR, NSMP CHARACTER ALG, BATCH, CONCT, CTRL, JOBD, METH C .. Array Arguments .. INTEGER IWORK(*) DOUBLE PRECISION DWORK(*), R(LDR, *), SV(*), U(LDU, *), $ Y(LDY, *)Arguments
Mode Parameters
METH CHARACTER*1 Specifies the subspace identification method to be used, as follows: = 'M': MOESP algorithm with past inputs and outputs; = 'N': N4SID algorithm. ALG CHARACTER*1 Specifies the algorithm for computing the triangular factor R, as follows: = 'C': Cholesky algorithm applied to the correlation matrix of the input-output data; = 'F': Fast QR algorithm; = 'Q': QR algorithm applied to the concatenated block Hankel matrices. JOBD CHARACTER*1 Specifies whether or not the matrices B and D should later be computed using the MOESP approach, as follows: = 'M': the matrices B and D should later be computed using the MOESP approach; = 'N': the matrices B and D should not be computed using the MOESP approach. This parameter is not relevant for METH = 'N'. BATCH CHARACTER*1 Specifies whether or not sequential data processing is to be used, and, for sequential processing, whether or not the current data block is the first block, an intermediate block, or the last block, as follows: = 'F': the first block in sequential data processing; = 'I': an intermediate block in sequential data processing; = 'L': the last block in sequential data processing; = 'O': one block only (non-sequential data processing). NOTE that when 100 cycles of sequential data processing are completed for BATCH = 'I', a warning is issued, to prevent for an infinite loop. CONCT CHARACTER*1 Specifies whether or not the successive data blocks in sequential data processing belong to a single experiment, as follows: = 'C': the current data block is a continuation of the previous data block and/or it will be continued by the next data block; = 'N': there is no connection between the current data block and the previous and/or the next ones. This parameter is not used if BATCH = 'O'. CTRL CHARACTER*1 Specifies whether or not the user's confirmation of the system order estimate is desired, as follows: = 'C': user's confirmation; = 'N': no confirmation. If CTRL = 'C', a reverse communication routine, IB01OY, is indirectly called (by SLICOT Library routine IB01OD), and, after inspecting the singular values and system order estimate, n, the user may accept n or set a new value. IB01OY is not called if CTRL = 'N'.Input/Output Parameters
NOBR (input) INTEGER The number of block rows, s, in the input and output block Hankel matrices to be processed. NOBR > 0. (In the MOESP theory, NOBR should be larger than n, the estimated dimension of state vector.) M (input) INTEGER The number of system inputs. M >= 0. When M = 0, no system inputs are processed. L (input) INTEGER The number of system outputs. L > 0. NSMP (input) INTEGER The number of rows of matrices U and Y (number of samples, t). (When sequential data processing is used, NSMP is the number of samples of the current data block.) NSMP >= 2*(M+L+1)*NOBR - 1, for non-sequential processing; NSMP >= 2*NOBR, for sequential processing. The total number of samples when calling the routine with BATCH = 'L' should be at least 2*(M+L+1)*NOBR - 1. The NSMP argument may vary from a cycle to another in sequential data processing, but NOBR, M, and L should be kept constant. For efficiency, it is advisable to use NSMP as large as possible. U (input) DOUBLE PRECISION array, dimension (LDU,M) The leading NSMP-by-M part of this array must contain the t-by-m input-data sequence matrix U, U = [u_1 u_2 ... u_m]. Column j of U contains the NSMP values of the j-th input component for consecutive time increments. If M = 0, this array is not referenced. LDU INTEGER The leading dimension of the array U. LDU >= NSMP, if M > 0; LDU >= 1, if M = 0. Y (input) DOUBLE PRECISION array, dimension (LDY,L) The leading NSMP-by-L part of this array must contain the t-by-l output-data sequence matrix Y, Y = [y_1 y_2 ... y_l]. Column j of Y contains the NSMP values of the j-th output component for consecutive time increments. LDY INTEGER The leading dimension of the array Y. LDY >= NSMP. N (output) INTEGER The estimated order of the system. If CTRL = 'C', the estimated order has been reset to a value specified by the user. R (output or input/output) DOUBLE PRECISION array, dimension ( LDR,2*(M+L)*NOBR ) On exit, if ALG = 'C' and BATCH = 'F' or 'I', the leading 2*(M+L)*NOBR-by-2*(M+L)*NOBR upper triangular part of this array contains the current upper triangular part of the correlation matrix in sequential data processing. If ALG = 'F' and BATCH = 'F' or 'I', the array R is not referenced. On exit, if INFO = 0, ALG = 'Q', and BATCH = 'F' or 'I', the leading 2*(M+L)*NOBR-by-2*(M+L)*NOBR upper triangular part of this array contains the current upper triangular factor R from the QR factorization of the concatenated block Hankel matrices. Denote R_ij, i,j = 1:4, the ij submatrix of R, partitioned by M*NOBR, M*NOBR, L*NOBR, and L*NOBR rows and columns. On exit, if INFO = 0 and BATCH = 'L' or 'O', the leading 2*(M+L)*NOBR-by-2*(M+L)*NOBR upper triangular part of this array contains the matrix S, the processed upper triangular factor R from the QR factorization of the concatenated block Hankel matrices, as required by other subroutines. Specifically, let S_ij, i,j = 1:4, be the ij submatrix of S, partitioned by M*NOBR, L*NOBR, M*NOBR, and L*NOBR rows and columns. The submatrix S_22 contains the matrix of left singular vectors needed subsequently. Useful information is stored in S_11 and in the block-column S_14 : S_44. For METH = 'M' and JOBD = 'M', the upper triangular part of S_31 contains the upper triangular factor in the QR factorization of the matrix R_1c = [ R_12' R_22' R_11' ]', and S_12 contains the corresponding leading part of the transformed matrix R_2c = [ R_13' R_23' R_14' ]'. For METH = 'N', the subarray S_41 : S_43 contains the transpose of the matrix contained in S_14 : S_34. The details of the contents of R need not be known if this routine is followed by SLICOT Library routine IB01BD. On entry, if ALG = 'C', or ALG = 'Q', and BATCH = 'I' or 'L', the leading 2*(M+L)*NOBR-by-2*(M+L)*NOBR upper triangular part of this array must contain the upper triangular matrix R computed at the previous call of this routine in sequential data processing. The array R need not be set on entry if ALG = 'F' or if BATCH = 'F' or 'O'. LDR INTEGER The leading dimension of the array R. LDR >= MAX( 2*(M+L)*NOBR, 3*M*NOBR ), for METH = 'M' and JOBD = 'M'; LDR >= 2*(M+L)*NOBR, for METH = 'M' and JOBD = 'N' or for METH = 'N'. SV (output) DOUBLE PRECISION array, dimension ( L*NOBR ) The singular values used to estimate the system order.Tolerances
RCOND DOUBLE PRECISION The tolerance to be used for estimating the rank of matrices. If the user sets RCOND > 0, the given value of RCOND is used as a lower bound for the reciprocal condition number; an m-by-n matrix whose estimated condition number is less than 1/RCOND is considered to be of full rank. If the user sets RCOND <= 0, then an implicitly computed, default tolerance, defined by RCONDEF = m*n*EPS, is used instead, where EPS is the relative machine precision (see LAPACK Library routine DLAMCH). This parameter is not used for METH = 'M'. TOL DOUBLE PRECISION Absolute tolerance used for determining an estimate of the system order. If TOL >= 0, the estimate is indicated by the index of the last singular value greater than or equal to TOL. (Singular values less than TOL are considered as zero.) When TOL = 0, an internally computed default value, TOL = NOBR*EPS*SV(1), is used, where SV(1) is the maximal singular value, and EPS is the relative machine precision (see LAPACK Library routine DLAMCH). When TOL < 0, the estimate is indicated by the index of the singular value that has the largest logarithmic gap to its successor.Workspace
IWORK INTEGER array, dimension (LIWORK) LIWORK >= (M+L)*NOBR, if METH = 'N'; LIWORK >= M+L, if METH = 'M' and ALG = 'F'; LIWORK >= 0, if METH = 'M' and ALG = 'C' or 'Q'. DWORK DOUBLE PRECISION array, dimension (LDWORK) On exit, if INFO = 0, DWORK(1) returns the optimal value of LDWORK, and, for METH = 'N', and BATCH = 'L' or 'O', DWORK(2) and DWORK(3) contain the reciprocal condition numbers of the triangular factors of the matrices U_f and r_1 [6]. On exit, if INFO = -23, DWORK(1) returns the minimum value of LDWORK. Let k = 0, if CONCT = 'N' and ALG = 'C' or 'Q'; k = 2*NOBR-1, if CONCT = 'C' and ALG = 'C' or 'Q'; k = 2*NOBR*(M+L+1), if CONCT = 'N' and ALG = 'F'; k = 2*NOBR*(M+L+2), if CONCT = 'C' and ALG = 'F'. The first (M+L)*k elements of DWORK should be preserved during successive calls of the routine with BATCH = 'F' or 'I', till the final call with BATCH = 'L'. LDWORK INTEGER The length of the array DWORK. LDWORK >= (4*NOBR-2)*(M+L), if ALG = 'C', BATCH = 'F' or 'I' and CONCT = 'C'; LDWORK >= 1, if ALG = 'C', BATCH = 'F' or 'I' and CONCT = 'N'; LDWORK >= max((4*NOBR-2)*(M+L), 5*L*NOBR), if METH = 'M', ALG = 'C', BATCH = 'L' and CONCT = 'C'; LDWORK >= max((2*M-1)*NOBR, (M+L)*NOBR, 5*L*NOBR), if METH = 'M', JOBD = 'M', ALG = 'C', BATCH = 'O', or (BATCH = 'L' and CONCT = 'N'); LDWORK >= 5*L*NOBR, if METH = 'M', JOBD = 'N', ALG = 'C', BATCH = 'O', or (BATCH = 'L' and CONCT = 'N'); LDWORK >= 5*(M+L)*NOBR+1, if METH = 'N', ALG = 'C', and BATCH = 'L' or 'O'; LDWORK >= (M+L)*2*NOBR*(M+L+3), if ALG = 'F', BATCH <> 'O' and CONCT = 'C'; LDWORK >= (M+L)*2*NOBR*(M+L+1), if ALG = 'F', BATCH = 'F', 'I' and CONCT = 'N'; LDWORK >= (M+L)*4*NOBR*(M+L+1)+(M+L)*2*NOBR, if ALG = 'F', BATCH = 'L' and CONCT = 'N', or BATCH = 'O'; LDWORK >= 4*(M+L)*NOBR, if ALG = 'Q', BATCH = 'F', and LDR >= NS = NSMP - 2*NOBR + 1; LDWORK >= max(4*(M+L)*NOBR, 5*L*NOBR), if METH = 'M', ALG = 'Q', BATCH = 'O', and LDR >= NS; LDWORK >= 5*(M+L)*NOBR+1, if METH = 'N', ALG = 'Q', BATCH = 'O', and LDR >= NS; LDWORK >= 6*(M+L)*NOBR, if ALG = 'Q', (BATCH = 'F' or 'O', and LDR < NS), or (BATCH = 'I' or 'L' and CONCT = 'N'); LDWORK >= 4*(NOBR+1)*(M+L)*NOBR, if ALG = 'Q', BATCH = 'I' or 'L' and CONCT = 'C'. The workspace used for ALG = 'Q' is LDRWRK*2*(M+L)*NOBR + 4*(M+L)*NOBR, where LDRWRK = LDWORK/(2*(M+L)*NOBR) - 2; recommended value LDRWRK = NS, assuming a large enough cache size. For good performance, LDWORK should be larger.Warning Indicator
IWARN INTEGER = 0: no warning; = 1: the number of 100 cycles in sequential data processing has been exhausted without signaling that the last block of data was get; the cycle counter was reinitialized; = 2: a fast algorithm was requested (ALG = 'C' or 'F'), but it failed, and the QR algorithm was then used (non-sequential data processing); = 3: all singular values were exactly zero, hence N = 0 (both input and output were identically zero); = 4: the least squares problems with coefficient matrix U_f, used for computing the weighted oblique projection (for METH = 'N'), have a rank-deficient coefficient matrix; = 5: the least squares problem with coefficient matrix r_1 [6], used for computing the weighted oblique projection (for METH = 'N'), has a rank-deficient coefficient matrix. NOTE: the values 4 and 5 of IWARN have no significance for the identification problem.Error Indicator
INFO INTEGER = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; = 1: a fast algorithm was requested (ALG = 'C', or 'F') in sequential data processing, but it failed; the routine can be repeatedly called again using the standard QR algorithm; = 2: the singular value decomposition (SVD) algorithm did not converge.Method
The procedure consists in three main steps, the first step being performed by one of the three algorithms included. 1.a) For non-sequential data processing using QR algorithm, a t x 2(m+l)s matrix H is constructed, where H = [ Uf' Up' Y' ], for METH = 'M', s+1,2s,t 1,s,t 1,2s,t H = [ U' Y' ], for METH = 'N', 1,2s,t 1,2s,t and Up , Uf , U , and Y are block Hankel 1,s,t s+1,2s,t 1,2s,t 1,2s,t matrices defined in terms of the input and output data [3]. A QR factorization is used to compress the data. The fast QR algorithm uses a QR factorization which exploits the block-Hankel structure. Actually, the Cholesky factor of H'*H is computed. 1.b) For sequential data processing using QR algorithm, the QR decomposition is done sequentially, by updating the upper triangular factor R. This is also performed internally if the workspace is not large enough to accommodate an entire batch. 1.c) For non-sequential or sequential data processing using Cholesky algorithm, the correlation matrix of input-output data is computed (sequentially, if requested), taking advantage of the block Hankel structure [7]. Then, the Cholesky factor of the correlation matrix is found, if possible. 2) A singular value decomposition (SVD) of a certain matrix is then computed, which reveals the order n of the system as the number of "non-zero" singular values. For the MOESP approach, this matrix is [ R_24' R_34' ]' := R(ms+1:(2m+l)s,(2m+l)s+1:2(m+l)s), where R is the upper triangular factor R constructed by SLICOT Library routine IB01MD. For the N4SID approach, a weighted oblique projection is computed from the upper triangular factor R and its SVD is then found. 3) The singular values are compared to the given, or default TOL, and the estimated order n is returned, possibly after user's confirmation.References
[1] Verhaegen M., and Dewilde, P. Subspace Model Identification. Part 1: The output-error state-space model identification class of algorithms. Int. J. Control, 56, pp. 1187-1210, 1992. [2] Verhaegen M. Subspace Model Identification. Part 3: Analysis of the ordinary output-error state-space model identification algorithm. Int. J. Control, 58, pp. 555-586, 1993. [3] Verhaegen M. Identification of the deterministic part of MIMO state space models given in innovations form from input-output data. Automatica, Vol.30, No.1, pp.61-74, 1994. [4] Van Overschee, P., and De Moor, B. N4SID: Subspace Algorithms for the Identification of Combined Deterministic-Stochastic Systems. Automatica, Vol.30, No.1, pp. 75-93, 1994. [5] Peternell, K., Scherrer, W. and Deistler, M. Statistical Analysis of Novel Subspace Identification Methods. Signal Processing, 52, pp. 161-177, 1996. [6] Sima, V. Subspace-based Algorithms for Multivariable System Identification. Studies in Informatics and Control, 5, pp. 335-344, 1996. [7] Sima, V. Cholesky or QR Factorization for Data Compression in Subspace-based Identification ? Proceedings of the Second NICONET Workshop on ``Numerical Control Software: SLICOT, a Useful Tool in Industry'', December 3, 1999, INRIA Rocquencourt, France, pp. 75-80, 1999.Numerical Aspects
The implemented method is numerically stable (when QR algorithm is used), reliable and efficient. The fast Cholesky or QR algorithms are more efficient, but the accuracy could diminish by forming the correlation matrix. The most time-consuming computational step is step 1: 2 The QR algorithm needs 0(t(2(m+l)s) ) floating point operations. 2 3 The Cholesky algorithm needs 0(2t(m+l) s)+0((2(m+l)s) ) floating point operations. 2 3 2 The fast QR algorithm needs 0(2t(m+l) s)+0(4(m+l) s ) floating point operations. 3 Step 2 of the algorithm requires 0(((m+l)s) ) floating point operations.Further Comments
For ALG = 'Q', BATCH = 'O' and LDR < NS, or BATCH <> 'O', the calculations could be rather inefficient if only minimal workspace (see argument LDWORK) is provided. It is advisable to provide as much workspace as possible. Almost optimal efficiency can be obtained for LDWORK = (NS+2)*(2*(M+L)*NOBR), assuming that the cache size is large enough to accommodate R, U, Y, and DWORK.Example
Program Text
* IB01AD EXAMPLE PROGRAM TEXT * Copyright (c) 2002-2017 NICONET e.V. * * .. Parameters .. INTEGER NIN, NOUT PARAMETER ( NIN = 5, NOUT = 6 ) INTEGER LDR, LDU, LDWORK, LDY, LIWORK, LMAX, MMAX, $ NOBRMX, NSMPMX PARAMETER ( LMAX = 5, MMAX = 5, NOBRMX = 20, NSMPMX = 2000, $ LDR = MAX( 2*( MMAX + LMAX )*NOBRMX, $ 3*MMAX*NOBRMX ), LDU = NSMPMX, $ LDWORK = MAX( 6*( MMAX + LMAX )*NOBRMX, $ ( MMAX + LMAX )*( 4*NOBRMX* $ ( MMAX + LMAX + 1 ) + 2*NOBRMX ), $ ( MMAX + LMAX )*4*NOBRMX* $ ( NOBRMX + 1 ) ), $ LDY = NSMPMX, LIWORK = ( MMAX + LMAX )*NOBRMX ) * .. Local Scalars .. LOGICAL NGIVEN CHARACTER ALG, BATCH, CONCT, CTRL, JOBD, METH INTEGER I, ICYCLE, II, INFO, IWARN, J, L, M, N, NCYCLE, $ NGIV, NOBR, NSAMPL, NSMP DOUBLE PRECISION RCOND, TOL * .. Local Arrays .. DOUBLE PRECISION DWORK(LDWORK), R(LDR, 2*(MMAX+LMAX)*NOBRMX), $ SV(LMAX*NOBRMX), U(LDU, MMAX), Y(LDY, LMAX) INTEGER IWORK(LIWORK) * .. External Functions .. LOGICAL LSAME EXTERNAL LSAME * .. External Subroutines .. EXTERNAL IB01AD * .. Intrinsic Functions .. INTRINSIC MAX * .. Executable Statements .. * WRITE ( NOUT, FMT = 99999 ) * Skip the heading in the data file and read the data. * If the value of N is positive, it will be taken as system order. READ ( NIN, FMT = '()' ) READ ( NIN, FMT = * ) NOBR, N, M, L, NSMP, RCOND, TOL, METH, ALG, $ JOBD, BATCH, CONCT, CTRL IF ( LSAME( BATCH, 'F' ) ) THEN READ ( NIN, FMT = * ) NCYCLE ELSE NCYCLE = 1 END IF NSAMPL = NCYCLE*NSMP * NGIVEN = N.GT.0 IF( NGIVEN ) $ NGIV = N IF ( NOBR.LE.0 .OR. NOBR.GT.NOBRMX ) THEN WRITE ( NOUT, FMT = 99997 ) NOBR ELSE IF ( M.LT.0 .OR. M.GT.MMAX ) THEN WRITE ( NOUT, FMT = 99996 ) M ELSE IF ( L.LE.0 .OR. L.GT.LMAX ) THEN WRITE ( NOUT, FMT = 99995 ) L ELSE IF ( NSMP.LT.0 .OR. NSMP.GT.NSMPMX .OR. $ ( NSMP.LT.2*( M + L + 1 )*NOBR - 1 .AND. $ LSAME( BATCH, 'O' ) ) .OR. $ ( NSAMPL.LT.2*( M + L + 1 )*NOBR - 1 .AND. $ LSAME( BATCH, 'L' ) ) .OR. $ NSMP.LT.2*NOBR .AND. ( LSAME( BATCH, 'F' ) .OR. $ LSAME( BATCH, 'I' ) ) ) THEN WRITE ( NOUT, FMT = 99994 ) NSMP ELSE IF ( NCYCLE.LE.0 .OR. NSAMPL.GT.NSMPMX ) THEN WRITE ( NOUT, FMT = 99993 ) NCYCLE ELSE * Read the matrices U and Y from the input file. IF ( M.GT.0 ) $ READ ( NIN, FMT = * ) $ ( ( U(I,J), J = 1, M ), I = 1, NSAMPL ) READ ( NIN, FMT = * ) ( ( Y(I,J), J = 1, L ), I = 1, NSAMPL ) * Compute the R factor from a QR (or Cholesky) factorization * of the Hankel-like matrix (or correlation matrix). DO 10 ICYCLE = 1, NCYCLE II = ( ICYCLE - 1 )*NSMP + 1 IF ( NCYCLE.GT.1 ) THEN IF ( ICYCLE.GT.1 ) BATCH = 'I' IF ( ICYCLE.EQ.NCYCLE ) BATCH = 'L' END IF CALL IB01AD( METH, ALG, JOBD, BATCH, CONCT, CTRL, NOBR, M, $ L, NSMP, U(II,1), LDU, Y(II,1), LDY, N, R, LDR, $ SV, RCOND, TOL, IWORK, DWORK, LDWORK, IWARN, $ INFO ) 10 CONTINUE IF ( INFO.NE.0 ) THEN WRITE ( NOUT, FMT = 99998 ) INFO ELSE IF ( IWARN.NE.0 ) $ WRITE ( NOUT, FMT = 99990 ) IWARN IF( NGIVEN ) $ N = NGIV WRITE ( NOUT, FMT = 99992 ) N WRITE ( NOUT, FMT = 99991 ) ( SV(I), I = 1,L*NOBR ) END IF END IF STOP 99999 FORMAT ( ' IB01AD EXAMPLE PROGRAM RESULTS', /1X) 99998 FORMAT ( ' INFO on exit from IB01AD = ',I2) 99997 FORMAT (/' NOBR is out of range.',/' NOBR = ', I5) 99996 FORMAT (/' M is out of range.',/' M = ', I5) 99995 FORMAT (/' L is out of range.',/' L = ', I5) 99994 FORMAT (/' NSMP is out of range.',/' NSMP = ', I5) 99993 FORMAT (/' NCYCLE is out of range.',/' NCYCLE = ', I5) 99992 FORMAT ( ' The order of the system is ', I5) 99991 FORMAT ( ' The singular values are ',/ (8(1X,F8.4))) 99990 FORMAT ( ' IWARN on exit from IB01AD = ',I2) ENDProgram Data
IB01AD EXAMPLE PROGRAM DATA 15 0 1 1 1000 0.0 -1.0 M C N O N N 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 6.41 3.41 3.41 3.41 6.41 3.41 3.41 3.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 3.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 6.41 3.41 3.41 3.41 3.41 3.41 3.41 3.41 6.41 6.41 6.41 6.41 6.41 6.41 4.766099 4.763659 4.839359 5.002979 5.017629 5.056699 5.154379 5.361949 5.425439 5.569519 5.681849 5.742899 5.803949 5.918729 5.821049 5.447419 5.061589 4.629349 4.267939 4.011519 3.850349 3.711159 3.569519 3.518239 3.652549 3.818609 3.862559 4.011519 4.353409 4.705049 5.083559 5.344859 5.274039 5.127519 4.761219 4.451089 4.221539 4.045709 3.874769 3.730689 3.662319 3.576849 3.542659 3.479169 3.454749 3.359509 3.298459 3.225199 3.200779 3.225199 3.227639 3.274039 3.457189 3.867449 4.321659 4.492599 4.431549 4.243519 4.050599 3.857679 3.730689 3.791739 3.921169 3.955359 3.847909 3.725809 3.611039 3.716039 4.092109 4.480389 4.814939 5.054259 5.303339 5.486489 5.672089 5.779529 5.799069 5.664759 5.291129 4.880879 4.558529 4.184909 3.889419 3.708719 3.623249 3.569519 3.718479 4.033499 4.412009 4.629349 4.558529 4.394919 4.180019 4.197119 4.431549 4.714819 4.961459 5.300899 5.567079 5.681849 5.545099 5.188569 4.883319 4.600049 4.270379 4.038389 3.838139 3.711159 3.591499 3.535329 3.486489 3.476729 3.425439 3.381489 3.369279 3.364389 3.347299 3.381489 3.420559 3.413229 3.452309 3.635459 4.038389 4.375379 4.727029 5.056699 5.298459 5.532889 5.466959 5.195899 4.885759 4.763659 4.875989 5.042049 5.283809 5.491379 5.596379 5.672089 5.772209 5.830819 5.933379 5.899189 5.935819 5.894309 5.918729 5.994429 5.957799 6.031059 6.062809 6.040829 6.096999 6.123859 6.162929 6.040829 5.845469 5.772209 5.799069 5.923609 5.928499 6.001759 6.001759 6.060369 5.882099 5.510909 5.322879 5.371719 5.454749 5.437649 5.159269 4.902859 4.587839 4.502369 4.595159 4.824709 5.064029 5.271599 5.466959 5.615919 5.528009 5.254499 4.883319 4.517019 4.197119 4.001759 3.806399 3.904079 3.923609 3.869889 3.806399 3.720929 3.818609 4.140949 4.529229 4.805179 5.086009 5.339969 5.532889 5.576849 5.667199 5.791739 5.850349 5.923609 5.921169 5.977339 5.740459 5.388809 5.000539 4.849129 4.944369 5.173919 5.369279 5.447419 5.603709 5.730689 5.850349 5.979779 5.991989 6.084789 5.940709 5.803949 5.791739 5.603709 5.264269 4.946809 4.619579 4.514579 4.433989 4.285029 4.121419 3.945589 3.984659 4.219099 4.546319 4.873549 5.154379 5.388809 5.613479 5.835699 5.884539 5.955359 5.762439 5.459629 5.061589 4.707499 4.458409 4.267939 4.053039 3.943149 3.825929 3.967569 4.280149 4.480389 4.492599 4.390039 4.197119 4.111649 3.982219 3.867449 3.767319 3.872329 4.236189 4.663539 4.971229 5.066469 4.902859 4.675749 4.392479 4.099439 4.114089 4.326539 4.643999 4.971229 5.159269 5.388809 5.576849 5.652549 5.803949 5.913839 5.886979 5.799069 5.730689 5.762439 5.813719 5.821049 5.928499 6.013969 5.764879 5.413229 5.098219 4.678189 4.372939 4.392479 4.590279 4.919949 5.017629 4.858899 4.675749 4.619579 4.834479 5.090889 5.376599 5.681849 5.823489 5.952919 6.062809 6.089669 6.075019 6.026179 5.994429 6.077459 5.857679 5.701389 5.730689 5.784419 5.823489 5.894309 5.762439 5.415679 4.961459 4.595159 4.331429 4.297239 4.582949 4.861339 5.173919 5.166589 4.919949 4.607369 4.370499 4.182469 4.038389 4.145839 4.431549 4.556089 4.480389 4.375379 4.370499 4.558529 4.858899 4.895529 4.741679 4.744129 4.875989 5.105539 5.239849 5.518239 5.652549 5.723369 5.855239 5.962679 5.984659 5.984659 6.055479 6.062809 6.055479 6.070129 5.784419 5.440099 5.056699 4.941929 5.010299 5.134849 5.313109 5.479169 5.623249 5.562199 5.330209 5.010299 4.665979 4.414459 4.201999 4.048159 4.079899 4.189789 4.131179 4.004199 3.916289 3.960239 4.199559 4.624469 4.883319 5.137289 5.379049 5.623249 5.762439 5.833259 5.686739 5.366839 5.225199 5.239849 5.354629 5.508469 5.596379 5.752669 5.874769 5.906519 5.894309 5.742899 5.447419 5.024959 4.883319 4.885759 4.893089 4.714819 4.451089 4.233749 4.043269 3.864999 3.757559 3.669639 3.593939 3.547539 3.506029 3.454749 3.398579 3.361949 3.339969 3.374159 3.520679 3.713599 3.757559 3.779529 3.696509 3.777089 3.886979 3.904079 3.850349 3.965129 4.282589 4.521899 4.714819 4.971229 5.220319 5.532889 5.652549 5.781979 5.955359 6.035939 6.118969 6.133629 6.153159 6.192229 6.143389 6.167809 5.991989 5.652549 5.459629 5.437649 5.339969 5.098219 4.785639 4.492599 4.236189 4.067689 3.933379 3.823489 3.730689 3.611039 3.564639 3.549989 3.557309 3.513359 3.515799 3.694059 4.072579 4.480389 4.705049 4.612259 4.385149 4.201999 4.026179 3.904079 3.774649 3.691619 3.845469 4.201999 4.585399 4.902859 5.256949 5.510909 5.640339 5.843029 5.974889 5.935819 5.821049 5.528009 5.171479 4.810059 4.453529 4.380269 4.565859 4.805179 5.125079 5.354629 5.589059 5.764879 5.923609 5.940709 5.857679 5.694059 5.486489 5.149499 4.844249 4.541439 4.267939 4.060369 3.960239 3.789299 3.642779 3.525569 3.498699 3.454749 3.408349 3.379049 3.376599 3.361949 3.359509 3.369279 3.398579 3.579289 3.948029 4.412009 4.585399 4.514579 4.343639 4.155599 3.984659 4.043269 4.307009 4.421779 4.353409 4.223979 4.053039 3.940709 3.838139 3.730689 3.652549 3.611039 3.564639 3.496259 3.462069 3.454749 3.425439 3.379049 3.432769 3.623249 3.974889 4.380269 4.714819 5.073799 5.369279 5.603709 5.745349 5.652549 5.401019 5.015189 4.709939 4.416899 4.236189 4.236189 4.248399 4.221539 4.297239 4.590279 4.893089 5.134849 5.427889 5.379049 5.364389 5.452309 5.567079 5.672089 5.769769 5.830819 5.923609 5.965129 6.057919 6.050599 6.072579 6.111649 6.070129 5.896749 5.755109 5.718479 5.821049 6.001759 6.001759 5.901629 5.557309 5.173919 4.800289 4.431549 4.194679 4.006639 3.850349 3.747789 3.642779 3.591499 3.569519 3.528009 3.537779 3.554869 3.493819 3.447419 3.440099 3.408349 3.410789 3.452309 3.681849 4.060369 4.441319 4.854019 5.154379 5.425439 5.596379 5.586619 5.354629 5.027399 4.863779 4.761219 4.570739 4.368059 4.397359 4.573189 4.841809 5.203219 5.452309 5.652549 5.855239 5.906519 5.952919 5.828369 5.791739 5.799069 5.813719 5.877209 5.955359 5.781979 5.518239 5.127519 4.763659 4.492599 4.233749 4.011519 3.855239 3.691619 3.635459 3.818609 4.155599 4.590279 4.988329 5.076239 4.907739 4.648889 4.377829 4.216649 4.287469 4.590279 4.846689 5.139729 5.388809 5.689179 5.884539 6.043269 6.170259 6.211769 6.250839 6.209329 6.013969 5.701389 5.469399 5.479169 5.557309 5.728249 5.882099 5.984659 5.901629 5.581729 5.371719 5.418119 5.510909 5.667199 5.791739 5.698949 5.484049 5.154379 4.980999 5.061589 5.195899 5.359509 5.615919 5.762439 5.857679 5.948029 5.835699 5.706269 5.498699 5.188569 5.117749 5.191009 5.315549 5.532889 5.444979 5.396139 5.274039 5.027399 4.744129 4.668419 4.651329 4.514579 4.267939 4.260609 4.263049 4.189789 4.277699 4.600049 4.932159 5.283809 5.528009 5.740459 5.874769 5.955359 5.991989 5.845469 5.528009 5.061589 4.734359 4.534109 4.534109 4.697729 4.744129 4.619579 4.643999 4.832039 5.132399 5.410789 5.625689 5.603709 5.315549 4.961459 4.619579 4.358289 4.155599 4.033499 3.886979 3.772209 3.640339 3.532889 3.435209 3.427889 3.422999 3.398579 3.603709 4.023729 4.451089 4.792969 4.902859 4.780759 4.590279 4.336309 4.145839 4.216649 4.433989 4.714819 5.098219 5.359509 5.569519 5.772209 5.921169 6.055479 5.962679 5.642779 5.435209 5.388809 5.537779 5.681849 5.701389 5.615919 5.667199 5.740459 5.803949 5.882099 5.950469 6.072579 6.148279 6.116529 6.177579 6.201999 6.206889 5.991989 5.564639 5.178799 4.998089 5.051819 5.232529 5.484049 5.686739 5.899189 5.869889 5.977339 6.053039 6.079899 6.128739 6.079899 6.167809 6.194679 6.236189 6.053039 5.652549 5.274039 4.858899 4.534109 4.455969 4.619579 4.866229 5.117749 5.166589 5.056699 5.002979 5.098219 5.325319 5.567079 5.466959 5.252059 4.946809 4.880879 4.980999 5.225199 5.459629 5.723369 5.791739 5.906519 5.991989 5.835699 5.528009 5.142169 4.775869 4.490159 4.236189 4.023729 3.886979 3.752669 3.681849 3.806399 4.145839 4.600049 5.002979 5.303339 5.552429 5.615919 5.523119 5.611039 5.713599 5.845469 5.899189 5.994429 6.092109 6.092109 6.143389 6.153159 6.233749 6.187349 6.013969 5.835699 5.774649 5.686739 5.537779 5.327759 5.054259 4.700169 4.394919 4.180019 4.043269 3.877209 3.752669 3.728249 3.869889 4.206889 4.355849 4.426669 4.453529 4.521899 4.392479 4.155599 3.965129 3.877209 3.970009 4.258169 4.421779 4.336309 4.299679 4.392479 4.675749 4.761219 4.658659 4.490159 4.307009 4.126299 3.972449 4.077459 4.372939 4.741679 5.088449 5.186129 5.037169 4.785639 4.563419 4.534109 4.705049 4.741679 4.648889 4.431549 4.238629 4.065249 3.943149 3.811279 3.691619 3.652549 3.825929 4.223979 4.424219 4.429109 4.319219 4.138509 3.965129 3.886979 3.801509 3.701389 3.640339 3.767319 4.150719 4.648889 4.990769 5.088449 5.022509 4.783199 4.685519 4.665979 4.707499 4.912619 5.195899 5.415679 5.623249 5.740459 5.899189 5.928499 6.050599 6.153159 5.965129 5.586619 5.381489 5.371719 5.486489 5.567079 5.821049 5.913839 5.994429 6.011519 5.999309 6.018849 5.821049 5.728249 5.740459 5.764879 5.882099 5.926049 5.750229 5.415679 4.995649 4.861339 4.902859 5.103099 5.364389 5.596379 5.752669 5.845469 5.928499 6.006639 5.840579 5.518239 5.173919 4.739239 4.458409 4.426669 4.602489 4.822269 5.183689 5.430329 5.652549 5.821049 5.706269 5.369279 5.027399 4.705049 4.414459 4.145839 3.965129 4.033499 4.372939 4.683079Program Results
IB01AD EXAMPLE PROGRAM RESULTS The order of the system is 4 The singular values are 69.8841 14.9963 3.6675 1.9677 0.3000 0.2078 0.1651 0.1373 0.1133 0.1059 0.0856 0.0784 0.0733 0.0678 0.0571