LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
cunt03.f
Go to the documentation of this file.
00001 *> \brief \b CUNT03
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *  Definition:
00009 *  ===========
00010 *
00011 *       SUBROUTINE CUNT03( RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK,
00012 *                          RWORK, RESULT, INFO )
00013 * 
00014 *       .. Scalar Arguments ..
00015 *       CHARACTER*( * )    RC
00016 *       INTEGER            INFO, K, LDU, LDV, LWORK, MU, MV, N
00017 *       REAL               RESULT
00018 *       ..
00019 *       .. Array Arguments ..
00020 *       REAL               RWORK( * )
00021 *       COMPLEX            U( LDU, * ), V( LDV, * ), WORK( * )
00022 *       ..
00023 *  
00024 *
00025 *> \par Purpose:
00026 *  =============
00027 *>
00028 *> \verbatim
00029 *>
00030 *> CUNT03 compares two unitary matrices U and V to see if their
00031 *> corresponding rows or columns span the same spaces.  The rows are
00032 *> checked if RC = 'R', and the columns are checked if RC = 'C'.
00033 *>
00034 *> RESULT is the maximum of
00035 *>
00036 *>    | V*V' - I | / ( MV ulp ), if RC = 'R', or
00037 *>
00038 *>    | V'*V - I | / ( MV ulp ), if RC = 'C',
00039 *>
00040 *> and the maximum over rows (or columns) 1 to K of
00041 *>
00042 *>    | U(i) - S*V(i) |/ ( N ulp )
00043 *>
00044 *> where abs(S) = 1 (chosen to minimize the expression), U(i) is the
00045 *> i-th row (column) of U, and V(i) is the i-th row (column) of V.
00046 *> \endverbatim
00047 *
00048 *  Arguments:
00049 *  ==========
00050 *
00051 *> \param[in] RC
00052 *> \verbatim
00053 *>          RC is CHARACTER*1
00054 *>          If RC = 'R' the rows of U and V are to be compared.
00055 *>          If RC = 'C' the columns of U and V are to be compared.
00056 *> \endverbatim
00057 *>
00058 *> \param[in] MU
00059 *> \verbatim
00060 *>          MU is INTEGER
00061 *>          The number of rows of U if RC = 'R', and the number of
00062 *>          columns if RC = 'C'.  If MU = 0 CUNT03 does nothing.
00063 *>          MU must be at least zero.
00064 *> \endverbatim
00065 *>
00066 *> \param[in] MV
00067 *> \verbatim
00068 *>          MV is INTEGER
00069 *>          The number of rows of V if RC = 'R', and the number of
00070 *>          columns if RC = 'C'.  If MV = 0 CUNT03 does nothing.
00071 *>          MV must be at least zero.
00072 *> \endverbatim
00073 *>
00074 *> \param[in] N
00075 *> \verbatim
00076 *>          N is INTEGER
00077 *>          If RC = 'R', the number of columns in the matrices U and V,
00078 *>          and if RC = 'C', the number of rows in U and V.  If N = 0
00079 *>          CUNT03 does nothing.  N must be at least zero.
00080 *> \endverbatim
00081 *>
00082 *> \param[in] K
00083 *> \verbatim
00084 *>          K is INTEGER
00085 *>          The number of rows or columns of U and V to compare.
00086 *>          0 <= K <= max(MU,MV).
00087 *> \endverbatim
00088 *>
00089 *> \param[in] U
00090 *> \verbatim
00091 *>          U is COMPLEX array, dimension (LDU,N)
00092 *>          The first matrix to compare.  If RC = 'R', U is MU by N, and
00093 *>          if RC = 'C', U is N by MU.
00094 *> \endverbatim
00095 *>
00096 *> \param[in] LDU
00097 *> \verbatim
00098 *>          LDU is INTEGER
00099 *>          The leading dimension of U.  If RC = 'R', LDU >= max(1,MU),
00100 *>          and if RC = 'C', LDU >= max(1,N).
00101 *> \endverbatim
00102 *>
00103 *> \param[in] V
00104 *> \verbatim
00105 *>          V is COMPLEX array, dimension (LDV,N)
00106 *>          The second matrix to compare.  If RC = 'R', V is MV by N, and
00107 *>          if RC = 'C', V is N by MV.
00108 *> \endverbatim
00109 *>
00110 *> \param[in] LDV
00111 *> \verbatim
00112 *>          LDV is INTEGER
00113 *>          The leading dimension of V.  If RC = 'R', LDV >= max(1,MV),
00114 *>          and if RC = 'C', LDV >= max(1,N).
00115 *> \endverbatim
00116 *>
00117 *> \param[out] WORK
00118 *> \verbatim
00119 *>          WORK is COMPLEX array, dimension (LWORK)
00120 *> \endverbatim
00121 *>
00122 *> \param[in] LWORK
00123 *> \verbatim
00124 *>          LWORK is INTEGER
00125 *>          The length of the array WORK.  For best performance, LWORK
00126 *>          should be at least N*N if RC = 'C' or M*M if RC = 'R', but
00127 *>          the tests will be done even if LWORK is 0.
00128 *> \endverbatim
00129 *>
00130 *> \param[out] RWORK
00131 *> \verbatim
00132 *>          RWORK is REAL array, dimension (max(MV,N))
00133 *> \endverbatim
00134 *>
00135 *> \param[out] RESULT
00136 *> \verbatim
00137 *>          RESULT is REAL
00138 *>          The value computed by the test described above.  RESULT is
00139 *>          limited to 1/ulp to avoid overflow.
00140 *> \endverbatim
00141 *>
00142 *> \param[out] INFO
00143 *> \verbatim
00144 *>          INFO is INTEGER
00145 *>          0  indicates a successful exit
00146 *>          -k indicates the k-th parameter had an illegal value
00147 *> \endverbatim
00148 *
00149 *  Authors:
00150 *  ========
00151 *
00152 *> \author Univ. of Tennessee 
00153 *> \author Univ. of California Berkeley 
00154 *> \author Univ. of Colorado Denver 
00155 *> \author NAG Ltd. 
00156 *
00157 *> \date November 2011
00158 *
00159 *> \ingroup complex_eig
00160 *
00161 *  =====================================================================
00162       SUBROUTINE CUNT03( RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK,
00163      $                   RWORK, RESULT, INFO )
00164 *
00165 *  -- LAPACK test routine (version 3.4.0) --
00166 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00167 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00168 *     November 2011
00169 *
00170 *     .. Scalar Arguments ..
00171       CHARACTER*( * )    RC
00172       INTEGER            INFO, K, LDU, LDV, LWORK, MU, MV, N
00173       REAL               RESULT
00174 *     ..
00175 *     .. Array Arguments ..
00176       REAL               RWORK( * )
00177       COMPLEX            U( LDU, * ), V( LDV, * ), WORK( * )
00178 *     ..
00179 *
00180 *  =====================================================================
00181 *
00182 *
00183 *     .. Parameters ..
00184       REAL               ZERO, ONE
00185       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
00186 *     ..
00187 *     .. Local Scalars ..
00188       INTEGER            I, IRC, J, LMX
00189       REAL               RES1, RES2, ULP
00190       COMPLEX            S, SU, SV
00191 *     ..
00192 *     .. External Functions ..
00193       LOGICAL            LSAME
00194       INTEGER            ICAMAX
00195       REAL               SLAMCH
00196       EXTERNAL           LSAME, ICAMAX, SLAMCH
00197 *     ..
00198 *     .. Intrinsic Functions ..
00199       INTRINSIC          ABS, CMPLX, MAX, MIN, REAL
00200 *     ..
00201 *     .. External Subroutines ..
00202       EXTERNAL           CUNT01, XERBLA
00203 *     ..
00204 *     .. Executable Statements ..
00205 *
00206 *     Check inputs
00207 *
00208       INFO = 0
00209       IF( LSAME( RC, 'R' ) ) THEN
00210          IRC = 0
00211       ELSE IF( LSAME( RC, 'C' ) ) THEN
00212          IRC = 1
00213       ELSE
00214          IRC = -1
00215       END IF
00216       IF( IRC.EQ.-1 ) THEN
00217          INFO = -1
00218       ELSE IF( MU.LT.0 ) THEN
00219          INFO = -2
00220       ELSE IF( MV.LT.0 ) THEN
00221          INFO = -3
00222       ELSE IF( N.LT.0 ) THEN
00223          INFO = -4
00224       ELSE IF( K.LT.0 .OR. K.GT.MAX( MU, MV ) ) THEN
00225          INFO = -5
00226       ELSE IF( ( IRC.EQ.0 .AND. LDU.LT.MAX( 1, MU ) ) .OR.
00227      $         ( IRC.EQ.1 .AND. LDU.LT.MAX( 1, N ) ) ) THEN
00228          INFO = -7
00229       ELSE IF( ( IRC.EQ.0 .AND. LDV.LT.MAX( 1, MV ) ) .OR.
00230      $         ( IRC.EQ.1 .AND. LDV.LT.MAX( 1, N ) ) ) THEN
00231          INFO = -9
00232       END IF
00233       IF( INFO.NE.0 ) THEN
00234          CALL XERBLA( 'CUNT03', -INFO )
00235          RETURN
00236       END IF
00237 *
00238 *     Initialize result
00239 *
00240       RESULT = ZERO
00241       IF( MU.EQ.0 .OR. MV.EQ.0 .OR. N.EQ.0 )
00242      $   RETURN
00243 *
00244 *     Machine constants
00245 *
00246       ULP = SLAMCH( 'Precision' )
00247 *
00248       IF( IRC.EQ.0 ) THEN
00249 *
00250 *        Compare rows
00251 *
00252          RES1 = ZERO
00253          DO 20 I = 1, K
00254             LMX = ICAMAX( N, U( I, 1 ), LDU )
00255             IF( V( I, LMX ).EQ.CMPLX( ZERO ) ) THEN
00256                SV = ONE
00257             ELSE
00258                SV = ABS( V( I, LMX ) ) / V( I, LMX )
00259             END IF
00260             IF( U( I, LMX ).EQ.CMPLX( ZERO ) ) THEN
00261                SU = ONE
00262             ELSE
00263                SU = ABS( U( I, LMX ) ) / U( I, LMX )
00264             END IF
00265             S = SV / SU
00266             DO 10 J = 1, N
00267                RES1 = MAX( RES1, ABS( U( I, J )-S*V( I, J ) ) )
00268    10       CONTINUE
00269    20    CONTINUE
00270          RES1 = RES1 / ( REAL( N )*ULP )
00271 *
00272 *        Compute orthogonality of rows of V.
00273 *
00274          CALL CUNT01( 'Rows', MV, N, V, LDV, WORK, LWORK, RWORK, RES2 )
00275 *
00276       ELSE
00277 *
00278 *        Compare columns
00279 *
00280          RES1 = ZERO
00281          DO 40 I = 1, K
00282             LMX = ICAMAX( N, U( 1, I ), 1 )
00283             IF( V( LMX, I ).EQ.CMPLX( ZERO ) ) THEN
00284                SV = ONE
00285             ELSE
00286                SV = ABS( V( LMX, I ) ) / V( LMX, I )
00287             END IF
00288             IF( U( LMX, I ).EQ.CMPLX( ZERO ) ) THEN
00289                SU = ONE
00290             ELSE
00291                SU = ABS( U( LMX, I ) ) / U( LMX, I )
00292             END IF
00293             S = SV / SU
00294             DO 30 J = 1, N
00295                RES1 = MAX( RES1, ABS( U( J, I )-S*V( J, I ) ) )
00296    30       CONTINUE
00297    40    CONTINUE
00298          RES1 = RES1 / ( REAL( N )*ULP )
00299 *
00300 *        Compute orthogonality of columns of V.
00301 *
00302          CALL CUNT01( 'Columns', N, MV, V, LDV, WORK, LWORK, RWORK,
00303      $                RES2 )
00304       END IF
00305 *
00306       RESULT = MIN( MAX( RES1, RES2 ), ONE / ULP )
00307       RETURN
00308 *
00309 *     End of CUNT03
00310 *
00311       END
 All Files Functions