OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
double
Reference BLAS » Level2

Functions

subroutine dgbmv (trans, m, n, kl, ku, alpha, a, lda, x, incx, beta, y, incy)
 DGBMV
subroutine dgemv (trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
 DGEMV
subroutine dger (m, n, alpha, x, incx, y, incy, a, lda)
 DGER
subroutine dsbmv (uplo, n, k, alpha, a, lda, x, incx, beta, y, incy)
 DSBMV
subroutine dspmv (uplo, n, alpha, ap, x, incx, beta, y, incy)
 DSPMV
subroutine dspr (uplo, n, alpha, x, incx, ap)
 DSPR
subroutine dspr2 (uplo, n, alpha, x, incx, y, incy, ap)
 DSPR2
subroutine dsymv (uplo, n, alpha, a, lda, x, incx, beta, y, incy)
 DSYMV
subroutine dsyr (uplo, n, alpha, x, incx, a, lda)
 DSYR
subroutine dsyr2 (uplo, n, alpha, x, incx, y, incy, a, lda)
 DSYR2
subroutine dtbmv (uplo, trans, diag, n, k, a, lda, x, incx)
 DTBMV
subroutine dtbsv (uplo, trans, diag, n, k, a, lda, x, incx)
 DTBSV
subroutine dtpmv (uplo, trans, diag, n, ap, x, incx)
 DTPMV
subroutine dtpsv (uplo, trans, diag, n, ap, x, incx)
 DTPSV
subroutine dtrmv (uplo, trans, diag, n, a, lda, x, incx)
 DTRMV

Detailed Description

This is the group of double LEVEL 2 BLAS routines.

Function Documentation

◆ dgbmv()

subroutine dgbmv ( character trans,
integer m,
integer n,
integer kl,
integer ku,
double precision alpha,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx,
double precision beta,
double precision, dimension(*) y,
integer incy )

DGBMV

Purpose:
!>
!> DGBMV  performs one of the matrix-vector operations
!>
!>    y := alpha*A*x + beta*y,   or   y := alpha*A**T*x + beta*y,
!>
!> where alpha and beta are scalars, x and y are vectors and A is an
!> m by n band matrix, with kl sub-diagonals and ku super-diagonals.
!>
Parameters
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the operation to be performed as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   y := alpha*A*x + beta*y.
!>
!>              TRANS = 'T' or 't'   y := alpha*A**T*x + beta*y.
!>
!>              TRANS = 'C' or 'c'   y := alpha*A**T*x + beta*y.
!>
[in]M
!>          M is INTEGER
!>           On entry, M specifies the number of rows of the matrix A.
!>           M must be at least zero.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the number of columns of the matrix A.
!>           N must be at least zero.
!>
[in]KL
!>          KL is INTEGER
!>           On entry, KL specifies the number of sub-diagonals of the
!>           matrix A. KL must satisfy  0 .le. KL.
!>
[in]KU
!>          KU is INTEGER
!>           On entry, KU specifies the number of super-diagonals of the
!>           matrix A. KU must satisfy  0 .le. KU.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry, the leading ( kl + ku + 1 ) by n part of the
!>           array A must contain the matrix of coefficients, supplied
!>           column by column, with the leading diagonal of the matrix in
!>           row ( ku + 1 ) of the array, the first super-diagonal
!>           starting at position 2 in row ku, the first sub-diagonal
!>           starting at position 1 in row ( ku + 2 ), and so on.
!>           Elements in the array A that do not correspond to elements
!>           in the band matrix (such as the top left ku by ku triangle)
!>           are not referenced.
!>           The following program segment will transfer a band matrix
!>           from conventional full matrix storage to band storage:
!>
!>                 DO 20, J = 1, N
!>                    K = KU + 1 - J
!>                    DO 10, I = MAX( 1, J - KU ), MIN( M, J + KL )
!>                       A( K + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           ( kl + ku + 1 ).
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
!>           and at least
!>           ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
!>           Before entry, the incremented array X must contain the
!>           vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]BETA
!>          BETA is DOUBLE PRECISION.
!>           On entry, BETA specifies the scalar beta. When BETA is
!>           supplied as zero then Y need not be set on input.
!>
[in,out]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
!>           and at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
!>           Before entry, the incremented array Y must contain the
!>           vector y. On exit, Y is overwritten by the updated vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 184 of file dgbmv.f.

185*
186* -- Reference BLAS level2 routine --
187* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
188* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
189*
190* .. Scalar Arguments ..
191 DOUBLE PRECISION ALPHA,BETA
192 INTEGER INCX,INCY,KL,KU,LDA,M,N
193 CHARACTER TRANS
194* ..
195* .. Array Arguments ..
196 DOUBLE PRECISION A(LDA,*),X(*),Y(*)
197* ..
198*
199* =====================================================================
200*
201* .. Parameters ..
202 DOUBLE PRECISION ONE,ZERO
203 parameter(one=1.0d+0,zero=0.0d+0)
204* ..
205* .. Local Scalars ..
206 DOUBLE PRECISION TEMP
207 INTEGER I,INFO,IX,IY,J,JX,JY,K,KUP1,KX,KY,LENX,LENY
208* ..
209* .. External Functions ..
210 LOGICAL LSAME
211 EXTERNAL lsame
212* ..
213* .. External Subroutines ..
214 EXTERNAL xerbla
215* ..
216* .. Intrinsic Functions ..
217 INTRINSIC max,min
218* ..
219*
220* Test the input parameters.
221*
222 info = 0
223 IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
224 + .NOT.lsame(trans,'C')) THEN
225 info = 1
226 ELSE IF (m.LT.0) THEN
227 info = 2
228 ELSE IF (n.LT.0) THEN
229 info = 3
230 ELSE IF (kl.LT.0) THEN
231 info = 4
232 ELSE IF (ku.LT.0) THEN
233 info = 5
234 ELSE IF (lda.LT. (kl+ku+1)) THEN
235 info = 8
236 ELSE IF (incx.EQ.0) THEN
237 info = 10
238 ELSE IF (incy.EQ.0) THEN
239 info = 13
240 END IF
241 IF (info.NE.0) THEN
242 CALL xerbla('DGBMV ',info)
243 RETURN
244 END IF
245*
246* Quick return if possible.
247*
248 IF ((m.EQ.0) .OR. (n.EQ.0) .OR.
249 + ((alpha.EQ.zero).AND. (beta.EQ.one))) RETURN
250*
251* Set LENX and LENY, the lengths of the vectors x and y, and set
252* up the start points in X and Y.
253*
254 IF (lsame(trans,'N')) THEN
255 lenx = n
256 leny = m
257 ELSE
258 lenx = m
259 leny = n
260 END IF
261 IF (incx.GT.0) THEN
262 kx = 1
263 ELSE
264 kx = 1 - (lenx-1)*incx
265 END IF
266 IF (incy.GT.0) THEN
267 ky = 1
268 ELSE
269 ky = 1 - (leny-1)*incy
270 END IF
271*
272* Start the operations. In this version the elements of A are
273* accessed sequentially with one pass through the band part of A.
274*
275* First form y := beta*y.
276*
277 IF (beta.NE.one) THEN
278 IF (incy.EQ.1) THEN
279 IF (beta.EQ.zero) THEN
280 DO 10 i = 1,leny
281 y(i) = zero
282 10 CONTINUE
283 ELSE
284 DO 20 i = 1,leny
285 y(i) = beta*y(i)
286 20 CONTINUE
287 END IF
288 ELSE
289 iy = ky
290 IF (beta.EQ.zero) THEN
291 DO 30 i = 1,leny
292 y(iy) = zero
293 iy = iy + incy
294 30 CONTINUE
295 ELSE
296 DO 40 i = 1,leny
297 y(iy) = beta*y(iy)
298 iy = iy + incy
299 40 CONTINUE
300 END IF
301 END IF
302 END IF
303 IF (alpha.EQ.zero) RETURN
304 kup1 = ku + 1
305 IF (lsame(trans,'N')) THEN
306*
307* Form y := alpha*A*x + y.
308*
309 jx = kx
310 IF (incy.EQ.1) THEN
311 DO 60 j = 1,n
312 temp = alpha*x(jx)
313 k = kup1 - j
314 DO 50 i = max(1,j-ku),min(m,j+kl)
315 y(i) = y(i) + temp*a(k+i,j)
316 50 CONTINUE
317 jx = jx + incx
318 60 CONTINUE
319 ELSE
320 DO 80 j = 1,n
321 temp = alpha*x(jx)
322 iy = ky
323 k = kup1 - j
324 DO 70 i = max(1,j-ku),min(m,j+kl)
325 y(iy) = y(iy) + temp*a(k+i,j)
326 iy = iy + incy
327 70 CONTINUE
328 jx = jx + incx
329 IF (j.GT.ku) ky = ky + incy
330 80 CONTINUE
331 END IF
332 ELSE
333*
334* Form y := alpha*A**T*x + y.
335*
336 jy = ky
337 IF (incx.EQ.1) THEN
338 DO 100 j = 1,n
339 temp = zero
340 k = kup1 - j
341 DO 90 i = max(1,j-ku),min(m,j+kl)
342 temp = temp + a(k+i,j)*x(i)
343 90 CONTINUE
344 y(jy) = y(jy) + alpha*temp
345 jy = jy + incy
346 100 CONTINUE
347 ELSE
348 DO 120 j = 1,n
349 temp = zero
350 ix = kx
351 k = kup1 - j
352 DO 110 i = max(1,j-ku),min(m,j+kl)
353 temp = temp + a(k+i,j)*x(ix)
354 ix = ix + incx
355 110 CONTINUE
356 y(jy) = y(jy) + alpha*temp
357 jy = jy + incy
358 IF (j.GT.ku) kx = kx + incx
359 120 CONTINUE
360 END IF
361 END IF
362*
363 RETURN
364*
365* End of DGBMV
366*
alpha
#define alpha
Definition eval.h:35
xerbla
subroutine xerbla(srname, info)
XERBLA
Definition xerbla.f:60
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:53
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21

◆ dgemv()

subroutine dgemv ( character trans,
integer m,
integer n,
double precision alpha,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx,
double precision beta,
double precision, dimension(*) y,
integer incy )

DGEMV

Purpose:
!>
!> DGEMV  performs one of the matrix-vector operations
!>
!>    y := alpha*A*x + beta*y,   or   y := alpha*A**T*x + beta*y,
!>
!> where alpha and beta are scalars, x and y are vectors and A is an
!> m by n matrix.
!>
Parameters
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the operation to be performed as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   y := alpha*A*x + beta*y.
!>
!>              TRANS = 'T' or 't'   y := alpha*A**T*x + beta*y.
!>
!>              TRANS = 'C' or 'c'   y := alpha*A**T*x + beta*y.
!>
[in]M
!>          M is INTEGER
!>           On entry, M specifies the number of rows of the matrix A.
!>           M must be at least zero.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the number of columns of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry, the leading m by n part of the array A must
!>           contain the matrix of coefficients.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           max( 1, m ).
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
!>           and at least
!>           ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
!>           Before entry, the incremented array X must contain the
!>           vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]BETA
!>          BETA is DOUBLE PRECISION.
!>           On entry, BETA specifies the scalar beta. When BETA is
!>           supplied as zero then Y need not be set on input.
!>
[in,out]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
!>           and at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
!>           Before entry with BETA non-zero, the incremented array Y
!>           must contain the vector y. On exit, Y is overwritten by the
!>           updated vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 155 of file dgemv.f.

156*
157* -- Reference BLAS level2 routine --
158* -- Reference BLAS 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 DOUBLE PRECISION ALPHA,BETA
163 INTEGER INCX,INCY,LDA,M,N
164 CHARACTER TRANS
165* ..
166* .. Array Arguments ..
167 DOUBLE PRECISION A(LDA,*),X(*),Y(*)
168* ..
169*
170* =====================================================================
171*
172* .. Parameters ..
173 DOUBLE PRECISION ONE,ZERO
174 parameter(one=1.0d+0,zero=0.0d+0)
175* ..
176* .. Local Scalars ..
177 DOUBLE PRECISION TEMP
178 INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY,LENX,LENY
179* ..
180* .. External Functions ..
181 LOGICAL LSAME
182 EXTERNAL lsame
183* ..
184* .. External Subroutines ..
185 EXTERNAL xerbla
186* ..
187* .. Intrinsic Functions ..
188 INTRINSIC max
189* ..
190*
191* Test the input parameters.
192*
193 info = 0
194 IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
195 + .NOT.lsame(trans,'C')) THEN
196 info = 1
197 ELSE IF (m.LT.0) THEN
198 info = 2
199 ELSE IF (n.LT.0) THEN
200 info = 3
201 ELSE IF (lda.LT.max(1,m)) THEN
202 info = 6
203 ELSE IF (incx.EQ.0) THEN
204 info = 8
205 ELSE IF (incy.EQ.0) THEN
206 info = 11
207 END IF
208 IF (info.NE.0) THEN
209 CALL xerbla('DGEMV ',info)
210 RETURN
211 END IF
212*
213* Quick return if possible.
214*
215 IF ((m.EQ.0) .OR. (n.EQ.0) .OR.
216 + ((alpha.EQ.zero).AND. (beta.EQ.one))) RETURN
217*
218* Set LENX and LENY, the lengths of the vectors x and y, and set
219* up the start points in X and Y.
220*
221 IF (lsame(trans,'N')) THEN
222 lenx = n
223 leny = m
224 ELSE
225 lenx = m
226 leny = n
227 END IF
228 IF (incx.GT.0) THEN
229 kx = 1
230 ELSE
231 kx = 1 - (lenx-1)*incx
232 END IF
233 IF (incy.GT.0) THEN
234 ky = 1
235 ELSE
236 ky = 1 - (leny-1)*incy
237 END IF
238*
239* Start the operations. In this version the elements of A are
240* accessed sequentially with one pass through A.
241*
242* First form y := beta*y.
243*
244 IF (beta.NE.one) THEN
245 IF (incy.EQ.1) THEN
246 IF (beta.EQ.zero) THEN
247 DO 10 i = 1,leny
248 y(i) = zero
249 10 CONTINUE
250 ELSE
251 DO 20 i = 1,leny
252 y(i) = beta*y(i)
253 20 CONTINUE
254 END IF
255 ELSE
256 iy = ky
257 IF (beta.EQ.zero) THEN
258 DO 30 i = 1,leny
259 y(iy) = zero
260 iy = iy + incy
261 30 CONTINUE
262 ELSE
263 DO 40 i = 1,leny
264 y(iy) = beta*y(iy)
265 iy = iy + incy
266 40 CONTINUE
267 END IF
268 END IF
269 END IF
270 IF (alpha.EQ.zero) RETURN
271 IF (lsame(trans,'N')) THEN
272*
273* Form y := alpha*A*x + y.
274*
275 jx = kx
276 IF (incy.EQ.1) THEN
277 DO 60 j = 1,n
278 temp = alpha*x(jx)
279 DO 50 i = 1,m
280 y(i) = y(i) + temp*a(i,j)
281 50 CONTINUE
282 jx = jx + incx
283 60 CONTINUE
284 ELSE
285 DO 80 j = 1,n
286 temp = alpha*x(jx)
287 iy = ky
288 DO 70 i = 1,m
289 y(iy) = y(iy) + temp*a(i,j)
290 iy = iy + incy
291 70 CONTINUE
292 jx = jx + incx
293 80 CONTINUE
294 END IF
295 ELSE
296*
297* Form y := alpha*A**T*x + y.
298*
299 jy = ky
300 IF (incx.EQ.1) THEN
301 DO 100 j = 1,n
302 temp = zero
303 DO 90 i = 1,m
304 temp = temp + a(i,j)*x(i)
305 90 CONTINUE
306 y(jy) = y(jy) + alpha*temp
307 jy = jy + incy
308 100 CONTINUE
309 ELSE
310 DO 120 j = 1,n
311 temp = zero
312 ix = kx
313 DO 110 i = 1,m
314 temp = temp + a(i,j)*x(ix)
315 ix = ix + incx
316 110 CONTINUE
317 y(jy) = y(jy) + alpha*temp
318 jy = jy + incy
319 120 CONTINUE
320 END IF
321 END IF
322*
323 RETURN
324*
325* End of DGEMV
326*

◆ dger()

subroutine dger ( integer m,
integer n,
double precision alpha,
double precision, dimension(*) x,
integer incx,
double precision, dimension(*) y,
integer incy,
double precision, dimension(lda,*) a,
integer lda )

DGER

Purpose:
!>
!> DGER   performs the rank 1 operation
!>
!>    A := alpha*x*y**T + A,
!>
!> where alpha is a scalar, x is an m element vector, y is an n element
!> vector and A is an m by n matrix.
!>
Parameters
[in]M
!>          M is INTEGER
!>           On entry, M specifies the number of rows of the matrix A.
!>           M must be at least zero.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the number of columns of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( m - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the m
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ).
!>           Before entry, the incremented array Y must contain the n
!>           element vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
[in,out]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry, the leading m by n part of the array A must
!>           contain the matrix of coefficients. On exit, A is
!>           overwritten by the updated matrix.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           max( 1, m ).
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 129 of file dger.f.

130*
131* -- Reference BLAS level2 routine --
132* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
133* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
134*
135* .. Scalar Arguments ..
136 DOUBLE PRECISION ALPHA
137 INTEGER INCX,INCY,LDA,M,N
138* ..
139* .. Array Arguments ..
140 DOUBLE PRECISION A(LDA,*),X(*),Y(*)
141* ..
142*
143* =====================================================================
144*
145* .. Parameters ..
146 DOUBLE PRECISION ZERO
147 parameter(zero=0.0d+0)
148* ..
149* .. Local Scalars ..
150 DOUBLE PRECISION TEMP
151 INTEGER I,INFO,IX,J,JY,KX
152* ..
153* .. External Subroutines ..
154 EXTERNAL xerbla
155* ..
156* .. Intrinsic Functions ..
157 INTRINSIC max
158* ..
159*
160* Test the input parameters.
161*
162 info = 0
163 IF (m.LT.0) THEN
164 info = 1
165 ELSE IF (n.LT.0) THEN
166 info = 2
167 ELSE IF (incx.EQ.0) THEN
168 info = 5
169 ELSE IF (incy.EQ.0) THEN
170 info = 7
171 ELSE IF (lda.LT.max(1,m)) THEN
172 info = 9
173 END IF
174 IF (info.NE.0) THEN
175 CALL xerbla('DGER ',info)
176 RETURN
177 END IF
178*
179* Quick return if possible.
180*
181 IF ((m.EQ.0) .OR. (n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
182*
183* Start the operations. In this version the elements of A are
184* accessed sequentially with one pass through A.
185*
186 IF (incy.GT.0) THEN
187 jy = 1
188 ELSE
189 jy = 1 - (n-1)*incy
190 END IF
191 IF (incx.EQ.1) THEN
192 DO 20 j = 1,n
193 IF (y(jy).NE.zero) THEN
194 temp = alpha*y(jy)
195 DO 10 i = 1,m
196 a(i,j) = a(i,j) + x(i)*temp
197 10 CONTINUE
198 END IF
199 jy = jy + incy
200 20 CONTINUE
201 ELSE
202 IF (incx.GT.0) THEN
203 kx = 1
204 ELSE
205 kx = 1 - (m-1)*incx
206 END IF
207 DO 40 j = 1,n
208 IF (y(jy).NE.zero) THEN
209 temp = alpha*y(jy)
210 ix = kx
211 DO 30 i = 1,m
212 a(i,j) = a(i,j) + x(ix)*temp
213 ix = ix + incx
214 30 CONTINUE
215 END IF
216 jy = jy + incy
217 40 CONTINUE
218 END IF
219*
220 RETURN
221*
222* End of DGER
223*

◆ dsbmv()

subroutine dsbmv ( character uplo,
integer n,
integer k,
double precision alpha,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx,
double precision beta,
double precision, dimension(*) y,
integer incy )

DSBMV

Purpose:
!>
!> DSBMV  performs the matrix-vector  operation
!>
!>    y := alpha*A*x + beta*y,
!>
!> where alpha and beta are scalars, x and y are n element vectors and
!> A is an n by n symmetric band matrix, with k super-diagonals.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the band matrix A is being supplied as
!>           follows:
!>
!>              UPLO = 'U' or 'u'   The upper triangular part of A is
!>                                  being supplied.
!>
!>              UPLO = 'L' or 'l'   The lower triangular part of A is
!>                                  being supplied.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]K
!>          K is INTEGER
!>           On entry, K specifies the number of super-diagonals of the
!>           matrix A. K must satisfy  0 .le. K.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
!>           by n part of the array A must contain the upper triangular
!>           band part of the symmetric matrix, supplied column by
!>           column, with the leading diagonal of the matrix in row
!>           ( k + 1 ) of the array, the first super-diagonal starting at
!>           position 2 in row k, and so on. The top left k by k triangle
!>           of the array A is not referenced.
!>           The following program segment will transfer the upper
!>           triangular part of a symmetric band matrix from conventional
!>           full matrix storage to band storage:
!>
!>                 DO 20, J = 1, N
!>                    M = K + 1 - J
!>                    DO 10, I = MAX( 1, J - K ), J
!>                       A( M + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
!>           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
!>           by n part of the array A must contain the lower triangular
!>           band part of the symmetric matrix, supplied column by
!>           column, with the leading diagonal of the matrix in row 1 of
!>           the array, the first sub-diagonal starting at position 1 in
!>           row 2, and so on. The bottom right k by k triangle of the
!>           array A is not referenced.
!>           The following program segment will transfer the lower
!>           triangular part of a symmetric band matrix from conventional
!>           full matrix storage to band storage:
!>
!>                 DO 20, J = 1, N
!>                    M = 1 - J
!>                    DO 10, I = J, MIN( N, J + K )
!>                       A( M + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           ( k + 1 ).
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the
!>           vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]BETA
!>          BETA is DOUBLE PRECISION.
!>           On entry, BETA specifies the scalar beta.
!>
[in,out]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ).
!>           Before entry, the incremented array Y must contain the
!>           vector y. On exit, Y is overwritten by the updated vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 183 of file dsbmv.f.

184*
185* -- Reference BLAS level2 routine --
186* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
187* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
188*
189* .. Scalar Arguments ..
190 DOUBLE PRECISION ALPHA,BETA
191 INTEGER INCX,INCY,K,LDA,N
192 CHARACTER UPLO
193* ..
194* .. Array Arguments ..
195 DOUBLE PRECISION A(LDA,*),X(*),Y(*)
196* ..
197*
198* =====================================================================
199*
200* .. Parameters ..
201 DOUBLE PRECISION ONE,ZERO
202 parameter(one=1.0d+0,zero=0.0d+0)
203* ..
204* .. Local Scalars ..
205 DOUBLE PRECISION TEMP1,TEMP2
206 INTEGER I,INFO,IX,IY,J,JX,JY,KPLUS1,KX,KY,L
207* ..
208* .. External Functions ..
209 LOGICAL LSAME
210 EXTERNAL lsame
211* ..
212* .. External Subroutines ..
213 EXTERNAL xerbla
214* ..
215* .. Intrinsic Functions ..
216 INTRINSIC max,min
217* ..
218*
219* Test the input parameters.
220*
221 info = 0
222 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
223 info = 1
224 ELSE IF (n.LT.0) THEN
225 info = 2
226 ELSE IF (k.LT.0) THEN
227 info = 3
228 ELSE IF (lda.LT. (k+1)) THEN
229 info = 6
230 ELSE IF (incx.EQ.0) THEN
231 info = 8
232 ELSE IF (incy.EQ.0) THEN
233 info = 11
234 END IF
235 IF (info.NE.0) THEN
236 CALL xerbla('DSBMV ',info)
237 RETURN
238 END IF
239*
240* Quick return if possible.
241*
242 IF ((n.EQ.0) .OR. ((alpha.EQ.zero).AND. (beta.EQ.one))) RETURN
243*
244* Set up the start points in X and Y.
245*
246 IF (incx.GT.0) THEN
247 kx = 1
248 ELSE
249 kx = 1 - (n-1)*incx
250 END IF
251 IF (incy.GT.0) THEN
252 ky = 1
253 ELSE
254 ky = 1 - (n-1)*incy
255 END IF
256*
257* Start the operations. In this version the elements of the array A
258* are accessed sequentially with one pass through A.
259*
260* First form y := beta*y.
261*
262 IF (beta.NE.one) THEN
263 IF (incy.EQ.1) THEN
264 IF (beta.EQ.zero) THEN
265 DO 10 i = 1,n
266 y(i) = zero
267 10 CONTINUE
268 ELSE
269 DO 20 i = 1,n
270 y(i) = beta*y(i)
271 20 CONTINUE
272 END IF
273 ELSE
274 iy = ky
275 IF (beta.EQ.zero) THEN
276 DO 30 i = 1,n
277 y(iy) = zero
278 iy = iy + incy
279 30 CONTINUE
280 ELSE
281 DO 40 i = 1,n
282 y(iy) = beta*y(iy)
283 iy = iy + incy
284 40 CONTINUE
285 END IF
286 END IF
287 END IF
288 IF (alpha.EQ.zero) RETURN
289 IF (lsame(uplo,'U')) THEN
290*
291* Form y when upper triangle of A is stored.
292*
293 kplus1 = k + 1
294 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
295 DO 60 j = 1,n
296 temp1 = alpha*x(j)
297 temp2 = zero
298 l = kplus1 - j
299 DO 50 i = max(1,j-k),j - 1
300 y(i) = y(i) + temp1*a(l+i,j)
301 temp2 = temp2 + a(l+i,j)*x(i)
302 50 CONTINUE
303 y(j) = y(j) + temp1*a(kplus1,j) + alpha*temp2
304 60 CONTINUE
305 ELSE
306 jx = kx
307 jy = ky
308 DO 80 j = 1,n
309 temp1 = alpha*x(jx)
310 temp2 = zero
311 ix = kx
312 iy = ky
313 l = kplus1 - j
314 DO 70 i = max(1,j-k),j - 1
315 y(iy) = y(iy) + temp1*a(l+i,j)
316 temp2 = temp2 + a(l+i,j)*x(ix)
317 ix = ix + incx
318 iy = iy + incy
319 70 CONTINUE
320 y(jy) = y(jy) + temp1*a(kplus1,j) + alpha*temp2
321 jx = jx + incx
322 jy = jy + incy
323 IF (j.GT.k) THEN
324 kx = kx + incx
325 ky = ky + incy
326 END IF
327 80 CONTINUE
328 END IF
329 ELSE
330*
331* Form y when lower triangle of A is stored.
332*
333 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
334 DO 100 j = 1,n
335 temp1 = alpha*x(j)
336 temp2 = zero
337 y(j) = y(j) + temp1*a(1,j)
338 l = 1 - j
339 DO 90 i = j + 1,min(n,j+k)
340 y(i) = y(i) + temp1*a(l+i,j)
341 temp2 = temp2 + a(l+i,j)*x(i)
342 90 CONTINUE
343 y(j) = y(j) + alpha*temp2
344 100 CONTINUE
345 ELSE
346 jx = kx
347 jy = ky
348 DO 120 j = 1,n
349 temp1 = alpha*x(jx)
350 temp2 = zero
351 y(jy) = y(jy) + temp1*a(1,j)
352 l = 1 - j
353 ix = jx
354 iy = jy
355 DO 110 i = j + 1,min(n,j+k)
356 ix = ix + incx
357 iy = iy + incy
358 y(iy) = y(iy) + temp1*a(l+i,j)
359 temp2 = temp2 + a(l+i,j)*x(ix)
360 110 CONTINUE
361 y(jy) = y(jy) + alpha*temp2
362 jx = jx + incx
363 jy = jy + incy
364 120 CONTINUE
365 END IF
366 END IF
367*
368 RETURN
369*
370* End of DSBMV
371*

◆ dspmv()

subroutine dspmv ( character uplo,
integer n,
double precision alpha,
double precision, dimension(*) ap,
double precision, dimension(*) x,
integer incx,
double precision beta,
double precision, dimension(*) y,
integer incy )

DSPMV

Purpose:
!>
!> DSPMV  performs the matrix-vector operation
!>
!>    y := alpha*A*x + beta*y,
!>
!> where alpha and beta are scalars, x and y are n element vectors and
!> A is an n by n symmetric matrix, supplied in packed form.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the matrix A is supplied in the packed
!>           array AP as follows:
!>
!>              UPLO = 'U' or 'u'   The upper triangular part of A is
!>                                  supplied in AP.
!>
!>              UPLO = 'L' or 'l'   The lower triangular part of A is
!>                                  supplied in AP.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]AP
!>          AP is DOUBLE PRECISION array, dimension at least
!>           ( ( n*( n + 1 ) )/2 ).
!>           Before entry with UPLO = 'U' or 'u', the array AP must
!>           contain the upper triangular part of the symmetric matrix
!>           packed sequentially, column by column, so that AP( 1 )
!>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
!>           and a( 2, 2 ) respectively, and so on.
!>           Before entry with UPLO = 'L' or 'l', the array AP must
!>           contain the lower triangular part of the symmetric matrix
!>           packed sequentially, column by column, so that AP( 1 )
!>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
!>           and a( 3, 1 ) respectively, and so on.
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]BETA
!>          BETA is DOUBLE PRECISION.
!>           On entry, BETA specifies the scalar beta. When BETA is
!>           supplied as zero then Y need not be set on input.
!>
[in,out]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ).
!>           Before entry, the incremented array Y must contain the n
!>           element vector y. On exit, Y is overwritten by the updated
!>           vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 146 of file dspmv.f.

