LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
Collaboration diagram for double:

Functions

subroutine dgbmv (TRANS, M, N, KL, KU, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
 DGBMV More...
 
subroutine dgemv (TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
 DGEMV More...
 
subroutine dger (M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
 DGER More...
 
subroutine dsbmv (UPLO, N, K, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
 DSBMV More...
 
subroutine dspmv (UPLO, N, ALPHA, AP, X, INCX, BETA, Y, INCY)
 DSPMV More...
 
subroutine dspr (UPLO, N, ALPHA, X, INCX, AP)
 DSPR More...
 
subroutine dspr2 (UPLO, N, ALPHA, X, INCX, Y, INCY, AP)
 DSPR2 More...
 
subroutine dsymv (UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
 DSYMV More...
 
subroutine dsyr (UPLO, N, ALPHA, X, INCX, A, LDA)
 DSYR More...
 
subroutine dsyr2 (UPLO, N, ALPHA, X, INCX, Y, INCY, A, LDA)
 DSYR2 More...
 
subroutine dtbmv (UPLO, TRANS, DIAG, N, K, A, LDA, X, INCX)
 DTBMV More...
 
subroutine dtbsv (UPLO, TRANS, DIAG, N, K, A, LDA, X, INCX)
 DTBSV More...
 
subroutine dtpmv (UPLO, TRANS, DIAG, N, AP, X, INCX)
 DTPMV More...
 
subroutine dtpsv (UPLO, TRANS, DIAG, N, AP, X, INCX)
 DTPSV More...
 
subroutine dtrmv (UPLO, TRANS, DIAG, N, A, LDA, X, INCX)
 DTRMV More...
 

Detailed Description

This is the group of double LEVEL 2 BLAS routines.

Function Documentation

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 of 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 of 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 of 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.
Date
November 2015
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 187 of file dgbmv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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 of 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.
Date
November 2015
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 158 of file dgemv.f.

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

Here is the call graph for this function:

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 of 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 of 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 of 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.
Date
November 2011
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 132 of file dger.f.

132 *
133 * -- Reference BLAS level2 routine (version 3.4.0) --
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 * November 2011
137 *
138 * .. Scalar Arguments ..
139  DOUBLE PRECISION alpha
140  INTEGER incx,incy,lda,m,n
141 * ..
142 * .. Array Arguments ..
143  DOUBLE PRECISION a(lda,*),x(*),y(*)
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,jy,kx
155 * ..
156 * .. External Subroutines ..
157  EXTERNAL xerbla
158 * ..
159 * .. Intrinsic Functions ..
160  INTRINSIC max
161 * ..
162 *
163 * Test the input parameters.
164 *
165  info = 0
166  IF (m.LT.0) THEN
167  info = 1
168  ELSE IF (n.LT.0) THEN
169  info = 2
170  ELSE IF (incx.EQ.0) THEN
171  info = 5
172  ELSE IF (incy.EQ.0) THEN
173  info = 7
174  ELSE IF (lda.LT.max(1,m)) THEN
175  info = 9
176  END IF
177  IF (info.NE.0) THEN
178  CALL xerbla('DGER ',info)
179  RETURN
180  END IF
181 *
182 * Quick return if possible.
183 *
184  IF ((m.EQ.0) .OR. (n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
185 *
186 * Start the operations. In this version the elements of A are
187 * accessed sequentially with one pass through A.
188 *
189  IF (incy.GT.0) THEN
190  jy = 1
191  ELSE
192  jy = 1 - (n-1)*incy
193  END IF
194  IF (incx.EQ.1) THEN
195  DO 20 j = 1,n
196  IF (y(jy).NE.zero) THEN
197  temp = alpha*y(jy)
198  DO 10 i = 1,m
199  a(i,j) = a(i,j) + x(i)*temp
200  10 CONTINUE
201  END IF
202  jy = jy + incy
203  20 CONTINUE
204  ELSE
205  IF (incx.GT.0) THEN
206  kx = 1
207  ELSE
208  kx = 1 - (m-1)*incx
209  END IF
210  DO 40 j = 1,n
211  IF (y(jy).NE.zero) THEN
212  temp = alpha*y(jy)
213  ix = kx
214  DO 30 i = 1,m
215  a(i,j) = a(i,j) + x(ix)*temp
216  ix = ix + incx
217  30 CONTINUE
218  END IF
219  jy = jy + incy
220  40 CONTINUE
221  END IF
222 *
223  RETURN
224 *
225 * End of DGER .
226 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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 of 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.
Date
November 2011
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 186 of file dsbmv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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 of 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.
Date
November 2011
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 149 of file dspmv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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.
Date
November 2011
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 dspr.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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 of 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.
Date
November 2011
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 144 of file dspr2.f.

144 *
145 * -- Reference BLAS level2 routine (version 3.4.0) --
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 * November 2011
149 *
150 * .. Scalar Arguments ..
151  DOUBLE PRECISION alpha
152  INTEGER incx,incy,n
153  CHARACTER uplo
154 * ..
155 * .. Array Arguments ..
156  DOUBLE PRECISION ap(*),x(*),y(*)
157 * ..
158 *
159 * =====================================================================
160 *
161 * .. Parameters ..
162  DOUBLE PRECISION zero
163  parameter(zero=0.0d+0)
164 * ..
165 * .. Local Scalars ..
166  DOUBLE PRECISION temp1,temp2
167  INTEGER i,info,ix,iy,j,jx,jy,k,kk,kx,ky
168 * ..
169 * .. External Functions ..
170  LOGICAL lsame
171  EXTERNAL lsame
172 * ..
173 * .. External Subroutines ..
174  EXTERNAL xerbla
175 * ..
176 *
177 * Test the input parameters.
178 *
179  info = 0
180  IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
181  info = 1
182  ELSE IF (n.LT.0) THEN
183  info = 2
184  ELSE IF (incx.EQ.0) THEN
185  info = 5
186  ELSE IF (incy.EQ.0) THEN
187  info = 7
188  END IF
189  IF (info.NE.0) THEN
190  CALL xerbla('DSPR2 ',info)
191  RETURN
192  END IF
193 *
194 * Quick return if possible.
195 *
196  IF ((n.EQ.0) .OR. (alpha.EQ.zero)) RETURN
197 *
198 * Set up the start points in X and Y if the increments are not both
199 * unity.
200 *
201  IF ((incx.NE.1) .OR. (incy.NE.1)) THEN
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  jx = kx
213  jy = ky
214  END IF
215 *
216 * Start the operations. In this version the elements of the array AP
217 * are accessed sequentially with one pass through AP.
218 *
219  kk = 1
220  IF (lsame(uplo,'U')) THEN
221 *
222 * Form A when upper triangle is stored in AP.
223 *
224  IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
225  DO 20 j = 1,n
226  IF ((x(j).NE.zero) .OR. (y(j).NE.zero)) THEN
227  temp1 = alpha*y(j)
228  temp2 = alpha*x(j)
229  k = kk
230  DO 10 i = 1,j
231  ap(k) = ap(k) + x(i)*temp1 + y(i)*temp2
232  k = k + 1
233  10 CONTINUE
234  END IF
235  kk = kk + j
236  20 CONTINUE
237  ELSE
238  DO 40 j = 1,n
239  IF ((x(jx).NE.zero) .OR. (y(jy).NE.zero)) THEN
240  temp1 = alpha*y(jy)
241  temp2 = alpha*x(jx)
242  ix = kx
243  iy = ky
244  DO 30 k = kk,kk + j - 1
245  ap(k) = ap(k) + x(ix)*temp1 + y(iy)*temp2
246  ix = ix + incx
247  iy = iy + incy
248  30 CONTINUE
249  END IF
250  jx = jx + incx
251  jy = jy + incy
252  kk = kk + j
253  40 CONTINUE
254  END IF
255  ELSE
256 *
257 * Form A when lower triangle is stored in AP.
258 *
259  IF ((incx.EQ.1) .AND. (incy.EQ.1)) THEN
260  DO 60 j = 1,n
261  IF ((x(j).NE.zero) .OR. (y(j).NE.zero)) THEN
262  temp1 = alpha*y(j)
263  temp2 = alpha*x(j)
264  k = kk
265  DO 50 i = j,n
266  ap(k) = ap(k) + x(i)*temp1 + y(i)*temp2
267  k = k + 1
268  50 CONTINUE
269  END IF
270  kk = kk + n - j + 1
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 k = kk,kk + n - j
280  ap(k) = ap(k) + 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  kk = kk + n - j + 1
288  80 CONTINUE
289  END IF
290  END IF
291 *
292  RETURN
293 *
294 * End of DSPR2 .
295 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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 of 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.
Date
November 2011
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 154 of file dsymv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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.
Date
November 2011
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 134 of file dsyr.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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 of 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.
Date
November 2011
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 149 of file dsyr2.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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
           tranformed 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.
Date
November 2011
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 188 of file dtbmv.f.

188 *
189 * -- Reference BLAS level2 routine (version 3.4.0) --
190 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
191 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
192 * November 2011
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('DTBMV ',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 sequentially with one pass through A.
265 *
266  IF (lsame(trans,'N')) THEN
267 *
268 * Form x := A*x.
269 *
270  IF (lsame(uplo,'U')) THEN
271  kplus1 = k + 1
272  IF (incx.EQ.1) THEN
273  DO 20 j = 1,n
274  IF (x(j).NE.zero) THEN
275  temp = x(j)
276  l = kplus1 - j
277  DO 10 i = max(1,j-k),j - 1
278  x(i) = x(i) + temp*a(l+i,j)
279  10 CONTINUE
280  IF (nounit) x(j) = x(j)*a(kplus1,j)
281  END IF
282  20 CONTINUE
283  ELSE
284  jx = kx
285  DO 40 j = 1,n
286  IF (x(jx).NE.zero) THEN
287  temp = x(jx)
288  ix = kx
289  l = kplus1 - j
290  DO 30 i = max(1,j-k),j - 1
291  x(ix) = x(ix) + temp*a(l+i,j)
292  ix = ix + incx
293  30 CONTINUE
294  IF (nounit) x(jx) = x(jx)*a(kplus1,j)
295  END IF
296  jx = jx + incx
297  IF (j.GT.k) kx = kx + incx
298  40 CONTINUE
299  END IF
300  ELSE
301  IF (incx.EQ.1) THEN
302  DO 60 j = n,1,-1
303  IF (x(j).NE.zero) THEN
304  temp = x(j)
305  l = 1 - j
306  DO 50 i = min(n,j+k),j + 1,-1
307  x(i) = x(i) + temp*a(l+i,j)
308  50 CONTINUE
309  IF (nounit) x(j) = x(j)*a(1,j)
310  END IF
311  60 CONTINUE
312  ELSE
313  kx = kx + (n-1)*incx
314  jx = kx
315  DO 80 j = n,1,-1
316  IF (x(jx).NE.zero) THEN
317  temp = x(jx)
318  ix = kx
319  l = 1 - j
320  DO 70 i = min(n,j+k),j + 1,-1
321  x(ix) = x(ix) + temp*a(l+i,j)
322  ix = ix - incx
323  70 CONTINUE
324  IF (nounit) x(jx) = x(jx)*a(1,j)
325  END IF
326  jx = jx - incx
327  IF ((n-j).GE.k) kx = kx - incx
328  80 CONTINUE
329  END IF
330  END IF
331  ELSE
332 *
333 * Form x := A**T*x.
334 *
335  IF (lsame(uplo,'U')) THEN
336  kplus1 = k + 1
337  IF (incx.EQ.1) THEN
338  DO 100 j = n,1,-1
339  temp = x(j)
340  l = kplus1 - j
341  IF (nounit) temp = temp*a(kplus1,j)
342  DO 90 i = j - 1,max(1,j-k),-1
343  temp = temp + a(l+i,j)*x(i)
344  90 CONTINUE
345  x(j) = temp
346  100 CONTINUE
347  ELSE
348  kx = kx + (n-1)*incx
349  jx = kx
350  DO 120 j = n,1,-1
351  temp = x(jx)
352  kx = kx - incx
353  ix = kx
354  l = kplus1 - j
355  IF (nounit) temp = temp*a(kplus1,j)
356  DO 110 i = j - 1,max(1,j-k),-1
357  temp = temp + a(l+i,j)*x(ix)
358  ix = ix - incx
359  110 CONTINUE
360  x(jx) = temp
361  jx = jx - incx
362  120 CONTINUE
363  END IF
364  ELSE
365  IF (incx.EQ.1) THEN
366  DO 140 j = 1,n
367  temp = x(j)
368  l = 1 - j
369  IF (nounit) temp = temp*a(1,j)
370  DO 130 i = j + 1,min(n,j+k)
371  temp = temp + a(l+i,j)*x(i)
372  130 CONTINUE
373  x(j) = temp
374  140 CONTINUE
375  ELSE
376  jx = kx
377  DO 160 j = 1,n
378  temp = x(jx)
379  kx = kx + incx
380  ix = kx
381  l = 1 - j
382  IF (nounit) temp = temp*a(1,j)
383  DO 150 i = j + 1,min(n,j+k)
384  temp = temp + a(l+i,j)*x(ix)
385  ix = ix + incx
386  150 CONTINUE
387  x(jx) = temp
388  jx = jx + incx
389  160 CONTINUE
390  END IF
391  END IF
392  END IF
393 *
394  RETURN
395 *
396 * End of DTBMV .
397 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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.
Date
November 2011
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 191 of file dtbsv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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
           tranformed 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.
Date
November 2011
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 144 of file dtpmv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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.
Date
November 2011
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 dtpsv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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 of 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 of 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
           tranformed 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.
Date
November 2011
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 149 of file dtrmv.f.

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

Here is the call graph for this function:

Here is the caller graph for this function: