LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ chetrs_rook()

 subroutine chetrs_rook ( character uplo, integer n, integer nrhs, complex, dimension( lda, * ) a, integer lda, integer, dimension( * ) ipiv, complex, dimension( ldb, * ) b, integer ldb, integer info )

CHETRS_ROOK computes the solution to a system of linear equations A * X = B for HE matrices using factorization obtained with one of the bounded diagonal pivoting methods (max 2 interchanges)

Purpose:
CHETRS_ROOK solves a system of linear equations A*X = B with a complex
Hermitian matrix A using the factorization A = U*D*U**H or
A = L*D*L**H computed by CHETRF_ROOK.
Parameters
 [in] UPLO UPLO is CHARACTER*1 Specifies whether the details of the factorization are stored as an upper or lower triangular matrix. = 'U': Upper triangular, form is A = U*D*U**H; = 'L': Lower triangular, form is A = L*D*L**H. [in] N N is INTEGER The order of the matrix A. N >= 0. [in] NRHS NRHS is INTEGER The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0. [in] A A is COMPLEX array, dimension (LDA,N) The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by CHETRF_ROOK. [in] LDA LDA is INTEGER The leading dimension of the array A. LDA >= max(1,N). [in] IPIV IPIV is INTEGER array, dimension (N) Details of the interchanges and the block structure of D as determined by CHETRF_ROOK. [in,out] B B is COMPLEX array, dimension (LDB,NRHS) On entry, the right hand side matrix B. On exit, the solution matrix X. [in] LDB LDB is INTEGER The leading dimension of the array B. LDB >= max(1,N). [out] INFO INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value
Contributors:
November 2013,  Igor Kozachenko,
Computer Science Division,
University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
School of Mathematics,
University of Manchester

Definition at line 134 of file chetrs_rook.f.