147*
148* -- Reference BLAS level2 routine --
149* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
150* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
151*
152* .. Scalar Arguments ..
153 DOUBLE PRECISION ALPHA,BETA
154 INTEGER INCX,INCY,N
155 CHARACTER UPLO
156* ..
157* .. Array Arguments ..
158 DOUBLE PRECISION AP(*),X(*),Y(*)
159* ..
160*
161* =====================================================================
162*
163* .. Parameters ..
164 DOUBLE PRECISION ONE,ZERO
165 parameter(one=1.0d+0,zero=0.0d+0)
166* ..
167* .. Local Scalars ..
168 DOUBLE PRECISION TEMP1,TEMP2
169 INTEGER I,INFO,IX,IY,J,JX,JY,K,KK,KX,KY
170* ..
171* .. External Functions ..
172 LOGICAL LSAME
173 EXTERNAL lsame
174* ..
175* .. External Subroutines ..
176 EXTERNAL xerbla
177* ..
178*
179* Test the input parameters.
180*
181 info = 0
182 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
183 info = 1
184 ELSE IF (n.LT.0) THEN
185 info = 2
186 ELSE IF (incx.EQ.0) THEN
187 info = 6
188 ELSE IF (incy.EQ.0) THEN
189 info = 9
190 END IF
191 IF (info.NE.0) THEN
192 CALL xerbla('DSPMV ',info)
193 RETURN
194 END IF
195*
196* Quick return if possible.
197*
198 IF ((n.EQ.0) .OR. ((alpha.EQ.zero).AND. (beta.EQ.one))) RETURN
199*
200* Set up the start points in X and Y.
201*
202 IF (incx.GT.0) THEN
203 kx = 1
204 ELSE
205 kx = 1 - (n-1)*incx
206 END IF
207 IF (incy.GT.0) THEN
208 ky = 1
209 ELSE
210 ky = 1 - (n-1)*incy
211 END IF
212*
213* Start the operations. In this version the elements of the array AP
214* are accessed sequentially with one pass through AP.
215*
216* First form y := beta*y.
217*
218 IF (beta.NE.one) THEN
219 IF (incy.EQ.1) THEN
220 IF (beta.EQ.zero) THEN
221 DO 10 i = 1,n
222 y(i) = zero
223 10 CONTINUE
224 ELSE
225 DO 20 i = 1,n
226 y(i) = beta*y(i)
227 20 CONTINUE
228 END IF
229 ELSE
230 iy = ky
231 IF (beta.EQ.zero) THEN
232 DO 30 i = 1,n
233 y(iy) = zero
234 iy = iy + incy
235 30 CONTINUE
236 ELSE
237 DO 40 i = 1,n
238 y(iy) = beta*y(iy)
239 iy = iy + incy
240 40 CONTINUE
241 END IF
242 END IF
243 END IF
244 IF (alpha.EQ.zero) RETURN
245 kk = 1
246 IF (lsame(uplo,'U')) THEN
247*
248* Form y when AP contains the upper triangle.
249*
250 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
251 DO 60 j = 1,n
252 temp1 = alpha*x(j)
253 temp2 = zero
254 k = kk
255 DO 50 i = 1,j - 1
256 y(i) = y(i) + temp1*ap(k)
257 temp2 = temp2 + ap(k)*x(i)
258 k = k + 1
259 50 CONTINUE
260 y(j) = y(j) + temp1*ap(kk+j-1) + alpha*temp2
261 kk = kk + j
262 60 CONTINUE
263 ELSE
264 jx = kx
265 jy = ky
266 DO 80 j = 1,n
267 temp1 = alpha*x(jx)
268 temp2 = zero
269 ix = kx
270 iy = ky
271 DO 70 k = kk,kk + j - 2
272 y(iy) = y(iy) + temp1*ap(k)
273 temp2 = temp2 + ap(k)*x(ix)
274 ix = ix + incx
275 iy = iy + incy
276 70 CONTINUE
277 y(jy) = y(jy) + temp1*ap(kk+j-1) + alpha*temp2
278 jx = jx + incx
279 jy = jy + incy
280 kk = kk + j
281 80 CONTINUE
282 END IF
283 ELSE
284*
285* Form y when AP contains the lower triangle.
286*
287 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
288 DO 100 j = 1,n
289 temp1 = alpha*x(j)
290 temp2 = zero
291 y(j) = y(j) + temp1*ap(kk)
292 k = kk + 1
293 DO 90 i = j + 1,n
294 y(i) = y(i) + temp1*ap(k)
295 temp2 = temp2 + ap(k)*x(i)
296 k = k + 1
297 90 CONTINUE
298 y(j) = y(j) + alpha*temp2
299 kk = kk + (n-j+1)
300 100 CONTINUE
301 ELSE
302 jx = kx
303 jy = ky
304 DO 120 j = 1,n
305 temp1 = alpha*x(jx)
306 temp2 = zero
307 y(jy) = y(jy) + temp1*ap(kk)
308 ix = jx
309 iy = jy
310 DO 110 k = kk + 1,kk + n - j
311 ix = ix + incx
312 iy = iy + incy
313 y(iy) = y(iy) + temp1*ap(k)
314 temp2 = temp2 + ap(k)*x(ix)
315 110 CONTINUE
316 y(jy) = y(jy) + alpha*temp2
317 jx = jx + incx
318 jy = jy + incy
319 kk = kk + (n-j+1)
320 120 CONTINUE
321 END IF
322 END IF
323*
324 RETURN
325*
326* End of DSPMV
327*

◆ dspr()

subroutine dspr ( character uplo,
integer n,
double precision alpha,
double precision, dimension(*) x,
integer incx,
double precision, dimension(*) ap )

DSPR

Purpose:
!>
!> DSPR    performs the symmetric rank 1 operation
!>
!>    A := alpha*x*x**T + A,
!>
!> where alpha is a real scalar, x is an n element vector and A is an
!> n by n symmetric matrix, supplied in packed form.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the matrix A is supplied in the packed
!>           array AP as follows:
!>
!>              UPLO = 'U' or 'u'   The upper triangular part of A is
!>                                  supplied in AP.
!>
!>              UPLO = 'L' or 'l'   The lower triangular part of A is
!>                                  supplied in AP.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in,out]AP
!>          AP is DOUBLE PRECISION array, dimension at least
!>           ( ( n*( n + 1 ) )/2 ).
!>           Before entry with  UPLO = 'U' or 'u', the array AP must
!>           contain the upper triangular part of the symmetric matrix
!>           packed sequentially, column by column, so that AP( 1 )
!>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
!>           and a( 2, 2 ) respectively, and so on. On exit, the array
!>           AP is overwritten by the upper triangular part of the
!>           updated matrix.
!>           Before entry with UPLO = 'L' or 'l', the array AP must
!>           contain the lower triangular part of the symmetric matrix
!>           packed sequentially, column by column, so that AP( 1 )
!>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
!>           and a( 3, 1 ) respectively, and so on. On exit, the array
!>           AP is overwritten by the lower triangular part of the
!>           updated matrix.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 126 of file dspr.f.

127*
128* -- Reference BLAS level2 routine --
129* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
130* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
131*
132* .. Scalar Arguments ..
133 DOUBLE PRECISION ALPHA
134 INTEGER INCX,N
135 CHARACTER UPLO
136* ..
137* .. Array Arguments ..
138 DOUBLE PRECISION AP(*),X(*)
139* ..
140*
141* =====================================================================
142*
143* .. Parameters ..
144 DOUBLE PRECISION ZERO
145 parameter(zero=0.0d+0)
146* ..
147* .. Local Scalars ..
148 DOUBLE PRECISION TEMP
149 INTEGER I,INFO,IX,J,JX,K,KK,KX
150* ..
151* .. External Functions ..
152 LOGICAL LSAME
153 EXTERNAL lsame
154* ..
155* .. External Subroutines ..
156 EXTERNAL xerbla
157* ..
158*
159* Test the input parameters.
160*
161 info = 0
162 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
163 info = 1
164 ELSE IF (n.LT.0) THEN
165 info = 2
166 ELSE IF (incx.EQ.0) THEN
167 info = 5
168 END IF
169 IF (info.NE.0) THEN
170 CALL xerbla('DSPR ',info)
171 RETURN
172 END IF
173*
174* Quick return if possible.
175*
176 IF ((n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
177*
178* Set the start point in X if the increment is not unity.
179*
180 IF (incx.LE.0) THEN
181 kx = 1 - (n-1)*incx
182 ELSE IF (incx.NE.1) THEN
183 kx = 1
184 END IF
185*
186* Start the operations. In this version the elements of the array AP
187* are accessed sequentially with one pass through AP.
188*
189 kk = 1
190 IF (lsame(uplo,'U')) THEN
191*
192* Form A when upper triangle is stored in AP.
193*
194 IF (incx.EQ.1) THEN
195 DO 20 j = 1,n
196 IF (x(j).NE.zero) THEN
197 temp = alpha*x(j)
198 k = kk
199 DO 10 i = 1,j
200 ap(k) = ap(k) + x(i)*temp
201 k = k + 1
202 10 CONTINUE
203 END IF
204 kk = kk + j
205 20 CONTINUE
206 ELSE
207 jx = kx
208 DO 40 j = 1,n
209 IF (x(jx).NE.zero) THEN
210 temp = alpha*x(jx)
211 ix = kx
212 DO 30 k = kk,kk + j - 1
213 ap(k) = ap(k) + x(ix)*temp
214 ix = ix + incx
215 30 CONTINUE
216 END IF
217 jx = jx + incx
218 kk = kk + j
219 40 CONTINUE
220 END IF
221 ELSE
222*
223* Form A when lower triangle is stored in AP.
224*
225 IF (incx.EQ.1) THEN
226 DO 60 j = 1,n
227 IF (x(j).NE.zero) THEN
228 temp = alpha*x(j)
229 k = kk
230 DO 50 i = j,n
231 ap(k) = ap(k) + x(i)*temp
232 k = k + 1
233 50 CONTINUE
234 END IF
235 kk = kk + n - j + 1
236 60 CONTINUE
237 ELSE
238 jx = kx
239 DO 80 j = 1,n
240 IF (x(jx).NE.zero) THEN
241 temp = alpha*x(jx)
242 ix = jx
243 DO 70 k = kk,kk + n - j
244 ap(k) = ap(k) + x(ix)*temp
245 ix = ix + incx
246 70 CONTINUE
247 END IF
248 jx = jx + incx
249 kk = kk + n - j + 1
250 80 CONTINUE
251 END IF
252 END IF
253*
254 RETURN
255*
256* End of DSPR
257*

◆ dspr2()

subroutine dspr2 ( character uplo,
integer n,
double precision alpha,
double precision, dimension(*) x,
integer incx,
double precision, dimension(*) y,
integer incy,
double precision, dimension(*) ap )

DSPR2

Purpose:
!>
!> DSPR2  performs the symmetric rank 2 operation
!>
!>    A := alpha*x*y**T + alpha*y*x**T + A,
!>
!> where alpha is a scalar, x and y are n element vectors and A is an
!> n by n symmetric matrix, supplied in packed form.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the matrix A is supplied in the packed
!>           array AP as follows:
!>
!>              UPLO = 'U' or 'u'   The upper triangular part of A is
!>                                  supplied in AP.
!>
!>              UPLO = 'L' or 'l'   The lower triangular part of A is
!>                                  supplied in AP.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ).
!>           Before entry, the incremented array Y must contain the n
!>           element vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
[in,out]AP
!>          AP is DOUBLE PRECISION array, dimension at least
!>           ( ( n*( n + 1 ) )/2 ).
!>           Before entry with  UPLO = 'U' or 'u', the array AP must
!>           contain the upper triangular part of the symmetric matrix
!>           packed sequentially, column by column, so that AP( 1 )
!>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
!>           and a( 2, 2 ) respectively, and so on. On exit, the array
!>           AP is overwritten by the upper triangular part of the
!>           updated matrix.
!>           Before entry with UPLO = 'L' or 'l', the array AP must
!>           contain the lower triangular part of the symmetric matrix
!>           packed sequentially, column by column, so that AP( 1 )
!>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
!>           and a( 3, 1 ) respectively, and so on. On exit, the array
!>           AP is overwritten by the lower triangular part of the
!>           updated matrix.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 141 of file dspr2.f.

142*
143* -- Reference BLAS level2 routine --
144* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
145* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
146*
147* .. Scalar Arguments ..
148 DOUBLE PRECISION ALPHA
149 INTEGER INCX,INCY,N
150 CHARACTER UPLO
151* ..
152* .. Array Arguments ..
153 DOUBLE PRECISION AP(*),X(*),Y(*)
154* ..
155*
156* =====================================================================
157*
158* .. Parameters ..
159 DOUBLE PRECISION ZERO
160 parameter(zero=0.0d+0)
161* ..
162* .. Local Scalars ..
163 DOUBLE PRECISION TEMP1,TEMP2
164 INTEGER I,INFO,IX,IY,J,JX,JY,K,KK,KX,KY
165* ..
166* .. External Functions ..
167 LOGICAL LSAME
168 EXTERNAL lsame
169* ..
170* .. External Subroutines ..
171 EXTERNAL xerbla
172* ..
173*
174* Test the input parameters.
175*
176 info = 0
177 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
178 info = 1
179 ELSE IF (n.LT.0) THEN
180 info = 2
181 ELSE IF (incx.EQ.0) THEN
182 info = 5
183 ELSE IF (incy.EQ.0) THEN
184 info = 7
185 END IF
186 IF (info.NE.0) THEN
187 CALL xerbla('DSPR2 ',info)
188 RETURN
189 END IF
190*
191* Quick return if possible.
192*
193 IF ((n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
194*
195* Set up the start points in X and Y if the increments are not both
196* unity.
197*
198 IF ((incx.NE.1) .OR. (incy.NE.1)) THEN
199 IF (incx.GT.0) THEN
200 kx = 1
201 ELSE
202 kx = 1 - (n-1)*incx
203 END IF
204 IF (incy.GT.0) THEN
205 ky = 1
206 ELSE
207 ky = 1 - (n-1)*incy
208 END IF
209 jx = kx
210 jy = ky
211 END IF
212*
213* Start the operations. In this version the elements of the array AP
214* are accessed sequentially with one pass through AP.
215*
216 kk = 1
217 IF (lsame(uplo,'U')) THEN
218*
219* Form A when upper triangle is stored in AP.
220*
221 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
222 DO 20 j = 1,n
223 IF ((x(j).NE.zero) .OR. (y(j).NE.zero)) THEN
224 temp1 = alpha*y(j)
225 temp2 = alpha*x(j)
226 k = kk
227 DO 10 i = 1,j
228 ap(k) = ap(k) + x(i)*temp1 + y(i)*temp2
229 k = k + 1
230 10 CONTINUE
231 END IF
232 kk = kk + j
233 20 CONTINUE
234 ELSE
235 DO 40 j = 1,n
236 IF ((x(jx).NE.zero) .OR. (y(jy).NE.zero)) THEN
237 temp1 = alpha*y(jy)
238 temp2 = alpha*x(jx)
239 ix = kx
240 iy = ky
241 DO 30 k = kk,kk + j - 1
242 ap(k) = ap(k) + x(ix)*temp1 + y(iy)*temp2
243 ix = ix + incx
244 iy = iy + incy
245 30 CONTINUE
246 END IF
247 jx = jx + incx
248 jy = jy + incy
249 kk = kk + j
250 40 CONTINUE
251 END IF
252 ELSE
253*
254* Form A when lower triangle is stored in AP.
255*
256 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
257 DO 60 j = 1,n
258 IF ((x(j).NE.zero) .OR. (y(j).NE.zero)) THEN
259 temp1 = alpha*y(j)
260 temp2 = alpha*x(j)
261 k = kk
262 DO 50 i = j,n
263 ap(k) = ap(k) + x(i)*temp1 + y(i)*temp2
264 k = k + 1
265 50 CONTINUE
266 END IF
267 kk = kk + n - j + 1
268 60 CONTINUE
269 ELSE
270 DO 80 j = 1,n
271 IF ((x(jx).NE.zero) .OR. (y(jy).NE.zero)) THEN
272 temp1 = alpha*y(jy)
273 temp2 = alpha*x(jx)
274 ix = jx
275 iy = jy
276 DO 70 k = kk,kk + n - j
277 ap(k) = ap(k) + x(ix)*temp1 + y(iy)*temp2
278 ix = ix + incx
279 iy = iy + incy
280 70 CONTINUE
281 END IF
282 jx = jx + incx
283 jy = jy + incy
284 kk = kk + n - j + 1
285 80 CONTINUE
286 END IF
287 END IF
288*
289 RETURN
290*
291* End of DSPR2
292*

◆ dsymv()

subroutine dsymv ( character uplo,
integer n,
double precision alpha,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx,
double precision beta,
double precision, dimension(*) y,
integer incy )

DSYMV

Purpose:
!>
!> DSYMV  performs the matrix-vector  operation
!>
!>    y := alpha*A*x + beta*y,
!>
!> where alpha and beta are scalars, x and y are n element vectors and
!> A is an n by n symmetric matrix.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the array A is to be referenced as
!>           follows:
!>
!>              UPLO = 'U' or 'u'   Only the upper triangular part of A
!>                                  is to be referenced.
!>
!>              UPLO = 'L' or 'l'   Only the lower triangular part of A
!>                                  is to be referenced.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with  UPLO = 'U' or 'u', the leading n by n
!>           upper triangular part of the array A must contain the upper
!>           triangular part of the symmetric matrix and the strictly
!>           lower triangular part of A is not referenced.
!>           Before entry with UPLO = 'L' or 'l', the leading n by n
!>           lower triangular part of the array A must contain the lower
!>           triangular part of the symmetric matrix and the strictly
!>           upper triangular part of A is not referenced.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           max( 1, n ).
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]BETA
!>          BETA is DOUBLE PRECISION.
!>           On entry, BETA specifies the scalar beta. When BETA is
!>           supplied as zero then Y need not be set on input.
!>
[in,out]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ).
!>           Before entry, the incremented array Y must contain the n
!>           element vector y. On exit, Y is overwritten by the updated
!>           vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 151 of file dsymv.f.

152*
153* -- Reference BLAS level2 routine --
154* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
155* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
156*
157* .. Scalar Arguments ..
158 DOUBLE PRECISION ALPHA,BETA
159 INTEGER INCX,INCY,LDA,N
160 CHARACTER UPLO
161* ..
162* .. Array Arguments ..
163 DOUBLE PRECISION A(LDA,*),X(*),Y(*)
164* ..
165*
166* =====================================================================
167*
168* .. Parameters ..
169 DOUBLE PRECISION ONE,ZERO
170 parameter(one=1.0d+0,zero=0.0d+0)
171* ..
172* .. Local Scalars ..
173 DOUBLE PRECISION TEMP1,TEMP2
174 INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY
175* ..
176* .. External Functions ..
177 LOGICAL LSAME
178 EXTERNAL lsame
179* ..
180* .. External Subroutines ..
181 EXTERNAL xerbla
182* ..
183* .. Intrinsic Functions ..
184 INTRINSIC max
185* ..
186*
187* Test the input parameters.
188*
189 info = 0
190 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
191 info = 1
192 ELSE IF (n.LT.0) THEN
193 info = 2
194 ELSE IF (lda.LT.max(1,n)) THEN
195 info = 5
196 ELSE IF (incx.EQ.0) THEN
197 info = 7
198 ELSE IF (incy.EQ.0) THEN
199 info = 10
200 END IF
201 IF (info.NE.0) THEN
202 CALL xerbla('DSYMV ',info)
203 RETURN
204 END IF
205*
206* Quick return if possible.
207*
208 IF ((n.EQ.0) .OR. ((alpha.EQ.zero).AND. (beta.EQ.one))) RETURN
209*
210* Set up the start points in X and Y.
211*
212 IF (incx.GT.0) THEN
213 kx = 1
214 ELSE
215 kx = 1 - (n-1)*incx
216 END IF
217 IF (incy.GT.0) THEN
218 ky = 1
219 ELSE
220 ky = 1 - (n-1)*incy
221 END IF
222*
223* Start the operations. In this version the elements of A are
224* accessed sequentially with one pass through the triangular part
225* of A.
226*
227* First form y := beta*y.
228*
229 IF (beta.NE.one) THEN
230 IF (incy.EQ.1) THEN
231 IF (beta.EQ.zero) THEN
232 DO 10 i = 1,n
233 y(i) = zero
234 10 CONTINUE
235 ELSE
236 DO 20 i = 1,n
237 y(i) = beta*y(i)
238 20 CONTINUE
239 END IF
240 ELSE
241 iy = ky
242 IF (beta.EQ.zero) THEN
243 DO 30 i = 1,n
244 y(iy) = zero
245 iy = iy + incy
246 30 CONTINUE
247 ELSE
248 DO 40 i = 1,n
249 y(iy) = beta*y(iy)
250 iy = iy + incy
251 40 CONTINUE
252 END IF
253 END IF
254 END IF
255 IF (alpha.EQ.zero) RETURN
256 IF (lsame(uplo,'U')) THEN
257*
258* Form y when A is stored in upper triangle.
259*
260 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
261 DO 60 j = 1,n
262 temp1 = alpha*x(j)
263 temp2 = zero
264 DO 50 i = 1,j - 1
265 y(i) = y(i) + temp1*a(i,j)
266 temp2 = temp2 + a(i,j)*x(i)
267 50 CONTINUE
268 y(j) = y(j) + temp1*a(j,j) + alpha*temp2
269 60 CONTINUE
270 ELSE
271 jx = kx
272 jy = ky
273 DO 80 j = 1,n
274 temp1 = alpha*x(jx)
275 temp2 = zero
276 ix = kx
277 iy = ky
278 DO 70 i = 1,j - 1
279 y(iy) = y(iy) + temp1*a(i,j)
280 temp2 = temp2 + a(i,j)*x(ix)
281 ix = ix + incx
282 iy = iy + incy
283 70 CONTINUE
284 y(jy) = y(jy) + temp1*a(j,j) + alpha*temp2
285 jx = jx + incx
286 jy = jy + incy
287 80 CONTINUE
288 END IF
289 ELSE
290*
291* Form y when A is stored in lower triangle.
292*
293 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
294 DO 100 j = 1,n
295 temp1 = alpha*x(j)
296 temp2 = zero
297 y(j) = y(j) + temp1*a(j,j)
298 DO 90 i = j + 1,n
299 y(i) = y(i) + temp1*a(i,j)
300 temp2 = temp2 + a(i,j)*x(i)
301 90 CONTINUE
302 y(j) = y(j) + alpha*temp2
303 100 CONTINUE
304 ELSE
305 jx = kx
306 jy = ky
307 DO 120 j = 1,n
308 temp1 = alpha*x(jx)
309 temp2 = zero
310 y(jy) = y(jy) + temp1*a(j,j)
311 ix = jx
312 iy = jy
313 DO 110 i = j + 1,n
314 ix = ix + incx
315 iy = iy + incy
316 y(iy) = y(iy) + temp1*a(i,j)
317 temp2 = temp2 + a(i,j)*x(ix)
318 110 CONTINUE
319 y(jy) = y(jy) + alpha*temp2
320 jx = jx + incx
321 jy = jy + incy
322 120 CONTINUE
323 END IF
324 END IF
325*
326 RETURN
327*
328* End of DSYMV
329*

◆ dsyr()

subroutine dsyr ( character uplo,
integer n,
double precision alpha,
double precision, dimension(*) x,
integer incx,
double precision, dimension(lda,*) a,
integer lda )

DSYR

Purpose:
!>
!> DSYR   performs the symmetric rank 1 operation
!>
!>    A := alpha*x*x**T + A,
!>
!> where alpha is a real scalar, x is an n element vector and A is an
!> n by n symmetric matrix.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the array A is to be referenced as
!>           follows:
!>
!>              UPLO = 'U' or 'u'   Only the upper triangular part of A
!>                                  is to be referenced.
!>
!>              UPLO = 'L' or 'l'   Only the lower triangular part of A
!>                                  is to be referenced.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in,out]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with  UPLO = 'U' or 'u', the leading n by n
!>           upper triangular part of the array A must contain the upper
!>           triangular part of the symmetric matrix and the strictly
!>           lower triangular part of A is not referenced. On exit, the
!>           upper triangular part of the array A is overwritten by the
!>           upper triangular part of the updated matrix.
!>           Before entry with UPLO = 'L' or 'l', the leading n by n
!>           lower triangular part of the array A must contain the lower
!>           triangular part of the symmetric matrix and the strictly
!>           upper triangular part of A is not referenced. On exit, the
!>           lower triangular part of the array A is overwritten by the
!>           lower triangular part of the updated matrix.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           max( 1, n ).
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 131 of file dsyr.f.

132*
133* -- Reference BLAS level2 routine --
134* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
135* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
136*
137* .. Scalar Arguments ..
138 DOUBLE PRECISION ALPHA
139 INTEGER INCX,LDA,N
140 CHARACTER UPLO
141* ..
142* .. Array Arguments ..
143 DOUBLE PRECISION A(LDA,*),X(*)
144* ..
145*
146* =====================================================================
147*
148* .. Parameters ..
149 DOUBLE PRECISION ZERO
150 parameter(zero=0.0d+0)
151* ..
152* .. Local Scalars ..
153 DOUBLE PRECISION TEMP
154 INTEGER I,INFO,IX,J,JX,KX
155* ..
156* .. External Functions ..
157 LOGICAL LSAME
158 EXTERNAL lsame
159* ..
160* .. External Subroutines ..
161 EXTERNAL xerbla
162* ..
163* .. Intrinsic Functions ..
164 INTRINSIC max
165* ..
166*
167* Test the input parameters.
168*
169 info = 0
170 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
171 info = 1
172 ELSE IF (n.LT.0) THEN
173 info = 2
174 ELSE IF (incx.EQ.0) THEN
175 info = 5
176 ELSE IF (lda.LT.max(1,n)) THEN
177 info = 7
178 END IF
179 IF (info.NE.0) THEN
180 CALL xerbla('DSYR ',info)
181 RETURN
182 END IF
183*
184* Quick return if possible.
185*
186 IF ((n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
187*
188* Set the start point in X if the increment is not unity.
189*
190 IF (incx.LE.0) THEN
191 kx = 1 - (n-1)*incx
192 ELSE IF (incx.NE.1) THEN
193 kx = 1
194 END IF
195*
196* Start the operations. In this version the elements of A are
197* accessed sequentially with one pass through the triangular part
198* of A.
199*
200 IF (lsame(uplo,'U')) THEN
201*
202* Form A when A is stored in upper triangle.
203*
204 IF (incx.EQ.1) THEN
205 DO 20 j = 1,n
206 IF (x(j).NE.zero) THEN
207 temp = alpha*x(j)
208 DO 10 i = 1,j
209 a(i,j) = a(i,j) + x(i)*temp
210 10 CONTINUE
211 END IF
212 20 CONTINUE
213 ELSE
214 jx = kx
215 DO 40 j = 1,n
216 IF (x(jx).NE.zero) THEN
217 temp = alpha*x(jx)
218 ix = kx
219 DO 30 i = 1,j
220 a(i,j) = a(i,j) + x(ix)*temp
221 ix = ix + incx
222 30 CONTINUE
223 END IF
224 jx = jx + incx
225 40 CONTINUE
226 END IF
227 ELSE
228*
229* Form A when A is stored in lower triangle.
230*
231 IF (incx.EQ.1) THEN
232 DO 60 j = 1,n
233 IF (x(j).NE.zero) THEN
234 temp = alpha*x(j)
235 DO 50 i = j,n
236 a(i,j) = a(i,j) + x(i)*temp
237 50 CONTINUE
238 END IF
239 60 CONTINUE
240 ELSE
241 jx = kx
242 DO 80 j = 1,n
243 IF (x(jx).NE.zero) THEN
244 temp = alpha*x(jx)
245 ix = jx
246 DO 70 i = j,n
247 a(i,j) = a(i,j) + x(ix)*temp
248 ix = ix + incx
249 70 CONTINUE
250 END IF
251 jx = jx + incx
252 80 CONTINUE
253 END IF
254 END IF
255*
256 RETURN
257*
258* End of DSYR
259*

◆ dsyr2()

subroutine dsyr2 ( character uplo,
integer n,
double precision alpha,
double precision, dimension(*) x,
integer incx,
double precision, dimension(*) y,
integer incy,
double precision, dimension(lda,*) a,
integer lda )

DSYR2

Purpose:
!>
!> DSYR2  performs the symmetric rank 2 operation
!>
!>    A := alpha*x*y**T + alpha*y*x**T + A,
!>
!> where alpha is a scalar, x and y are n element vectors and A is an n
!> by n symmetric matrix.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the upper or lower
!>           triangular part of the array A is to be referenced as
!>           follows:
!>
!>              UPLO = 'U' or 'u'   Only the upper triangular part of A
!>                                  is to be referenced.
!>
!>              UPLO = 'L' or 'l'   Only the lower triangular part of A
!>                                  is to be referenced.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION.
!>           On entry, ALPHA specifies the scalar alpha.
!>
[in]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
[in]Y
!>          Y is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCY ) ).
!>           Before entry, the incremented array Y must contain the n
!>           element vector y.
!>
[in]INCY
!>          INCY is INTEGER
!>           On entry, INCY specifies the increment for the elements of
!>           Y. INCY must not be zero.
!>
[in,out]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with  UPLO = 'U' or 'u', the leading n by n
!>           upper triangular part of the array A must contain the upper
!>           triangular part of the symmetric matrix and the strictly
!>           lower triangular part of A is not referenced. On exit, the
!>           upper triangular part of the array A is overwritten by the
!>           upper triangular part of the updated matrix.
!>           Before entry with UPLO = 'L' or 'l', the leading n by n
!>           lower triangular part of the array A must contain the lower
!>           triangular part of the symmetric matrix and the strictly
!>           upper triangular part of A is not referenced. On exit, the
!>           lower triangular part of the array A is overwritten by the
!>           lower triangular part of the updated matrix.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           max( 1, n ).
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 146 of file dsyr2.f.

147*
148* -- Reference BLAS level2 routine --
149* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
150* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
151*
152* .. Scalar Arguments ..
153 DOUBLE PRECISION ALPHA
154 INTEGER INCX,INCY,LDA,N
155 CHARACTER UPLO
156* ..
157* .. Array Arguments ..
158 DOUBLE PRECISION A(LDA,*),X(*),Y(*)
159* ..
160*
161* =====================================================================
162*
163* .. Parameters ..
164 DOUBLE PRECISION ZERO
165 parameter(zero=0.0d+0)
166* ..
167* .. Local Scalars ..
168 DOUBLE PRECISION TEMP1,TEMP2
169 INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY
170* ..
171* .. External Functions ..
172 LOGICAL LSAME
173 EXTERNAL lsame
174* ..
175* .. External Subroutines ..
176 EXTERNAL xerbla
177* ..
178* .. Intrinsic Functions ..
179 INTRINSIC max
180* ..
181*
182* Test the input parameters.
183*
184 info = 0
185 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
186 info = 1
187 ELSE IF (n.LT.0) THEN
188 info = 2
189 ELSE IF (incx.EQ.0) THEN
190 info = 5
191 ELSE IF (incy.EQ.0) THEN
192 info = 7
193 ELSE IF (lda.LT.max(1,n)) THEN
194 info = 9
195 END IF
196 IF (info.NE.0) THEN
197 CALL xerbla('DSYR2 ',info)
198 RETURN
199 END IF
200*
201* Quick return if possible.
202*
203 IF ((n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
204*
205* Set up the start points in X and Y if the increments are not both
206* unity.
207*
208 IF ((incx.NE.1) .OR. (incy.NE.1)) THEN
209 IF (incx.GT.0) THEN
210 kx = 1
211 ELSE
212 kx = 1 - (n-1)*incx
213 END IF
214 IF (incy.GT.0) THEN
215 ky = 1
216 ELSE
217 ky = 1 - (n-1)*incy
218 END IF
219 jx = kx
220 jy = ky
221 END IF
222*
223* Start the operations. In this version the elements of A are
224* accessed sequentially with one pass through the triangular part
225* of A.
226*
227 IF (lsame(uplo,'U')) THEN
228*
229* Form A when A is stored in the upper triangle.
230*
231 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
232 DO 20 j = 1,n
233 IF ((x(j).NE.zero) .OR. (y(j).NE.zero)) THEN
234 temp1 = alpha*y(j)
235 temp2 = alpha*x(j)
236 DO 10 i = 1,j
237 a(i,j) = a(i,j) + x(i)*temp1 + y(i)*temp2
238 10 CONTINUE
239 END IF
240 20 CONTINUE
241 ELSE
242 DO 40 j = 1,n
243 IF ((x(jx).NE.zero) .OR. (y(jy).NE.zero)) THEN
244 temp1 = alpha*y(jy)
245 temp2 = alpha*x(jx)
246 ix = kx
247 iy = ky
248 DO 30 i = 1,j
249 a(i,j) = a(i,j) + x(ix)*temp1 + y(iy)*temp2
250 ix = ix + incx
251 iy = iy + incy
252 30 CONTINUE
253 END IF
254 jx = jx + incx
255 jy = jy + incy
256 40 CONTINUE
257 END IF
258 ELSE
259*
260* Form A when A is stored in the lower triangle.
261*
262 IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
263 DO 60 j = 1,n
264 IF ((x(j).NE.zero) .OR. (y(j).NE.zero)) THEN
265 temp1 = alpha*y(j)
266 temp2 = alpha*x(j)
267 DO 50 i = j,n
268 a(i,j) = a(i,j) + x(i)*temp1 + y(i)*temp2
269 50 CONTINUE
270 END IF
271 60 CONTINUE
272 ELSE
273 DO 80 j = 1,n
274 IF ((x(jx).NE.zero) .OR. (y(jy).NE.zero)) THEN
275 temp1 = alpha*y(jy)
276 temp2 = alpha*x(jx)
277 ix = jx
278 iy = jy
279 DO 70 i = j,n
280 a(i,j) = a(i,j) + x(ix)*temp1 + y(iy)*temp2
281 ix = ix + incx
282 iy = iy + incy
283 70 CONTINUE
284 END IF
285 jx = jx + incx
286 jy = jy + incy
287 80 CONTINUE
288 END IF
289 END IF
290*
291 RETURN
292*
293* End of DSYR2
294*

◆ dtbmv()

subroutine dtbmv ( character uplo,
character trans,
character diag,
integer n,
integer k,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx )

DTBMV

Purpose:
!>
!> DTBMV  performs one of the matrix-vector operations
!>
!>    x := A*x,   or   x := A**T*x,
!>
!> where x is an n element vector and  A is an n by n unit, or non-unit,
!> upper or lower triangular band matrix, with ( k + 1 ) diagonals.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the matrix is an upper or
!>           lower triangular matrix as follows:
!>
!>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
!>
!>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the operation to be performed as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   x := A*x.
!>
!>              TRANS = 'T' or 't'   x := A**T*x.
!>
!>              TRANS = 'C' or 'c'   x := A**T*x.
!>
[in]DIAG
!>          DIAG is CHARACTER*1
!>           On entry, DIAG specifies whether or not A is unit
!>           triangular as follows:
!>
!>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
!>
!>              DIAG = 'N' or 'n'   A is not assumed to be unit
!>                                  triangular.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]K
!>          K is INTEGER
!>           On entry with UPLO = 'U' or 'u', K specifies the number of
!>           super-diagonals of the matrix A.
!>           On entry with UPLO = 'L' or 'l', K specifies the number of
!>           sub-diagonals of the matrix A.
!>           K must satisfy  0 .le. K.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
!>           by n part of the array A must contain the upper triangular
!>           band part of the matrix of coefficients, supplied column by
!>           column, with the leading diagonal of the matrix in row
!>           ( k + 1 ) of the array, the first super-diagonal starting at
!>           position 2 in row k, and so on. The top left k by k triangle
!>           of the array A is not referenced.
!>           The following program segment will transfer an upper
!>           triangular band matrix from conventional full matrix storage
!>           to band storage:
!>
!>                 DO 20, J = 1, N
!>                    M = K + 1 - J
!>                    DO 10, I = MAX( 1, J - K ), J
!>                       A( M + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
!>           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
!>           by n part of the array A must contain the lower triangular
!>           band part of the matrix of coefficients, supplied column by
!>           column, with the leading diagonal of the matrix in row 1 of
!>           the array, the first sub-diagonal starting at position 1 in
!>           row 2, and so on. The bottom right k by k triangle of the
!>           array A is not referenced.
!>           The following program segment will transfer a lower
!>           triangular band matrix from conventional full matrix storage
!>           to band storage:
!>
!>                 DO 20, J = 1, N
!>                    M = 1 - J
!>                    DO 10, I = J, MIN( N, J + K )
!>                       A( M + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
!>           Note that when DIAG = 'U' or 'u' the elements of the array A
!>           corresponding to the diagonal elements of the matrix are not
!>           referenced, but are assumed to be unity.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           ( k + 1 ).
!>
[in,out]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x. On exit, X is overwritten with the
!>           transformed vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 185 of file dtbmv.f.

186*
187* -- Reference BLAS level2 routine --
188* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
189* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
190*
191* .. Scalar Arguments ..
192 INTEGER INCX,K,LDA,N
193 CHARACTER DIAG,TRANS,UPLO
194* ..
195* .. Array Arguments ..
196 DOUBLE PRECISION A(LDA,*),X(*)
197* ..
198*
199* =====================================================================
200*
201* .. Parameters ..
202 DOUBLE PRECISION ZERO
203 parameter(zero=0.0d+0)
204* ..
205* .. Local Scalars ..
206 DOUBLE PRECISION TEMP
207 INTEGER I,INFO,IX,J,JX,KPLUS1,KX,L
208 LOGICAL NOUNIT
209* ..
210* .. External Functions ..
211 LOGICAL LSAME
212 EXTERNAL lsame
213* ..
214* .. External Subroutines ..
215 EXTERNAL xerbla
216* ..
217* .. Intrinsic Functions ..
218 INTRINSIC max,min
219* ..
220*
221* Test the input parameters.
222*
223 info = 0
224 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
225 info = 1
226 ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
227 + .NOT.lsame(trans,'C')) THEN
228 info = 2
229 ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
230 info = 3
231 ELSE IF (n.LT.0) THEN
232 info = 4
233 ELSE IF (k.LT.0) THEN
234 info = 5
235 ELSE IF (lda.LT. (k+1)) THEN
236 info = 7
237 ELSE IF (incx.EQ.0) THEN
238 info = 9
239 END IF
240 IF (info.NE.0) THEN
241 CALL xerbla('DTBMV ',info)
242 RETURN
243 END IF
244*
245* Quick return if possible.
246*
247 IF (n.EQ.0) RETURN
248*
249 nounit = lsame(diag,'N')
250*
251* Set up the start point in X if the increment is not unity. This
252* will be ( N - 1 )*INCX too small for descending loops.
253*
254 IF (incx.LE.0) THEN
255 kx = 1 - (n-1)*incx
256 ELSE IF (incx.NE.1) THEN
257 kx = 1
258 END IF
259*
260* Start the operations. In this version the elements of A are
261* accessed sequentially with one pass through A.
262*
263 IF (lsame(trans,'N')) THEN
264*
265* Form x := A*x.
266*
267 IF (lsame(uplo,'U')) THEN
268 kplus1 = k + 1
269 IF (incx.EQ.1) THEN
270 DO 20 j = 1,n
271 IF (x(j).NE.zero) THEN
272 temp = x(j)
273 l = kplus1 - j
274 DO 10 i = max(1,j-k),j - 1
275 x(i) = x(i) + temp*a(l+i,j)
276 10 CONTINUE
277 IF (nounit) x(j) = x(j)*a(kplus1,j)
278 END IF
279 20 CONTINUE
280 ELSE
281 jx = kx
282 DO 40 j = 1,n
283 IF (x(jx).NE.zero) THEN
284 temp = x(jx)
285 ix = kx
286 l = kplus1 - j
287 DO 30 i = max(1,j-k),j - 1
288 x(ix) = x(ix) + temp*a(l+i,j)
289 ix = ix + incx
290 30 CONTINUE
291 IF (nounit) x(jx) = x(jx)*a(kplus1,j)
292 END IF
293 jx = jx + incx
294 IF (j.GT.k) kx = kx + incx
295 40 CONTINUE
296 END IF
297 ELSE
298 IF (incx.EQ.1) THEN
299 DO 60 j = n,1,-1
300 IF (x(j).NE.zero) THEN
301 temp = x(j)
302 l = 1 - j
303 DO 50 i = min(n,j+k),j + 1,-1
304 x(i) = x(i) + temp*a(l+i,j)
305 50 CONTINUE
306 IF (nounit) x(j) = x(j)*a(1,j)
307 END IF
308 60 CONTINUE
309 ELSE
310 kx = kx + (n-1)*incx
311 jx = kx
312 DO 80 j = n,1,-1
313 IF (x(jx).NE.zero) THEN
314 temp = x(jx)
315 ix = kx
316 l = 1 - j
317 DO 70 i = min(n,j+k),j + 1,-1
318 x(ix) = x(ix) + temp*a(l+i,j)
319 ix = ix - incx
320 70 CONTINUE
321 IF (nounit) x(jx) = x(jx)*a(1,j)
322 END IF
323 jx = jx - incx
324 IF ((n-j).GE.k) kx = kx - incx
325 80 CONTINUE
326 END IF
327 END IF
328 ELSE
329*
330* Form x := A**T*x.
331*
332 IF (lsame(uplo,'U')) THEN
333 kplus1 = k + 1
334 IF (incx.EQ.1) THEN
335 DO 100 j = n,1,-1
336 temp = x(j)
337 l = kplus1 - j
338 IF (nounit) temp = temp*a(kplus1,j)
339 DO 90 i = j - 1,max(1,j-k),-1
340 temp = temp + a(l+i,j)*x(i)
341 90 CONTINUE
342 x(j) = temp
343 100 CONTINUE
344 ELSE
345 kx = kx + (n-1)*incx
346 jx = kx
347 DO 120 j = n,1,-1
348 temp = x(jx)
349 kx = kx - incx
350 ix = kx
351 l = kplus1 - j
352 IF (nounit) temp = temp*a(kplus1,j)
353 DO 110 i = j - 1,max(1,j-k),-1
354 temp = temp + a(l+i,j)*x(ix)
355 ix = ix - incx
356 110 CONTINUE
357 x(jx) = temp
358 jx = jx - incx
359 120 CONTINUE
360 END IF
361 ELSE
362 IF (incx.EQ.1) THEN
363 DO 140 j = 1,n
364 temp = x(j)
365 l = 1 - j
366 IF (nounit) temp = temp*a(1,j)
367 DO 130 i = j + 1,min(n,j+k)
368 temp = temp + a(l+i,j)*x(i)
369 130 CONTINUE
370 x(j) = temp
371 140 CONTINUE
372 ELSE
373 jx = kx
374 DO 160 j = 1,n
375 temp = x(jx)
376 kx = kx + incx
377 ix = kx
378 l = 1 - j
379 IF (nounit) temp = temp*a(1,j)
380 DO 150 i = j + 1,min(n,j+k)
381 temp = temp + a(l+i,j)*x(ix)
382 ix = ix + incx
383 150 CONTINUE
384 x(jx) = temp
385 jx = jx + incx
386 160 CONTINUE
387 END IF
388 END IF
389 END IF
390*
391 RETURN
392*
393* End of DTBMV
394*

◆ dtbsv()

subroutine dtbsv ( character uplo,
character trans,
character diag,
integer n,
integer k,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx )

DTBSV

Purpose:
!>
!> DTBSV  solves one of the systems of equations
!>
!>    A*x = b,   or   A**T*x = b,
!>
!> where b and x are n element vectors and A is an n by n unit, or
!> non-unit, upper or lower triangular band matrix, with ( k + 1 )
!> diagonals.
!>
!> No test for singularity or near-singularity is included in this
!> routine. Such tests must be performed before calling this routine.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the matrix is an upper or
!>           lower triangular matrix as follows:
!>
!>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
!>
!>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the equations to be solved as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   A*x = b.
!>
!>              TRANS = 'T' or 't'   A**T*x = b.
!>
!>              TRANS = 'C' or 'c'   A**T*x = b.
!>
[in]DIAG
!>          DIAG is CHARACTER*1
!>           On entry, DIAG specifies whether or not A is unit
!>           triangular as follows:
!>
!>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
!>
!>              DIAG = 'N' or 'n'   A is not assumed to be unit
!>                                  triangular.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]K
!>          K is INTEGER
!>           On entry with UPLO = 'U' or 'u', K specifies the number of
!>           super-diagonals of the matrix A.
!>           On entry with UPLO = 'L' or 'l', K specifies the number of
!>           sub-diagonals of the matrix A.
!>           K must satisfy  0 .le. K.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
!>           by n part of the array A must contain the upper triangular
!>           band part of the matrix of coefficients, supplied column by
!>           column, with the leading diagonal of the matrix in row
!>           ( k + 1 ) of the array, the first super-diagonal starting at
!>           position 2 in row k, and so on. The top left k by k triangle
!>           of the array A is not referenced.
!>           The following program segment will transfer an upper
!>           triangular band matrix from conventional full matrix storage
!>           to band storage:
!>
!>                 DO 20, J = 1, N
!>                    M = K + 1 - J
!>                    DO 10, I = MAX( 1, J - K ), J
!>                       A( M + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
!>           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
!>           by n part of the array A must contain the lower triangular
!>           band part of the matrix of coefficients, supplied column by
!>           column, with the leading diagonal of the matrix in row 1 of
!>           the array, the first sub-diagonal starting at position 1 in
!>           row 2, and so on. The bottom right k by k triangle of the
!>           array A is not referenced.
!>           The following program segment will transfer a lower
!>           triangular band matrix from conventional full matrix storage
!>           to band storage:
!>
!>                 DO 20, J = 1, N
!>                    M = 1 - J
!>                    DO 10, I = J, MIN( N, J + K )
!>                       A( M + I, J ) = matrix( I, J )
!>              10    CONTINUE
!>              20 CONTINUE
!>
!>           Note that when DIAG = 'U' or 'u' the elements of the array A
!>           corresponding to the diagonal elements of the matrix are not
!>           referenced, but are assumed to be unity.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           ( k + 1 ).
!>
[in,out]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element right-hand side vector b. On exit, X is overwritten
!>           with the solution vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 188 of file dtbsv.f.

189*
190* -- Reference BLAS level2 routine --
191* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
192* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
193*
194* .. Scalar Arguments ..
195 INTEGER INCX,K,LDA,N
196 CHARACTER DIAG,TRANS,UPLO
197* ..
198* .. Array Arguments ..
199 DOUBLE PRECISION A(LDA,*),X(*)
200* ..
201*
202* =====================================================================
203*
204* .. Parameters ..
205 DOUBLE PRECISION ZERO
206 parameter(zero=0.0d+0)
207* ..
208* .. Local Scalars ..
209 DOUBLE PRECISION TEMP
210 INTEGER I,INFO,IX,J,JX,KPLUS1,KX,L
211 LOGICAL NOUNIT
212* ..
213* .. External Functions ..
214 LOGICAL LSAME
215 EXTERNAL lsame
216* ..
217* .. External Subroutines ..
218 EXTERNAL xerbla
219* ..
220* .. Intrinsic Functions ..
221 INTRINSIC max,min
222* ..
223*
224* Test the input parameters.
225*
226 info = 0
227 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
228 info = 1
229 ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
230 + .NOT.lsame(trans,'C')) THEN
231 info = 2
232 ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
233 info = 3
234 ELSE IF (n.LT.0) THEN
235 info = 4
236 ELSE IF (k.LT.0) THEN
237 info = 5
238 ELSE IF (lda.LT. (k+1)) THEN
239 info = 7
240 ELSE IF (incx.EQ.0) THEN
241 info = 9
242 END IF
243 IF (info.NE.0) THEN
244 CALL xerbla('DTBSV ',info)
245 RETURN
246 END IF
247*
248* Quick return if possible.
249*
250 IF (n.EQ.0) RETURN
251*
252 nounit = lsame(diag,'N')
253*
254* Set up the start point in X if the increment is not unity. This
255* will be ( N - 1 )*INCX too small for descending loops.
256*
257 IF (incx.LE.0) THEN
258 kx = 1 - (n-1)*incx
259 ELSE IF (incx.NE.1) THEN
260 kx = 1
261 END IF
262*
263* Start the operations. In this version the elements of A are
264* accessed by sequentially with one pass through A.
265*
266 IF (lsame(trans,'N')) THEN
267*
268* Form x := inv( A )*x.
269*
270 IF (lsame(uplo,'U')) THEN
271 kplus1 = k + 1
272 IF (incx.EQ.1) THEN
273 DO 20 j = n,1,-1
274 IF (x(j).NE.zero) THEN
275 l = kplus1 - j
276 IF (nounit) x(j) = x(j)/a(kplus1,j)
277 temp = x(j)
278 DO 10 i = j - 1,max(1,j-k),-1
279 x(i) = x(i) - temp*a(l+i,j)
280 10 CONTINUE
281 END IF
282 20 CONTINUE
283 ELSE
284 kx = kx + (n-1)*incx
285 jx = kx
286 DO 40 j = n,1,-1
287 kx = kx - incx
288 IF (x(jx).NE.zero) THEN
289 ix = kx
290 l = kplus1 - j
291 IF (nounit) x(jx) = x(jx)/a(kplus1,j)
292 temp = x(jx)
293 DO 30 i = j - 1,max(1,j-k),-1
294 x(ix) = x(ix) - temp*a(l+i,j)
295 ix = ix - incx
296 30 CONTINUE
297 END IF
298 jx = jx - incx
299 40 CONTINUE
300 END IF
301 ELSE
302 IF (incx.EQ.1) THEN
303 DO 60 j = 1,n
304 IF (x(j).NE.zero) THEN
305 l = 1 - j
306 IF (nounit) x(j) = x(j)/a(1,j)
307 temp = x(j)
308 DO 50 i = j + 1,min(n,j+k)
309 x(i) = x(i) - temp*a(l+i,j)
310 50 CONTINUE
311 END IF
312 60 CONTINUE
313 ELSE
314 jx = kx
315 DO 80 j = 1,n
316 kx = kx + incx
317 IF (x(jx).NE.zero) THEN
318 ix = kx
319 l = 1 - j
320 IF (nounit) x(jx) = x(jx)/a(1,j)
321 temp = x(jx)
322 DO 70 i = j + 1,min(n,j+k)
323 x(ix) = x(ix) - temp*a(l+i,j)
324 ix = ix + incx
325 70 CONTINUE
326 END IF
327 jx = jx + incx
328 80 CONTINUE
329 END IF
330 END IF
331 ELSE
332*
333* Form x := inv( A**T)*x.
334*
335 IF (lsame(uplo,'U')) THEN
336 kplus1 = k + 1
337 IF (incx.EQ.1) THEN
338 DO 100 j = 1,n
339 temp = x(j)
340 l = kplus1 - j
341 DO 90 i = max(1,j-k),j - 1
342 temp = temp - a(l+i,j)*x(i)
343 90 CONTINUE
344 IF (nounit) temp = temp/a(kplus1,j)
345 x(j) = temp
346 100 CONTINUE
347 ELSE
348 jx = kx
349 DO 120 j = 1,n
350 temp = x(jx)
351 ix = kx
352 l = kplus1 - j
353 DO 110 i = max(1,j-k),j - 1
354 temp = temp - a(l+i,j)*x(ix)
355 ix = ix + incx
356 110 CONTINUE
357 IF (nounit) temp = temp/a(kplus1,j)
358 x(jx) = temp
359 jx = jx + incx
360 IF (j.GT.k) kx = kx + incx
361 120 CONTINUE
362 END IF
363 ELSE
364 IF (incx.EQ.1) THEN
365 DO 140 j = n,1,-1
366 temp = x(j)
367 l = 1 - j
368 DO 130 i = min(n,j+k),j + 1,-1
369 temp = temp - a(l+i,j)*x(i)
370 130 CONTINUE
371 IF (nounit) temp = temp/a(1,j)
372 x(j) = temp
373 140 CONTINUE
374 ELSE
375 kx = kx + (n-1)*incx
376 jx = kx
377 DO 160 j = n,1,-1
378 temp = x(jx)
379 ix = kx
380 l = 1 - j
381 DO 150 i = min(n,j+k),j + 1,-1
382 temp = temp - a(l+i,j)*x(ix)
383 ix = ix - incx
384 150 CONTINUE
385 IF (nounit) temp = temp/a(1,j)
386 x(jx) = temp
387 jx = jx - incx
388 IF ((n-j).GE.k) kx = kx - incx
389 160 CONTINUE
390 END IF
391 END IF
392 END IF
393*
394 RETURN
395*
396* End of DTBSV
397*

◆ dtpmv()

subroutine dtpmv ( character uplo,
character trans,
character diag,
integer n,
double precision, dimension(*) ap,
double precision, dimension(*) x,
integer incx )

DTPMV

Purpose:
!>
!> DTPMV  performs one of the matrix-vector operations
!>
!>    x := A*x,   or   x := A**T*x,
!>
!> where x is an n element vector and  A is an n by n unit, or non-unit,
!> upper or lower triangular matrix, supplied in packed form.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the matrix is an upper or
!>           lower triangular matrix as follows:
!>
!>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
!>
!>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the operation to be performed as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   x := A*x.
!>
!>              TRANS = 'T' or 't'   x := A**T*x.
!>
!>              TRANS = 'C' or 'c'   x := A**T*x.
!>
[in]DIAG
!>          DIAG is CHARACTER*1
!>           On entry, DIAG specifies whether or not A is unit
!>           triangular as follows:
!>
!>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
!>
!>              DIAG = 'N' or 'n'   A is not assumed to be unit
!>                                  triangular.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]AP
!>          AP is DOUBLE PRECISION array, dimension at least
!>           ( ( n*( n + 1 ) )/2 ).
!>           Before entry with  UPLO = 'U' or 'u', the array AP must
!>           contain the upper triangular matrix packed sequentially,
!>           column by column, so that AP( 1 ) contains a( 1, 1 ),
!>           AP( 2 ) and AP( 3 ) contain a( 1, 2 ) and a( 2, 2 )
!>           respectively, and so on.
!>           Before entry with UPLO = 'L' or 'l', the array AP must
!>           contain the lower triangular matrix packed sequentially,
!>           column by column, so that AP( 1 ) contains a( 1, 1 ),
!>           AP( 2 ) and AP( 3 ) contain a( 2, 1 ) and a( 3, 1 )
!>           respectively, and so on.
!>           Note that when  DIAG = 'U' or 'u', the diagonal elements of
!>           A are not referenced, but are assumed to be unity.
!>
[in,out]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x. On exit, X is overwritten with the
!>           transformed vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 141 of file dtpmv.f.

142*
143* -- Reference BLAS level2 routine --
144* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
145* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
146*
147* .. Scalar Arguments ..
148 INTEGER INCX,N
149 CHARACTER DIAG,TRANS,UPLO
150* ..
151* .. Array Arguments ..
152 DOUBLE PRECISION AP(*),X(*)
153* ..
154*
155* =====================================================================
156*
157* .. Parameters ..
158 DOUBLE PRECISION ZERO
159 parameter(zero=0.0d+0)
160* ..
161* .. Local Scalars ..
162 DOUBLE PRECISION TEMP
163 INTEGER I,INFO,IX,J,JX,K,KK,KX
164 LOGICAL NOUNIT
165* ..
166* .. External Functions ..
167 LOGICAL LSAME
168 EXTERNAL lsame
169* ..
170* .. External Subroutines ..
171 EXTERNAL xerbla
172* ..
173*
174* Test the input parameters.
175*
176 info = 0
177 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
178 info = 1
179 ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
180 + .NOT.lsame(trans,'C')) THEN
181 info = 2
182 ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
183 info = 3
184 ELSE IF (n.LT.0) THEN
185 info = 4
186 ELSE IF (incx.EQ.0) THEN
187 info = 7
188 END IF
189 IF (info.NE.0) THEN
190 CALL xerbla('DTPMV ',info)
191 RETURN
192 END IF
193*
194* Quick return if possible.
195*
196 IF (n.EQ.0) RETURN
197*
198 nounit = lsame(diag,'N')
199*
200* Set up the start point in X if the increment is not unity. This
201* will be ( N - 1 )*INCX too small for descending loops.
202*
203 IF (incx.LE.0) THEN
204 kx = 1 - (n-1)*incx
205 ELSE IF (incx.NE.1) THEN
206 kx = 1
207 END IF
208*
209* Start the operations. In this version the elements of AP are
210* accessed sequentially with one pass through AP.
211*
212 IF (lsame(trans,'N')) THEN
213*
214* Form x:= A*x.
215*
216 IF (lsame(uplo,'U')) THEN
217 kk = 1
218 IF (incx.EQ.1) THEN
219 DO 20 j = 1,n
220 IF (x(j).NE.zero) THEN
221 temp = x(j)
222 k = kk
223 DO 10 i = 1,j - 1
224 x(i) = x(i) + temp*ap(k)
225 k = k + 1
226 10 CONTINUE
227 IF (nounit) x(j) = x(j)*ap(kk+j-1)
228 END IF
229 kk = kk + j
230 20 CONTINUE
231 ELSE
232 jx = kx
233 DO 40 j = 1,n
234 IF (x(jx).NE.zero) THEN
235 temp = x(jx)
236 ix = kx
237 DO 30 k = kk,kk + j - 2
238 x(ix) = x(ix) + temp*ap(k)
239 ix = ix + incx
240 30 CONTINUE
241 IF (nounit) x(jx) = x(jx)*ap(kk+j-1)
242 END IF
243 jx = jx + incx
244 kk = kk + j
245 40 CONTINUE
246 END IF
247 ELSE
248 kk = (n* (n+1))/2
249 IF (incx.EQ.1) THEN
250 DO 60 j = n,1,-1
251 IF (x(j).NE.zero) THEN
252 temp = x(j)
253 k = kk
254 DO 50 i = n,j + 1,-1
255 x(i) = x(i) + temp*ap(k)
256 k = k - 1
257 50 CONTINUE
258 IF (nounit) x(j) = x(j)*ap(kk-n+j)
259 END IF
260 kk = kk - (n-j+1)
261 60 CONTINUE
262 ELSE
263 kx = kx + (n-1)*incx
264 jx = kx
265 DO 80 j = n,1,-1
266 IF (x(jx).NE.zero) THEN
267 temp = x(jx)
268 ix = kx
269 DO 70 k = kk,kk - (n- (j+1)),-1
270 x(ix) = x(ix) + temp*ap(k)
271 ix = ix - incx
272 70 CONTINUE
273 IF (nounit) x(jx) = x(jx)*ap(kk-n+j)
274 END IF
275 jx = jx - incx
276 kk = kk - (n-j+1)
277 80 CONTINUE
278 END IF
279 END IF
280 ELSE
281*
282* Form x := A**T*x.
283*
284 IF (lsame(uplo,'U')) THEN
285 kk = (n* (n+1))/2
286 IF (incx.EQ.1) THEN
287 DO 100 j = n,1,-1
288 temp = x(j)
289 IF (nounit) temp = temp*ap(kk)
290 k = kk - 1
291 DO 90 i = j - 1,1,-1
292 temp = temp + ap(k)*x(i)
293 k = k - 1
294 90 CONTINUE
295 x(j) = temp
296 kk = kk - j
297 100 CONTINUE
298 ELSE
299 jx = kx + (n-1)*incx
300 DO 120 j = n,1,-1
301 temp = x(jx)
302 ix = jx
303 IF (nounit) temp = temp*ap(kk)
304 DO 110 k = kk - 1,kk - j + 1,-1
305 ix = ix - incx
306 temp = temp + ap(k)*x(ix)
307 110 CONTINUE
308 x(jx) = temp
309 jx = jx - incx
310 kk = kk - j
311 120 CONTINUE
312 END IF
313 ELSE
314 kk = 1
315 IF (incx.EQ.1) THEN
316 DO 140 j = 1,n
317 temp = x(j)
318 IF (nounit) temp = temp*ap(kk)
319 k = kk + 1
320 DO 130 i = j + 1,n
321 temp = temp + ap(k)*x(i)
322 k = k + 1
323 130 CONTINUE
324 x(j) = temp
325 kk = kk + (n-j+1)
326 140 CONTINUE
327 ELSE
328 jx = kx
329 DO 160 j = 1,n
330 temp = x(jx)
331 ix = jx
332 IF (nounit) temp = temp*ap(kk)
333 DO 150 k = kk + 1,kk + n - j
334 ix = ix + incx
335 temp = temp + ap(k)*x(ix)
336 150 CONTINUE
337 x(jx) = temp
338 jx = jx + incx
339 kk = kk + (n-j+1)
340 160 CONTINUE
341 END IF
342 END IF
343 END IF
344*
345 RETURN
346*
347* End of DTPMV
348*

◆ dtpsv()

subroutine dtpsv ( character uplo,
character trans,
character diag,
integer n,
double precision, dimension(*) ap,
double precision, dimension(*) x,
integer incx )

DTPSV

Purpose:
!>
!> DTPSV  solves one of the systems of equations
!>
!>    A*x = b,   or   A**T*x = b,
!>
!> where b and x are n element vectors and A is an n by n unit, or
!> non-unit, upper or lower triangular matrix, supplied in packed form.
!>
!> No test for singularity or near-singularity is included in this
!> routine. Such tests must be performed before calling this routine.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the matrix is an upper or
!>           lower triangular matrix as follows:
!>
!>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
!>
!>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the equations to be solved as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   A*x = b.
!>
!>              TRANS = 'T' or 't'   A**T*x = b.
!>
!>              TRANS = 'C' or 'c'   A**T*x = b.
!>
[in]DIAG
!>          DIAG is CHARACTER*1
!>           On entry, DIAG specifies whether or not A is unit
!>           triangular as follows:
!>
!>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
!>
!>              DIAG = 'N' or 'n'   A is not assumed to be unit
!>                                  triangular.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]AP
!>          AP is DOUBLE PRECISION array, dimension at least
!>           ( ( n*( n + 1 ) )/2 ).
!>           Before entry with  UPLO = 'U' or 'u', the array AP must
!>           contain the upper triangular matrix packed sequentially,
!>           column by column, so that AP( 1 ) contains a( 1, 1 ),
!>           AP( 2 ) and AP( 3 ) contain a( 1, 2 ) and a( 2, 2 )
!>           respectively, and so on.
!>           Before entry with UPLO = 'L' or 'l', the array AP must
!>           contain the lower triangular matrix packed sequentially,
!>           column by column, so that AP( 1 ) contains a( 1, 1 ),
!>           AP( 2 ) and AP( 3 ) contain a( 2, 1 ) and a( 3, 1 )
!>           respectively, and so on.
!>           Note that when  DIAG = 'U' or 'u', the diagonal elements of
!>           A are not referenced, but are assumed to be unity.
!>
[in,out]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element right-hand side vector b. On exit, X is overwritten
!>           with the solution vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 143 of file dtpsv.f.

144*
145* -- Reference BLAS level2 routine --
146* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
147* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
148*
149* .. Scalar Arguments ..
150 INTEGER INCX,N
151 CHARACTER DIAG,TRANS,UPLO
152* ..
153* .. Array Arguments ..
154 DOUBLE PRECISION AP(*),X(*)
155* ..
156*
157* =====================================================================
158*
159* .. Parameters ..
160 DOUBLE PRECISION ZERO
161 parameter(zero=0.0d+0)
162* ..
163* .. Local Scalars ..
164 DOUBLE PRECISION TEMP
165 INTEGER I,INFO,IX,J,JX,K,KK,KX
166 LOGICAL NOUNIT
167* ..
168* .. External Functions ..
169 LOGICAL LSAME
170 EXTERNAL lsame
171* ..
172* .. External Subroutines ..
173 EXTERNAL xerbla
174* ..
175*
176* Test the input parameters.
177*
178 info = 0
179 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
180 info = 1
181 ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
182 + .NOT.lsame(trans,'C')) THEN
183 info = 2
184 ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
185 info = 3
186 ELSE IF (n.LT.0) THEN
187 info = 4
188 ELSE IF (incx.EQ.0) THEN
189 info = 7
190 END IF
191 IF (info.NE.0) THEN
192 CALL xerbla('DTPSV ',info)
193 RETURN
194 END IF
195*
196* Quick return if possible.
197*
198 IF (n.EQ.0) RETURN
199*
200 nounit = lsame(diag,'N')
201*
202* Set up the start point in X if the increment is not unity. This
203* will be ( N - 1 )*INCX too small for descending loops.
204*
205 IF (incx.LE.0) THEN
206 kx = 1 - (n-1)*incx
207 ELSE IF (incx.NE.1) THEN
208 kx = 1
209 END IF
210*
211* Start the operations. In this version the elements of AP are
212* accessed sequentially with one pass through AP.
213*
214 IF (lsame(trans,'N')) THEN
215*
216* Form x := inv( A )*x.
217*
218 IF (lsame(uplo,'U')) THEN
219 kk = (n* (n+1))/2
220 IF (incx.EQ.1) THEN
221 DO 20 j = n,1,-1
222 IF (x(j).NE.zero) THEN
223 IF (nounit) x(j) = x(j)/ap(kk)
224 temp = x(j)
225 k = kk - 1
226 DO 10 i = j - 1,1,-1
227 x(i) = x(i) - temp*ap(k)
228 k = k - 1
229 10 CONTINUE
230 END IF
231 kk = kk - j
232 20 CONTINUE
233 ELSE
234 jx = kx + (n-1)*incx
235 DO 40 j = n,1,-1
236 IF (x(jx).NE.zero) THEN
237 IF (nounit) x(jx) = x(jx)/ap(kk)
238 temp = x(jx)
239 ix = jx
240 DO 30 k = kk - 1,kk - j + 1,-1
241 ix = ix - incx
242 x(ix) = x(ix) - temp*ap(k)
243 30 CONTINUE
244 END IF
245 jx = jx - incx
246 kk = kk - j
247 40 CONTINUE
248 END IF
249 ELSE
250 kk = 1
251 IF (incx.EQ.1) THEN
252 DO 60 j = 1,n
253 IF (x(j).NE.zero) THEN
254 IF (nounit) x(j) = x(j)/ap(kk)
255 temp = x(j)
256 k = kk + 1
257 DO 50 i = j + 1,n
258 x(i) = x(i) - temp*ap(k)
259 k = k + 1
260 50 CONTINUE
261 END IF
262 kk = kk + (n-j+1)
263 60 CONTINUE
264 ELSE
265 jx = kx
266 DO 80 j = 1,n
267 IF (x(jx).NE.zero) THEN
268 IF (nounit) x(jx) = x(jx)/ap(kk)
269 temp = x(jx)
270 ix = jx
271 DO 70 k = kk + 1,kk + n - j
272 ix = ix + incx
273 x(ix) = x(ix) - temp*ap(k)
274 70 CONTINUE
275 END IF
276 jx = jx + incx
277 kk = kk + (n-j+1)
278 80 CONTINUE
279 END IF
280 END IF
281 ELSE
282*
283* Form x := inv( A**T )*x.
284*
285 IF (lsame(uplo,'U')) THEN
286 kk = 1
287 IF (incx.EQ.1) THEN
288 DO 100 j = 1,n
289 temp = x(j)
290 k = kk
291 DO 90 i = 1,j - 1
292 temp = temp - ap(k)*x(i)
293 k = k + 1
294 90 CONTINUE
295 IF (nounit) temp = temp/ap(kk+j-1)
296 x(j) = temp
297 kk = kk + j
298 100 CONTINUE
299 ELSE
300 jx = kx
301 DO 120 j = 1,n
302 temp = x(jx)
303 ix = kx
304 DO 110 k = kk,kk + j - 2
305 temp = temp - ap(k)*x(ix)
306 ix = ix + incx
307 110 CONTINUE
308 IF (nounit) temp = temp/ap(kk+j-1)
309 x(jx) = temp
310 jx = jx + incx
311 kk = kk + j
312 120 CONTINUE
313 END IF
314 ELSE
315 kk = (n* (n+1))/2
316 IF (incx.EQ.1) THEN
317 DO 140 j = n,1,-1
318 temp = x(j)
319 k = kk
320 DO 130 i = n,j + 1,-1
321 temp = temp - ap(k)*x(i)
322 k = k - 1
323 130 CONTINUE
324 IF (nounit) temp = temp/ap(kk-n+j)
325 x(j) = temp
326 kk = kk - (n-j+1)
327 140 CONTINUE
328 ELSE
329 kx = kx + (n-1)*incx
330 jx = kx
331 DO 160 j = n,1,-1
332 temp = x(jx)
333 ix = kx
334 DO 150 k = kk,kk - (n- (j+1)),-1
335 temp = temp - ap(k)*x(ix)
336 ix = ix - incx
337 150 CONTINUE
338 IF (nounit) temp = temp/ap(kk-n+j)
339 x(jx) = temp
340 jx = jx - incx
341 kk = kk - (n-j+1)
342 160 CONTINUE
343 END IF
344 END IF
345 END IF
346*
347 RETURN
348*
349* End of DTPSV
350*

◆ dtrmv()

subroutine dtrmv ( character uplo,
character trans,
character diag,
integer n,
double precision, dimension(lda,*) a,
integer lda,
double precision, dimension(*) x,
integer incx )

DTRMV

Purpose:
!>
!> DTRMV  performs one of the matrix-vector operations
!>
!>    x := A*x,   or   x := A**T*x,
!>
!> where x is an n element vector and  A is an n by n unit, or non-unit,
!> upper or lower triangular matrix.
!>
Parameters
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the matrix is an upper or
!>           lower triangular matrix as follows:
!>
!>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
!>
!>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS specifies the operation to be performed as
!>           follows:
!>
!>              TRANS = 'N' or 'n'   x := A*x.
!>
!>              TRANS = 'T' or 't'   x := A**T*x.
!>
!>              TRANS = 'C' or 'c'   x := A**T*x.
!>
[in]DIAG
!>          DIAG is CHARACTER*1
!>           On entry, DIAG specifies whether or not A is unit
!>           triangular as follows:
!>
!>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
!>
!>              DIAG = 'N' or 'n'   A is not assumed to be unit
!>                                  triangular.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the order of the matrix A.
!>           N must be at least zero.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension ( LDA, N )
!>           Before entry with  UPLO = 'U' or 'u', the leading n by n
!>           upper triangular part of the array A must contain the upper
!>           triangular matrix and the strictly lower triangular part of
!>           A is not referenced.
!>           Before entry with UPLO = 'L' or 'l', the leading n by n
!>           lower triangular part of the array A must contain the lower
!>           triangular matrix and the strictly upper triangular part of
!>           A is not referenced.
!>           Note that when  DIAG = 'U' or 'u', the diagonal elements of
!>           A are not referenced either, but are assumed to be unity.
!>
[in]LDA
!>          LDA is INTEGER
!>           On entry, LDA specifies the first dimension of A as declared
!>           in the calling (sub) program. LDA must be at least
!>           max( 1, n ).
!>
[in,out]X
!>          X is DOUBLE PRECISION array, dimension at least
!>           ( 1 + ( n - 1 )*abs( INCX ) ).
!>           Before entry, the incremented array X must contain the n
!>           element vector x. On exit, X is overwritten with the
!>           transformed vector x.
!>
[in]INCX
!>          INCX is INTEGER
!>           On entry, INCX specifies the increment for the elements of
!>           X. INCX must not be zero.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  Level 2 Blas routine.
!>  The vector and matrix arguments are not referenced when N = 0, or M = 0
!>
!>  -- Written on 22-October-1986.
!>     Jack Dongarra, Argonne National Lab.
!>     Jeremy Du Croz, Nag Central Office.
!>     Sven Hammarling, Nag Central Office.
!>     Richard Hanson, Sandia National Labs.
!>

Definition at line 146 of file dtrmv.f.

147*
148* -- Reference BLAS level2 routine --
149* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
150* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
151*
152* .. Scalar Arguments ..
153 INTEGER INCX,LDA,N
154 CHARACTER DIAG,TRANS,UPLO
155* ..
156* .. Array Arguments ..
157 DOUBLE PRECISION A(LDA,*),X(*)
158* ..
159*
160* =====================================================================
161*
162* .. Parameters ..
163 DOUBLE PRECISION ZERO
164 parameter(zero=0.0d+0)
165* ..
166* .. Local Scalars ..
167 DOUBLE PRECISION TEMP
168 INTEGER I,INFO,IX,J,JX,KX
169 LOGICAL NOUNIT
170* ..
171* .. External Functions ..
172 LOGICAL LSAME
173 EXTERNAL lsame
174* ..
175* .. External Subroutines ..
176 EXTERNAL xerbla
177* ..
178* .. Intrinsic Functions ..
179 INTRINSIC max
180* ..
181*
182* Test the input parameters.
183*
184 info = 0
185 IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
186 info = 1
187 ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
188 + .NOT.lsame(trans,'C')) THEN
189 info = 2
190 ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
191 info = 3
192 ELSE IF (n.LT.0) THEN
193 info = 4
194 ELSE IF (lda.LT.max(1,n)) THEN
195 info = 6
196 ELSE IF (incx.EQ.0) THEN
197 info = 8
198 END IF
199 IF (info.NE.0) THEN
200 CALL xerbla('DTRMV ',info)
201 RETURN
202 END IF
203*
204* Quick return if possible.
205*
206 IF (n.EQ.0) RETURN
207*
208 nounit = lsame(diag,'N')
209*
210* Set up the start point in X if the increment is not unity. This
211* will be ( N - 1 )*INCX too small for descending loops.
212*
213 IF (incx.LE.0) THEN
214 kx = 1 - (n-1)*incx
215 ELSE IF (incx.NE.1) THEN
216 kx = 1
217 END IF
218*
219* Start the operations. In this version the elements of A are
220* accessed sequentially with one pass through A.
221*
222 IF (lsame(trans,'N')) THEN
223*
224* Form x := A*x.
225*
226 IF (lsame(uplo,'U')) THEN
227 IF (incx.EQ.1) THEN
228 DO 20 j = 1,n
229 IF (x(j).NE.zero) THEN
230 temp = x(j)
231 DO 10 i = 1,j - 1
232 x(i) = x(i) + temp*a(i,j)
233 10 CONTINUE
234 IF (nounit) x(j) = x(j)*a(j,j)
235 END IF
236 20 CONTINUE
237 ELSE
238 jx = kx
239 DO 40 j = 1,n
240 IF (x(jx).NE.zero) THEN
241 temp = x(jx)
242 ix = kx
243 DO 30 i = 1,j - 1
244 x(ix) = x(ix) + temp*a(i,j)
245 ix = ix + incx
246 30 CONTINUE
247 IF (nounit) x(jx) = x(jx)*a(j,j)
248 END IF
249 jx = jx + incx
250 40 CONTINUE
251 END IF
252 ELSE
253 IF (incx.EQ.1) THEN
254 DO 60 j = n,1,-1
255 IF (x(j).NE.zero) THEN
256 temp = x(j)
257 DO 50 i = n,j + 1,-1
258 x(i) = x(i) + temp*a(i,j)
259 50 CONTINUE
260 IF (nounit) x(j) = x(j)*a(j,j)
261 END IF
262 60 CONTINUE
263 ELSE
264 kx = kx + (n-1)*incx
265 jx = kx
266 DO 80 j = n,1,-1
267 IF (x(jx).NE.zero) THEN
268 temp = x(jx)
269 ix = kx
270 DO 70 i = n,j + 1,-1
271 x(ix) = x(ix) + temp*a(i,j)
272 ix = ix - incx
273 70 CONTINUE
274 IF (nounit) x(jx) = x(jx)*a(j,j)
275 END IF
276 jx = jx - incx
277 80 CONTINUE
278 END IF
279 END IF
280 ELSE
281*
282* Form x := A**T*x.
283*
284 IF (lsame(uplo,'U')) THEN
285 IF (incx.EQ.1) THEN
286 DO 100 j = n,1,-1
287 temp = x(j)
288 IF (nounit) temp = temp*a(j,j)
289 DO 90 i = j - 1,1,-1
290 temp = temp + a(i,j)*x(i)
291 90 CONTINUE
292 x(j) = temp
293 100 CONTINUE
294 ELSE
295 jx = kx + (n-1)*incx
296 DO 120 j = n,1,-1
297 temp = x(jx)
298 ix = jx
299 IF (nounit) temp = temp*a(j,j)
300 DO 110 i = j - 1,1,-1
301 ix = ix - incx
302 temp = temp + a(i,j)*x(ix)
303 110 CONTINUE
304 x(jx) = temp
305 jx = jx - incx
306 120 CONTINUE
307 END IF
308 ELSE
309 IF (incx.EQ.1) THEN
310 DO 140 j = 1,n
311 temp = x(j)
312 IF (nounit) temp = temp*a(j,j)
313 DO 130 i = j + 1,n
314 temp = temp + a(i,j)*x(i)
315 130 CONTINUE
316 x(j) = temp
317 140 CONTINUE
318 ELSE
319 jx = kx
320 DO 160 j = 1,n
321 temp = x(jx)
322 ix = jx
323 IF (nounit) temp = temp*a(j,j)
324 DO 150 i = j + 1,n
325 ix = ix + incx
326 temp = temp + a(i,j)*x(ix)
327 150 CONTINUE
328 x(jx) = temp
329 jx = jx + incx
330 160 CONTINUE
331 END IF
332 END IF
333 END IF
334*
335 RETURN
336*
337* End of DTRMV
338*