LAPACK  3.8.0 LAPACK: Linear Algebra PACKage

## ◆ dsytrs_rook()

 subroutine dsytrs_rook ( character UPLO, integer N, integer NRHS, double precision, dimension( lda, * ) A, integer LDA, integer, dimension( * ) IPIV, double precision, dimension( ldb, * ) B, integer LDB, integer INFO )

DSYTRS_ROOK

Purpose:
``` DSYTRS_ROOK solves a system of linear equations A*X = B with
a real symmetric matrix A using the factorization A = U*D*U**T or
A = L*D*L**T computed by DSYTRF_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**T; = 'L': Lower triangular, form is A = L*D*L**T.``` [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 DOUBLE PRECISION array, dimension (LDA,N) The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by DSYTRF_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 DSYTRF_ROOK.``` [in,out] B ``` B is DOUBLE PRECISION 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```
Date
April 2012
Contributors:
```   April 2012, 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 138 of file dsytrs_rook.f.

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