136*
137* -- LAPACK computational routine --
138* -- LAPACK is a software package provided by Univ. of Tennessee, --
139* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
140*
141* .. Scalar Arguments ..
142 CHARACTER UPLO
143 INTEGER INFO, LDA, LDB, N, NRHS
144* ..
145* .. Array Arguments ..
146 INTEGER IPIV( * )
147 COMPLEX A( LDA, * ), B( LDB, * )
148* ..
149*
150* =====================================================================
151*
152* .. Parameters ..
153 COMPLEX ONE
154 parameter( one = ( 1.0e+0, 0.0e+0 ) )
155* ..
156* .. Local Scalars ..
157 LOGICAL UPPER
158 INTEGER J, K, KP
159 REAL S
160 COMPLEX AK, AKM1, AKM1K, BK, BKM1, DENOM
161* ..
162* .. External Functions ..
163 LOGICAL LSAME
164 EXTERNAL lsame
165* ..
166* .. External Subroutines ..
167 EXTERNAL cgemv, cgeru, clacgv, csscal, cswap, xerbla
168* ..
169* .. Intrinsic Functions ..
170 INTRINSIC conjg, max, real
171* ..
172* .. Executable Statements ..
173*
174 info = 0
175 upper = lsame( uplo, 'U' )
176 IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
177 info = -1
178 ELSE IF( n.LT.0 ) THEN
179 info = -2
180 ELSE IF( nrhs.LT.0 ) THEN
181 info = -3
182 ELSE IF( lda.LT.max( 1, n ) ) THEN
183 info = -5
184 ELSE IF( ldb.LT.max( 1, n ) ) THEN
185 info = -8
186 END IF
187 IF( info.NE.0 ) THEN
188 CALL xerbla( 'CHETRS_ROOK', -info )
189 RETURN
190 END IF
191*
192* Quick return if possible
193*
194 IF( n.EQ.0 .OR. nrhs.EQ.0 )
195 \$ RETURN
196*
197 IF( upper ) THEN
198*
199* Solve A*X = B, where A = U*D*U**H.
200*
201* First solve U*D*X = B, overwriting B with X.
202*
203* K is the main loop index, decreasing from N to 1 in steps of
204* 1 or 2, depending on the size of the diagonal blocks.
205*
206 k = n
207 10 CONTINUE
208*
209* If K < 1, exit from loop.
210*
211 IF( k.LT.1 )
212 \$ GO TO 30
213*
214 IF( ipiv( k ).GT.0 ) THEN
215*
216* 1 x 1 diagonal block
217*
218* Interchange rows K and IPIV(K).
219*
220 kp = ipiv( k )
221 IF( kp.NE.k )
222 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
223*
224* Multiply by inv(U(K)), where U(K) is the transformation
225* stored in column K of A.
226*
227 CALL cgeru( k-1, nrhs, -one, a( 1, k ), 1, b( k, 1 ), ldb,
228 \$ b( 1, 1 ), ldb )
229*
230* Multiply by the inverse of the diagonal block.
231*
232 s = real( one ) / real( a( k, k ) )
233 CALL csscal( nrhs, s, b( k, 1 ), ldb )
234 k = k - 1
235 ELSE
236*
237* 2 x 2 diagonal block
238*
239* Interchange rows K and -IPIV(K), then K-1 and -IPIV(K-1)
240*
241 kp = -ipiv( k )
242 IF( kp.NE.k )
243 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
244*
245 kp = -ipiv( k-1)
246 IF( kp.NE.k-1 )
247 \$ CALL cswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ), ldb )
248*
249* Multiply by inv(U(K)), where U(K) is the transformation
250* stored in columns K-1 and K of A.
251*
252 CALL cgeru( k-2, nrhs, -one, a( 1, k ), 1, b( k, 1 ), ldb,
253 \$ b( 1, 1 ), ldb )
254 CALL cgeru( k-2, nrhs, -one, a( 1, k-1 ), 1, b( k-1, 1 ),
255 \$ ldb, b( 1, 1 ), ldb )
256*
257* Multiply by the inverse of the diagonal block.
258*
259 akm1k = a( k-1, k )
260 akm1 = a( k-1, k-1 ) / akm1k
261 ak = a( k, k ) / conjg( akm1k )
262 denom = akm1*ak - one
263 DO 20 j = 1, nrhs
264 bkm1 = b( k-1, j ) / akm1k
265 bk = b( k, j ) / conjg( akm1k )
266 b( k-1, j ) = ( ak*bkm1-bk ) / denom
267 b( k, j ) = ( akm1*bk-bkm1 ) / denom
268 20 CONTINUE
269 k = k - 2
270 END IF
271*
272 GO TO 10
273 30 CONTINUE
274*
275* Next solve U**H *X = B, overwriting B with X.
276*
277* K is the main loop index, increasing from 1 to N in steps of
278* 1 or 2, depending on the size of the diagonal blocks.
279*
280 k = 1
281 40 CONTINUE
282*
283* If K > N, exit from loop.
284*
285 IF( k.GT.n )
286 \$ GO TO 50
287*
288 IF( ipiv( k ).GT.0 ) THEN
289*
290* 1 x 1 diagonal block
291*
292* Multiply by inv(U**H(K)), where U(K) is the transformation
293* stored in column K of A.
294*
295 IF( k.GT.1 ) THEN
296 CALL clacgv( nrhs, b( k, 1 ), ldb )
297 CALL cgemv( 'Conjugate transpose', k-1, nrhs, -one, b,
298 \$ ldb, a( 1, k ), 1, one, b( k, 1 ), ldb )
299 CALL clacgv( nrhs, b( k, 1 ), ldb )
300 END IF
301*
302* Interchange rows K and IPIV(K).
303*
304 kp = ipiv( k )
305 IF( kp.NE.k )
306 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
307 k = k + 1
308 ELSE
309*
310* 2 x 2 diagonal block
311*
312* Multiply by inv(U**H(K+1)), where U(K+1) is the transformation
313* stored in columns K and K+1 of A.
314*
315 IF( k.GT.1 ) THEN
316 CALL clacgv( nrhs, b( k, 1 ), ldb )
317 CALL cgemv( 'Conjugate transpose', k-1, nrhs, -one, b,
318 \$ ldb, a( 1, k ), 1, one, b( k, 1 ), ldb )
319 CALL clacgv( nrhs, b( k, 1 ), ldb )
320*
321 CALL clacgv( nrhs, b( k+1, 1 ), ldb )
322 CALL cgemv( 'Conjugate transpose', k-1, nrhs, -one, b,
323 \$ ldb, a( 1, k+1 ), 1, one, b( k+1, 1 ), ldb )
324 CALL clacgv( nrhs, b( k+1, 1 ), ldb )
325 END IF
326*
327* Interchange rows K and -IPIV(K), then K+1 and -IPIV(K+1)
328*
329 kp = -ipiv( k )
330 IF( kp.NE.k )
331 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
332*
333 kp = -ipiv( k+1 )
334 IF( kp.NE.k+1 )
335 \$ CALL cswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ), ldb )
336*
337 k = k + 2
338 END IF
339*
340 GO TO 40
341 50 CONTINUE
342*
343 ELSE
344*
345* Solve A*X = B, where A = L*D*L**H.
346*
347* First solve L*D*X = B, overwriting B with X.
348*
349* K is the main loop index, increasing from 1 to N in steps of
350* 1 or 2, depending on the size of the diagonal blocks.
351*
352 k = 1
353 60 CONTINUE
354*
355* If K > N, exit from loop.
356*
357 IF( k.GT.n )
358 \$ GO TO 80
359*
360 IF( ipiv( k ).GT.0 ) THEN
361*
362* 1 x 1 diagonal block
363*
364* Interchange rows K and IPIV(K).
365*
366 kp = ipiv( k )
367 IF( kp.NE.k )
368 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
369*
370* Multiply by inv(L(K)), where L(K) is the transformation
371* stored in column K of A.
372*
373 IF( k.LT.n )
374 \$ CALL cgeru( n-k, nrhs, -one, a( k+1, k ), 1, b( k, 1 ),
375 \$ ldb, b( k+1, 1 ), ldb )
376*
377* Multiply by the inverse of the diagonal block.
378*
379 s = real( one ) / real( a( k, k ) )
380 CALL csscal( nrhs, s, b( k, 1 ), ldb )
381 k = k + 1
382 ELSE
383*
384* 2 x 2 diagonal block
385*
386* Interchange rows K and -IPIV(K), then K+1 and -IPIV(K+1)
387*
388 kp = -ipiv( k )
389 IF( kp.NE.k )
390 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
391*
392 kp = -ipiv( k+1 )
393 IF( kp.NE.k+1 )
394 \$ CALL cswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ), ldb )
395*
396* Multiply by inv(L(K)), where L(K) is the transformation
397* stored in columns K and K+1 of A.
398*
399 IF( k.LT.n-1 ) THEN
400 CALL cgeru( n-k-1, nrhs, -one, a( k+2, k ), 1, b( k, 1 ),
401 \$ ldb, b( k+2, 1 ), ldb )
402 CALL cgeru( n-k-1, nrhs, -one, a( k+2, k+1 ), 1,
403 \$ b( k+1, 1 ), ldb, b( k+2, 1 ), ldb )
404 END IF
405*
406* Multiply by the inverse of the diagonal block.
407*
408 akm1k = a( k+1, k )
409 akm1 = a( k, k ) / conjg( akm1k )
410 ak = a( k+1, k+1 ) / akm1k
411 denom = akm1*ak - one
412 DO 70 j = 1, nrhs
413 bkm1 = b( k, j ) / conjg( akm1k )
414 bk = b( k+1, j ) / akm1k
415 b( k, j ) = ( ak*bkm1-bk ) / denom
416 b( k+1, j ) = ( akm1*bk-bkm1 ) / denom
417 70 CONTINUE
418 k = k + 2
419 END IF
420*
421 GO TO 60
422 80 CONTINUE
423*
424* Next solve L**H *X = B, overwriting B with X.
425*
426* K is the main loop index, decreasing from N to 1 in steps of
427* 1 or 2, depending on the size of the diagonal blocks.
428*
429 k = n
430 90 CONTINUE
431*
432* If K < 1, exit from loop.
433*
434 IF( k.LT.1 )
435 \$ GO TO 100
436*
437 IF( ipiv( k ).GT.0 ) THEN
438*
439* 1 x 1 diagonal block
440*
441* Multiply by inv(L**H(K)), where L(K) is the transformation
442* stored in column K of A.
443*
444 IF( k.LT.n ) THEN
445 CALL clacgv( nrhs, b( k, 1 ), ldb )
446 CALL cgemv( 'Conjugate transpose', n-k, nrhs, -one,
447 \$ b( k+1, 1 ), ldb, a( k+1, k ), 1, one,
448 \$ b( k, 1 ), ldb )
449 CALL clacgv( nrhs, b( k, 1 ), ldb )
450 END IF
451*
452* Interchange rows K and IPIV(K).
453*
454 kp = ipiv( k )
455 IF( kp.NE.k )
456 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
457 k = k - 1
458 ELSE
459*
460* 2 x 2 diagonal block
461*
462* Multiply by inv(L**H(K-1)), where L(K-1) is the transformation
463* stored in columns K-1 and K of A.
464*
465 IF( k.LT.n ) THEN
466 CALL clacgv( nrhs, b( k, 1 ), ldb )
467 CALL cgemv( 'Conjugate transpose', n-k, nrhs, -one,
468 \$ b( k+1, 1 ), ldb, a( k+1, k ), 1, one,
469 \$ b( k, 1 ), ldb )
470 CALL clacgv( nrhs, b( k, 1 ), ldb )
471*
472 CALL clacgv( nrhs, b( k-1, 1 ), ldb )
473 CALL cgemv( 'Conjugate transpose', n-k, nrhs, -one,
474 \$ b( k+1, 1 ), ldb, a( k+1, k-1 ), 1, one,
475 \$ b( k-1, 1 ), ldb )
476 CALL clacgv( nrhs, b( k-1, 1 ), ldb )
477 END IF
478*
479* Interchange rows K and -IPIV(K), then K-1 and -IPIV(K-1)
480*
481 kp = -ipiv( k )
482 IF( kp.NE.k )
483 \$ CALL cswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
484*
485 kp = -ipiv( k-1 )
486 IF( kp.NE.k-1 )
487 \$ CALL cswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ), ldb )
488*
489 k = k - 2
490 END IF
491*
492 GO TO 90
493 100 CONTINUE
494 END IF
495*
496 RETURN
497*
498* End of CHETRS_ROOK
499*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine cgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
CGEMV
Definition cgemv.f:160
subroutine cgeru(m, n, alpha, x, incx, y, incy, a, lda)
CGERU
Definition cgeru.f:130
subroutine clacgv(n, x, incx)
CLACGV conjugates a complex vector.
Definition clacgv.f:74
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine csscal(n, sa, cx, incx)
CSSCAL
Definition csscal.f:78
subroutine cswap(n, cx, incx, cy, incy)
CSWAP
Definition cswap.f:81
Here is the call graph for this function:
Here is the caller graph for this function: