LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zunbdb5.f
Go to the documentation of this file.
1*> \brief \b ZUNBDB5
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download ZUNBDB5 + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zunbdb5.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zunbdb5.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zunbdb5.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE ZUNBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
22* LDQ2, WORK, LWORK, INFO )
23*
24* .. Scalar Arguments ..
25* INTEGER INCX1, INCX2, INFO, LDQ1, LDQ2, LWORK, M1, M2,
26* $ N
27* ..
28* .. Array Arguments ..
29* COMPLEX*16 Q1(LDQ1,*), Q2(LDQ2,*), WORK(*), X1(*), X2(*)
30* ..
31*
32*
33*> \par Purpose:
34* =============
35*>
36*>\verbatim
37*>
38*> ZUNBDB5 orthogonalizes the column vector
39*> X = [ X1 ]
40*> [ X2 ]
41*> with respect to the columns of
42*> Q = [ Q1 ] .
43*> [ Q2 ]
44*> The columns of Q must be orthonormal.
45*>
46*> If the projection is zero according to Kahan's "twice is enough"
47*> criterion, then some other vector from the orthogonal complement
48*> is returned. This vector is chosen in an arbitrary but deterministic
49*> way.
50*>
51*>\endverbatim
52*
53* Arguments:
54* ==========
55*
56*> \param[in] M1
57*> \verbatim
58*> M1 is INTEGER
59*> The dimension of X1 and the number of rows in Q1. 0 <= M1.
60*> \endverbatim
61*>
62*> \param[in] M2
63*> \verbatim
64*> M2 is INTEGER
65*> The dimension of X2 and the number of rows in Q2. 0 <= M2.
66*> \endverbatim
67*>
68*> \param[in] N
69*> \verbatim
70*> N is INTEGER
71*> The number of columns in Q1 and Q2. 0 <= N.
72*> \endverbatim
73*>
74*> \param[in,out] X1
75*> \verbatim
76*> X1 is COMPLEX*16 array, dimension (M1)
77*> On entry, the top part of the vector to be orthogonalized.
78*> On exit, the top part of the projected vector.
79*> \endverbatim
80*>
81*> \param[in] INCX1
82*> \verbatim
83*> INCX1 is INTEGER
84*> Increment for entries of X1.
85*> \endverbatim
86*>
87*> \param[in,out] X2
88*> \verbatim
89*> X2 is COMPLEX*16 array, dimension (M2)
90*> On entry, the bottom part of the vector to be
91*> orthogonalized. On exit, the bottom part of the projected
92*> vector.
93*> \endverbatim
94*>
95*> \param[in] INCX2
96*> \verbatim
97*> INCX2 is INTEGER
98*> Increment for entries of X2.
99*> \endverbatim
100*>
101*> \param[in] Q1
102*> \verbatim
103*> Q1 is COMPLEX*16 array, dimension (LDQ1, N)
104*> The top part of the orthonormal basis matrix.
105*> \endverbatim
106*>
107*> \param[in] LDQ1
108*> \verbatim
109*> LDQ1 is INTEGER
110*> The leading dimension of Q1. LDQ1 >= M1.
111*> \endverbatim
112*>
113*> \param[in] Q2
114*> \verbatim
115*> Q2 is COMPLEX*16 array, dimension (LDQ2, N)
116*> The bottom part of the orthonormal basis matrix.
117*> \endverbatim
118*>
119*> \param[in] LDQ2
120*> \verbatim
121*> LDQ2 is INTEGER
122*> The leading dimension of Q2. LDQ2 >= M2.
123*> \endverbatim
124*>
125*> \param[out] WORK
126*> \verbatim
127*> WORK is COMPLEX*16 array, dimension (LWORK)
128*> \endverbatim
129*>
130*> \param[in] LWORK
131*> \verbatim
132*> LWORK is INTEGER
133*> The dimension of the array WORK. LWORK >= N.
134*> \endverbatim
135*>
136*> \param[out] INFO
137*> \verbatim
138*> INFO is INTEGER
139*> = 0: successful exit.
140*> < 0: if INFO = -i, the i-th argument had an illegal value.
141*> \endverbatim
142*
143* Authors:
144* ========
145*
146*> \author Univ. of Tennessee
147*> \author Univ. of California Berkeley
148*> \author Univ. of Colorado Denver
149*> \author NAG Ltd.
150*
151*> \ingroup unbdb5
152*
153* =====================================================================
154 SUBROUTINE zunbdb5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
155 $ LDQ2, WORK, LWORK, INFO )
156*
157* -- LAPACK computational routine --
158* -- LAPACK is a software package provided by Univ. of Tennessee, --
159* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
160*
161* .. Scalar Arguments ..
162 INTEGER INCX1, INCX2, INFO, LDQ1, LDQ2, LWORK, M1, M2,
163 $ n
164* ..
165* .. Array Arguments ..
166 COMPLEX*16 Q1(LDQ1,*), Q2(LDQ2,*), WORK(*), X1(*), X2(*)
167* ..
168*
169* =====================================================================
170*
171* .. Parameters ..
172 DOUBLE PRECISION REALZERO
173 parameter( realzero = 0.0d0 )
174 COMPLEX*16 ONE, ZERO
175 parameter( one = (1.0d0,0.0d0), zero = (0.0d0,0.0d0) )
176* ..
177* .. Local Scalars ..
178 INTEGER CHILDINFO, I, J
179 DOUBLE PRECISION EPS, NORM, SCL, SSQ
180* ..
181* .. External Subroutines ..
182 EXTERNAL zlassq, zunbdb6, zscal, xerbla
183* ..
184* .. External Functions ..
185 DOUBLE PRECISION DLAMCH, DZNRM2
186 EXTERNAL dlamch, dznrm2
187* ..
188* .. Intrinsic Function ..
189 INTRINSIC max
190* ..
191* .. Executable Statements ..
192*
193* Test input arguments
194*
195 info = 0
196 IF( m1 .LT. 0 ) THEN
197 info = -1
198 ELSE IF( m2 .LT. 0 ) THEN
199 info = -2
200 ELSE IF( n .LT. 0 ) THEN
201 info = -3
202 ELSE IF( incx1 .LT. 1 ) THEN
203 info = -5
204 ELSE IF( incx2 .LT. 1 ) THEN
205 info = -7
206 ELSE IF( ldq1 .LT. max( 1, m1 ) ) THEN
207 info = -9
208 ELSE IF( ldq2 .LT. max( 1, m2 ) ) THEN
209 info = -11
210 ELSE IF( lwork .LT. n ) THEN
211 info = -13
212 END IF
213*
214 IF( info .NE. 0 ) THEN
215 CALL xerbla( 'ZUNBDB5', -info )
216 RETURN
217 END IF
218*
219 eps = dlamch( 'Precision' )
220*
221* Project X onto the orthogonal complement of Q if X is nonzero
222*
223 scl = realzero
224 ssq = realzero
225 CALL zlassq( m1, x1, incx1, scl, ssq )
226 CALL zlassq( m2, x2, incx2, scl, ssq )
227 norm = scl * sqrt( ssq )
228*
229 IF( norm .GT. n * eps ) THEN
230* Scale vector to unit norm to avoid problems in the caller code.
231* Computing the reciprocal is undesirable but
232* * xLASCL cannot be used because of the vector increments and
233* * the round-off error has a negligible impact on
234* orthogonalization.
235 CALL zscal( m1, one / norm, x1, incx1 )
236 CALL zscal( m2, one / norm, x2, incx2 )
237 CALL zunbdb6( m1, m2, n, x1, incx1, x2, incx2, q1, ldq1, q2,
238 $ ldq2, work, lwork, childinfo )
239*
240* If the projection is nonzero, then return
241*
242 IF( dznrm2(m1,x1,incx1) .NE. realzero
243 $ .OR. dznrm2(m2,x2,incx2) .NE. realzero ) THEN
244 RETURN
245 END IF
246 END IF
247*
248* Project each standard basis vector e_1,...,e_M1 in turn, stopping
249* when a nonzero projection is found
250*
251 DO i = 1, m1
252 DO j = 1, m1
253 x1(j) = zero
254 END DO
255 x1(i) = one
256 DO j = 1, m2
257 x2(j) = zero
258 END DO
259 CALL zunbdb6( m1, m2, n, x1, incx1, x2, incx2, q1, ldq1, q2,
260 $ ldq2, work, lwork, childinfo )
261 IF( dznrm2(m1,x1,incx1) .NE. realzero
262 $ .OR. dznrm2(m2,x2,incx2) .NE. realzero ) THEN
263 RETURN
264 END IF
265 END DO
266*
267* Project each standard basis vector e_(M1+1),...,e_(M1+M2) in turn,
268* stopping when a nonzero projection is found
269*
270 DO i = 1, m2
271 DO j = 1, m1
272 x1(j) = zero
273 END DO
274 DO j = 1, m2
275 x2(j) = zero
276 END DO
277 x2(i) = one
278 CALL zunbdb6( m1, m2, n, x1, incx1, x2, incx2, q1, ldq1, q2,
279 $ ldq2, work, lwork, childinfo )
280 IF( dznrm2(m1,x1,incx1) .NE. realzero
281 $ .OR. dznrm2(m2,x2,incx2) .NE. realzero ) THEN
282 RETURN
283 END IF
284 END DO
285*
286 RETURN
287*
288* End of ZUNBDB5
289*
290 END
291
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
zlassq
subroutine zlassq(n, x, incx, scale, sumsq)
ZLASSQ updates a sum of squares represented in scaled form.
Definition zlassq.f90:124
zscal
subroutine zscal(n, za, zx, incx)
ZSCAL
Definition zscal.f:78
zunbdb5
subroutine zunbdb5(m1, m2, n, x1, incx1, x2, incx2, q1, ldq1, q2, ldq2, work, lwork, info)
ZUNBDB5
Definition zunbdb5.f:156
zunbdb6
subroutine zunbdb6(m1, m2, n, x1, incx1, x2, incx2, q1, ldq1, q2, ldq2, work, lwork, info)
ZUNBDB6
Definition zunbdb6.f:159