LAPACK 3.3.1 Linear Algebra PACKage

# sbdt02.f

Go to the documentation of this file.
```00001       SUBROUTINE SBDT02( M, N, B, LDB, C, LDC, U, LDU, WORK, RESID )
00002 *
00003 *  -- LAPACK test routine (version 3.1) --
00004 *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
00005 *     November 2006
00006 *
00007 *     .. Scalar Arguments ..
00008       INTEGER            LDB, LDC, LDU, M, N
00009       REAL               RESID
00010 *     ..
00011 *     .. Array Arguments ..
00012       REAL               B( LDB, * ), C( LDC, * ), U( LDU, * ),
00013      \$                   WORK( * )
00014 *     ..
00015 *
00016 *  Purpose
00017 *  =======
00018 *
00019 *  SBDT02 tests the change of basis C = U' * B by computing the residual
00020 *
00021 *     RESID = norm( B - U * C ) / ( max(m,n) * norm(B) * EPS ),
00022 *
00023 *  where B and C are M by N matrices, U is an M by M orthogonal matrix,
00024 *  and EPS is the machine precision.
00025 *
00026 *  Arguments
00027 *  =========
00028 *
00029 *  M       (input) INTEGER
00030 *          The number of rows of the matrices B and C and the order of
00031 *          the matrix Q.
00032 *
00033 *  N       (input) INTEGER
00034 *          The number of columns of the matrices B and C.
00035 *
00036 *  B       (input) REAL array, dimension (LDB,N)
00037 *          The m by n matrix B.
00038 *
00039 *  LDB     (input) INTEGER
00040 *          The leading dimension of the array B.  LDB >= max(1,M).
00041 *
00042 *  C       (input) REAL array, dimension (LDC,N)
00043 *          The m by n matrix C, assumed to contain U' * B.
00044 *
00045 *  LDC     (input) INTEGER
00046 *          The leading dimension of the array C.  LDC >= max(1,M).
00047 *
00048 *  U       (input) REAL array, dimension (LDU,M)
00049 *          The m by m orthogonal matrix U.
00050 *
00051 *  LDU     (input) INTEGER
00052 *          The leading dimension of the array U.  LDU >= max(1,M).
00053 *
00054 *  WORK    (workspace) REAL array, dimension (M)
00055 *
00056 *  RESID   (output) REAL
00057 *          RESID = norm( B - U * C ) / ( max(m,n) * norm(B) * EPS ),
00058 *
00059 * ======================================================================
00060 *
00061 *     .. Parameters ..
00062       REAL               ZERO, ONE
00063       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
00064 *     ..
00065 *     .. Local Scalars ..
00066       INTEGER            J
00067       REAL               BNORM, EPS, REALMN
00068 *     ..
00069 *     .. External Functions ..
00070       REAL               SASUM, SLAMCH, SLANGE
00071       EXTERNAL           SASUM, SLAMCH, SLANGE
00072 *     ..
00073 *     .. External Subroutines ..
00074       EXTERNAL           SCOPY, SGEMV
00075 *     ..
00076 *     .. Intrinsic Functions ..
00077       INTRINSIC          MAX, MIN, REAL
00078 *     ..
00079 *     .. Executable Statements ..
00080 *
00081 *     Quick return if possible
00082 *
00083       RESID = ZERO
00084       IF( M.LE.0 .OR. N.LE.0 )
00085      \$   RETURN
00086       REALMN = REAL( MAX( M, N ) )
00087       EPS = SLAMCH( 'Precision' )
00088 *
00089 *     Compute norm( B - U * C )
00090 *
00091       DO 10 J = 1, N
00092          CALL SCOPY( M, B( 1, J ), 1, WORK, 1 )
00093          CALL SGEMV( 'No transpose', M, M, -ONE, U, LDU, C( 1, J ), 1,
00094      \$               ONE, WORK, 1 )
00095          RESID = MAX( RESID, SASUM( M, WORK, 1 ) )
00096    10 CONTINUE
00097 *
00098 *     Compute norm of B.
00099 *
00100       BNORM = SLANGE( '1', M, N, B, LDB, WORK )
00101 *
00102       IF( BNORM.LE.ZERO ) THEN
00103          IF( RESID.NE.ZERO )
00104      \$      RESID = ONE / EPS
00105       ELSE
00106          IF( BNORM.GE.RESID ) THEN
00107             RESID = ( RESID / BNORM ) / ( REALMN*EPS )
00108          ELSE
00109             IF( BNORM.LT.ONE ) THEN
00110                RESID = ( MIN( RESID, REALMN*BNORM ) / BNORM ) /
00111      \$                 ( REALMN*EPS )
00112             ELSE
00113                RESID = MIN( RESID / BNORM, REALMN ) / ( REALMN*EPS )
00114             END IF
00115          END IF
00116       END IF
00117       RETURN
00118 *
00119 *     End of SBDT02
00120 *
00121       END
```