1 SUBROUTINE DORT03( RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK,
2 $ RESULT, INFO )
3 *
4 * -- LAPACK test routine (version 3.1) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2006
7 *
8 * .. Scalar Arguments ..
9 CHARACTER*( * ) RC
10 INTEGER INFO, K, LDU, LDV, LWORK, MU, MV, N
11 DOUBLE PRECISION RESULT
12 * ..
13 * .. Array Arguments ..
14 DOUBLE PRECISION U( LDU, * ), V( LDV, * ), WORK( * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * DORT03 compares two orthogonal matrices U and V to see if their
21 * corresponding rows or columns span the same spaces. The rows are
22 * checked if RC = 'R', and the columns are checked if RC = 'C'.
23 *
24 * RESULT is the maximum of
25 *
26 * | V*V' - I | / ( MV ulp ), if RC = 'R', or
27 *
28 * | V'*V - I | / ( MV ulp ), if RC = 'C',
29 *
30 * and the maximum over rows (or columns) 1 to K of
31 *
32 * | U(i) - S*V(i) |/ ( N ulp )
33 *
34 * where S is +-1 (chosen to minimize the expression), U(i) is the i-th
35 * row (column) of U, and V(i) is the i-th row (column) of V.
36 *
37 * Arguments
38 * ==========
39 *
40 * RC (input) CHARACTER*1
41 * If RC = 'R' the rows of U and V are to be compared.
42 * If RC = 'C' the columns of U and V are to be compared.
43 *
44 * MU (input) INTEGER
45 * The number of rows of U if RC = 'R', and the number of
46 * columns if RC = 'C'. If MU = 0 DORT03 does nothing.
47 * MU must be at least zero.
48 *
49 * MV (input) INTEGER
50 * The number of rows of V if RC = 'R', and the number of
51 * columns if RC = 'C'. If MV = 0 DORT03 does nothing.
52 * MV must be at least zero.
53 *
54 * N (input) INTEGER
55 * If RC = 'R', the number of columns in the matrices U and V,
56 * and if RC = 'C', the number of rows in U and V. If N = 0
57 * DORT03 does nothing. N must be at least zero.
58 *
59 * K (input) INTEGER
60 * The number of rows or columns of U and V to compare.
61 * 0 <= K <= max(MU,MV).
62 *
63 * U (input) DOUBLE PRECISION array, dimension (LDU,N)
64 * The first matrix to compare. If RC = 'R', U is MU by N, and
65 * if RC = 'C', U is N by MU.
66 *
67 * LDU (input) INTEGER
68 * The leading dimension of U. If RC = 'R', LDU >= max(1,MU),
69 * and if RC = 'C', LDU >= max(1,N).
70 *
71 * V (input) DOUBLE PRECISION array, dimension (LDV,N)
72 * The second matrix to compare. If RC = 'R', V is MV by N, and
73 * if RC = 'C', V is N by MV.
74 *
75 * LDV (input) INTEGER
76 * The leading dimension of V. If RC = 'R', LDV >= max(1,MV),
77 * and if RC = 'C', LDV >= max(1,N).
78 *
79 * WORK (workspace) DOUBLE PRECISION array, dimension (LWORK)
80 *
81 * LWORK (input) INTEGER
82 * The length of the array WORK. For best performance, LWORK
83 * should be at least N*N if RC = 'C' or M*M if RC = 'R', but
84 * the tests will be done even if LWORK is 0.
85 *
86 * RESULT (output) DOUBLE PRECISION
87 * The value computed by the test described above. RESULT is
88 * limited to 1/ulp to avoid overflow.
89 *
90 * INFO (output) INTEGER
91 * 0 indicates a successful exit
92 * -k indicates the k-th parameter had an illegal value
93 *
94 * =====================================================================
95 *
96 * .. Parameters ..
97 DOUBLE PRECISION ZERO, ONE
98 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
99 * ..
100 * .. Local Scalars ..
101 INTEGER I, IRC, J, LMX
102 DOUBLE PRECISION RES1, RES2, S, ULP
103 * ..
104 * .. External Functions ..
105 LOGICAL LSAME
106 INTEGER IDAMAX
107 DOUBLE PRECISION DLAMCH
108 EXTERNAL LSAME, IDAMAX, DLAMCH
109 * ..
110 * .. Intrinsic Functions ..
111 INTRINSIC ABS, DBLE, MAX, MIN, SIGN
112 * ..
113 * .. External Subroutines ..
114 EXTERNAL DORT01, XERBLA
115 * ..
116 * .. Executable Statements ..
117 *
118 * Check inputs
119 *
120 INFO = 0
121 IF( LSAME( RC, 'R' ) ) THEN
122 IRC = 0
123 ELSE IF( LSAME( RC, 'C' ) ) THEN
124 IRC = 1
125 ELSE
126 IRC = -1
127 END IF
128 IF( IRC.EQ.-1 ) THEN
129 INFO = -1
130 ELSE IF( MU.LT.0 ) THEN
131 INFO = -2
132 ELSE IF( MV.LT.0 ) THEN
133 INFO = -3
134 ELSE IF( N.LT.0 ) THEN
135 INFO = -4
136 ELSE IF( K.LT.0 .OR. K.GT.MAX( MU, MV ) ) THEN
137 INFO = -5
138 ELSE IF( ( IRC.EQ.0 .AND. LDU.LT.MAX( 1, MU ) ) .OR.
139 $ ( IRC.EQ.1 .AND. LDU.LT.MAX( 1, N ) ) ) THEN
140 INFO = -7
141 ELSE IF( ( IRC.EQ.0 .AND. LDV.LT.MAX( 1, MV ) ) .OR.
142 $ ( IRC.EQ.1 .AND. LDV.LT.MAX( 1, N ) ) ) THEN
143 INFO = -9
144 END IF
145 IF( INFO.NE.0 ) THEN
146 CALL XERBLA( 'DORT03', -INFO )
147 RETURN
148 END IF
149 *
150 * Initialize result
151 *
152 RESULT = ZERO
153 IF( MU.EQ.0 .OR. MV.EQ.0 .OR. N.EQ.0 )
154 $ RETURN
155 *
156 * Machine constants
157 *
158 ULP = DLAMCH( 'Precision' )
159 *
160 IF( IRC.EQ.0 ) THEN
161 *
162 * Compare rows
163 *
164 RES1 = ZERO
165 DO 20 I = 1, K
166 LMX = IDAMAX( N, U( I, 1 ), LDU )
167 S = SIGN( ONE, U( I, LMX ) )*SIGN( ONE, V( I, LMX ) )
168 DO 10 J = 1, N
169 RES1 = MAX( RES1, ABS( U( I, J )-S*V( I, J ) ) )
170 10 CONTINUE
171 20 CONTINUE
172 RES1 = RES1 / ( DBLE( N )*ULP )
173 *
174 * Compute orthogonality of rows of V.
175 *
176 CALL DORT01( 'Rows', MV, N, V, LDV, WORK, LWORK, RES2 )
177 *
178 ELSE
179 *
180 * Compare columns
181 *
182 RES1 = ZERO
183 DO 40 I = 1, K
184 LMX = IDAMAX( N, U( 1, I ), 1 )
185 S = SIGN( ONE, U( LMX, I ) )*SIGN( ONE, V( LMX, I ) )
186 DO 30 J = 1, N
187 RES1 = MAX( RES1, ABS( U( J, I )-S*V( J, I ) ) )
188 30 CONTINUE
189 40 CONTINUE
190 RES1 = RES1 / ( DBLE( N )*ULP )
191 *
192 * Compute orthogonality of columns of V.
193 *
194 CALL DORT01( 'Columns', N, MV, V, LDV, WORK, LWORK, RES2 )
195 END IF
196 *
197 RESULT = MIN( MAX( RES1, RES2 ), ONE / ULP )
198 RETURN
199 *
200 * End of DORT03
201 *
202 END
2 $ RESULT, INFO )
3 *
4 * -- LAPACK test routine (version 3.1) --
5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
6 * November 2006
7 *
8 * .. Scalar Arguments ..
9 CHARACTER*( * ) RC
10 INTEGER INFO, K, LDU, LDV, LWORK, MU, MV, N
11 DOUBLE PRECISION RESULT
12 * ..
13 * .. Array Arguments ..
14 DOUBLE PRECISION U( LDU, * ), V( LDV, * ), WORK( * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * DORT03 compares two orthogonal matrices U and V to see if their
21 * corresponding rows or columns span the same spaces. The rows are
22 * checked if RC = 'R', and the columns are checked if RC = 'C'.
23 *
24 * RESULT is the maximum of
25 *
26 * | V*V' - I | / ( MV ulp ), if RC = 'R', or
27 *
28 * | V'*V - I | / ( MV ulp ), if RC = 'C',
29 *
30 * and the maximum over rows (or columns) 1 to K of
31 *
32 * | U(i) - S*V(i) |/ ( N ulp )
33 *
34 * where S is +-1 (chosen to minimize the expression), U(i) is the i-th
35 * row (column) of U, and V(i) is the i-th row (column) of V.
36 *
37 * Arguments
38 * ==========
39 *
40 * RC (input) CHARACTER*1
41 * If RC = 'R' the rows of U and V are to be compared.
42 * If RC = 'C' the columns of U and V are to be compared.
43 *
44 * MU (input) INTEGER
45 * The number of rows of U if RC = 'R', and the number of
46 * columns if RC = 'C'. If MU = 0 DORT03 does nothing.
47 * MU must be at least zero.
48 *
49 * MV (input) INTEGER
50 * The number of rows of V if RC = 'R', and the number of
51 * columns if RC = 'C'. If MV = 0 DORT03 does nothing.
52 * MV must be at least zero.
53 *
54 * N (input) INTEGER
55 * If RC = 'R', the number of columns in the matrices U and V,
56 * and if RC = 'C', the number of rows in U and V. If N = 0
57 * DORT03 does nothing. N must be at least zero.
58 *
59 * K (input) INTEGER
60 * The number of rows or columns of U and V to compare.
61 * 0 <= K <= max(MU,MV).
62 *
63 * U (input) DOUBLE PRECISION array, dimension (LDU,N)
64 * The first matrix to compare. If RC = 'R', U is MU by N, and
65 * if RC = 'C', U is N by MU.
66 *
67 * LDU (input) INTEGER
68 * The leading dimension of U. If RC = 'R', LDU >= max(1,MU),
69 * and if RC = 'C', LDU >= max(1,N).
70 *
71 * V (input) DOUBLE PRECISION array, dimension (LDV,N)
72 * The second matrix to compare. If RC = 'R', V is MV by N, and
73 * if RC = 'C', V is N by MV.
74 *
75 * LDV (input) INTEGER
76 * The leading dimension of V. If RC = 'R', LDV >= max(1,MV),
77 * and if RC = 'C', LDV >= max(1,N).
78 *
79 * WORK (workspace) DOUBLE PRECISION array, dimension (LWORK)
80 *
81 * LWORK (input) INTEGER
82 * The length of the array WORK. For best performance, LWORK
83 * should be at least N*N if RC = 'C' or M*M if RC = 'R', but
84 * the tests will be done even if LWORK is 0.
85 *
86 * RESULT (output) DOUBLE PRECISION
87 * The value computed by the test described above. RESULT is
88 * limited to 1/ulp to avoid overflow.
89 *
90 * INFO (output) INTEGER
91 * 0 indicates a successful exit
92 * -k indicates the k-th parameter had an illegal value
93 *
94 * =====================================================================
95 *
96 * .. Parameters ..
97 DOUBLE PRECISION ZERO, ONE
98 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
99 * ..
100 * .. Local Scalars ..
101 INTEGER I, IRC, J, LMX
102 DOUBLE PRECISION RES1, RES2, S, ULP
103 * ..
104 * .. External Functions ..
105 LOGICAL LSAME
106 INTEGER IDAMAX
107 DOUBLE PRECISION DLAMCH
108 EXTERNAL LSAME, IDAMAX, DLAMCH
109 * ..
110 * .. Intrinsic Functions ..
111 INTRINSIC ABS, DBLE, MAX, MIN, SIGN
112 * ..
113 * .. External Subroutines ..
114 EXTERNAL DORT01, XERBLA
115 * ..
116 * .. Executable Statements ..
117 *
118 * Check inputs
119 *
120 INFO = 0
121 IF( LSAME( RC, 'R' ) ) THEN
122 IRC = 0
123 ELSE IF( LSAME( RC, 'C' ) ) THEN
124 IRC = 1
125 ELSE
126 IRC = -1
127 END IF
128 IF( IRC.EQ.-1 ) THEN
129 INFO = -1
130 ELSE IF( MU.LT.0 ) THEN
131 INFO = -2
132 ELSE IF( MV.LT.0 ) THEN
133 INFO = -3
134 ELSE IF( N.LT.0 ) THEN
135 INFO = -4
136 ELSE IF( K.LT.0 .OR. K.GT.MAX( MU, MV ) ) THEN
137 INFO = -5
138 ELSE IF( ( IRC.EQ.0 .AND. LDU.LT.MAX( 1, MU ) ) .OR.
139 $ ( IRC.EQ.1 .AND. LDU.LT.MAX( 1, N ) ) ) THEN
140 INFO = -7
141 ELSE IF( ( IRC.EQ.0 .AND. LDV.LT.MAX( 1, MV ) ) .OR.
142 $ ( IRC.EQ.1 .AND. LDV.LT.MAX( 1, N ) ) ) THEN
143 INFO = -9
144 END IF
145 IF( INFO.NE.0 ) THEN
146 CALL XERBLA( 'DORT03', -INFO )
147 RETURN
148 END IF
149 *
150 * Initialize result
151 *
152 RESULT = ZERO
153 IF( MU.EQ.0 .OR. MV.EQ.0 .OR. N.EQ.0 )
154 $ RETURN
155 *
156 * Machine constants
157 *
158 ULP = DLAMCH( 'Precision' )
159 *
160 IF( IRC.EQ.0 ) THEN
161 *
162 * Compare rows
163 *
164 RES1 = ZERO
165 DO 20 I = 1, K
166 LMX = IDAMAX( N, U( I, 1 ), LDU )
167 S = SIGN( ONE, U( I, LMX ) )*SIGN( ONE, V( I, LMX ) )
168 DO 10 J = 1, N
169 RES1 = MAX( RES1, ABS( U( I, J )-S*V( I, J ) ) )
170 10 CONTINUE
171 20 CONTINUE
172 RES1 = RES1 / ( DBLE( N )*ULP )
173 *
174 * Compute orthogonality of rows of V.
175 *
176 CALL DORT01( 'Rows', MV, N, V, LDV, WORK, LWORK, RES2 )
177 *
178 ELSE
179 *
180 * Compare columns
181 *
182 RES1 = ZERO
183 DO 40 I = 1, K
184 LMX = IDAMAX( N, U( 1, I ), 1 )
185 S = SIGN( ONE, U( LMX, I ) )*SIGN( ONE, V( LMX, I ) )
186 DO 30 J = 1, N
187 RES1 = MAX( RES1, ABS( U( J, I )-S*V( J, I ) ) )
188 30 CONTINUE
189 40 CONTINUE
190 RES1 = RES1 / ( DBLE( N )*ULP )
191 *
192 * Compute orthogonality of columns of V.
193 *
194 CALL DORT01( 'Columns', N, MV, V, LDV, WORK, LWORK, RES2 )
195 END IF
196 *
197 RESULT = MIN( MAX( RES1, RES2 ), ONE / ULP )
198 RETURN
199 *
200 * End of DORT03
201 *
202 END