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

## ◆ ssptrs()

 subroutine ssptrs ( character uplo, integer n, integer nrhs, real, dimension( * ) ap, integer, dimension( * ) ipiv, real, dimension( ldb, * ) b, integer ldb, integer info )

SSPTRS

Purpose:
``` SSPTRS solves a system of linear equations A*X = B with a real
symmetric matrix A stored in packed format using the factorization
A = U*D*U**T or A = L*D*L**T computed by SSPTRF.```
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] AP ``` AP is REAL array, dimension (N*(N+1)/2) The block diagonal matrix D and the multipliers used to obtain the factor U or L as computed by SSPTRF, stored as a packed triangular matrix.``` [in] IPIV ``` IPIV is INTEGER array, dimension (N) Details of the interchanges and the block structure of D as determined by SSPTRF.``` [in,out] B ``` B is REAL 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```

Definition at line 114 of file ssptrs.f.

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