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

◆ dsbev_2stage()

subroutine dsbev_2stage ( character  jobz,
character  uplo,
integer  n,
integer  kd,
double precision, dimension( ldab, * )  ab,
integer  ldab,
double precision, dimension( * )  w,
double precision, dimension( ldz, * )  z,
integer  ldz,
double precision, dimension( * )  work,
integer  lwork,
integer  info 
)

DSBEV_2STAGE computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices

Download DSBEV_2STAGE + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DSBEV_2STAGE computes all the eigenvalues and, optionally, eigenvectors of
 a real symmetric band matrix A using the 2stage technique for
 the reduction to tridiagonal.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
                  Not available in this release.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]KD
          KD is INTEGER
          The number of superdiagonals of the matrix A if UPLO = 'U',
          or the number of subdiagonals if UPLO = 'L'.  KD >= 0.
[in,out]AB
          AB is DOUBLE PRECISION array, dimension (LDAB, N)
          On entry, the upper or lower triangle of the symmetric band
          matrix A, stored in the first KD+1 rows of the array.  The
          j-th column of A is stored in the j-th column of the array AB
          as follows:
          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).

          On exit, AB is overwritten by values generated during the
          reduction to tridiagonal form.  If UPLO = 'U', the first
          superdiagonal and the diagonal of the tridiagonal matrix T
          are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
          the diagonal and first subdiagonal of T are returned in the
          first two rows of AB.
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KD + 1.
[out]W
          W is DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]Z
          Z is DOUBLE PRECISION array, dimension (LDZ, N)
          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
          eigenvectors of the matrix A, with the i-th column of Z
          holding the eigenvector associated with W(i).
          If JOBZ = 'N', then Z is not referenced.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= max(1,N).
[out]WORK
          WORK is DOUBLE PRECISION array, dimension LWORK
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of the array WORK. LWORK >= 1, when N <= 1;
          otherwise
          If JOBZ = 'N' and N > 1, LWORK must be queried.
                                   LWORK = MAX(1, dimension) where
                                   dimension = (2KD+1)*N + KD*NTHREADS + N
                                   where KD is the size of the band.
                                   NTHREADS is the number of threads used when
                                   openMP compilation is enabled, otherwise =1.
          If JOBZ = 'V' and N > 1, LWORK must be queried. Not yet available.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, the algorithm failed to converge; i
                off-diagonal elements of an intermediate tridiagonal
                form did not converge to zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  All details about the 2stage techniques are available in:

  Azzam Haidar, Hatem Ltaief, and Jack Dongarra.
  Parallel reduction to condensed forms for symmetric eigenvalue problems
  using aggregated fine-grained and memory-aware kernels. In Proceedings
  of 2011 International Conference for High Performance Computing,
  Networking, Storage and Analysis (SC '11), New York, NY, USA,
  Article 8 , 11 pages.
  http://doi.acm.org/10.1145/2063384.2063394

  A. Haidar, J. Kurzak, P. Luszczek, 2013.
  An improved parallel singular value algorithm and its implementation
  for multicore hardware, In Proceedings of 2013 International Conference
  for High Performance Computing, Networking, Storage and Analysis (SC '13).
  Denver, Colorado, USA, 2013.
  Article 90, 12 pages.
  http://doi.acm.org/10.1145/2503210.2503292

  A. Haidar, R. Solca, S. Tomov, T. Schulthess and J. Dongarra.
  A novel hybrid CPU-GPU generalized eigensolver for electronic structure
  calculations based on fine-grained memory aware tasks.
  International Journal of High Performance Computing Applications.
  Volume 28 Issue 2, Pages 196-209, May 2014.
  http://hpc.sagepub.com/content/28/2/196

Definition at line 202 of file dsbev_2stage.f.

204*
205 IMPLICIT NONE
206*
207* -- LAPACK driver routine --
208* -- LAPACK is a software package provided by Univ. of Tennessee, --
209* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
210*
211* .. Scalar Arguments ..
212 CHARACTER JOBZ, UPLO
213 INTEGER INFO, KD, LDAB, LDZ, N, LWORK
214* ..
215* .. Array Arguments ..
216 DOUBLE PRECISION AB( LDAB, * ), W( * ), WORK( * ), Z( LDZ, * )
217* ..
218*
219* =====================================================================
220*
221* .. Parameters ..
222 DOUBLE PRECISION ZERO, ONE
223 parameter( zero = 0.0d0, one = 1.0d0 )
224* ..
225* .. Local Scalars ..
226 LOGICAL LOWER, WANTZ, LQUERY
227 INTEGER IINFO, IMAX, INDE, INDWRK, ISCALE,
228 $ LLWORK, LWMIN, LHTRD, LWTRD, IB, INDHOUS
229 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
230 $ SMLNUM
231* ..
232* .. External Functions ..
233 LOGICAL LSAME
234 INTEGER ILAENV2STAGE
235 DOUBLE PRECISION DLAMCH, DLANSB
236 EXTERNAL lsame, dlamch, dlansb, ilaenv2stage
237* ..
238* .. External Subroutines ..
239 EXTERNAL dlascl, dscal, dsteqr, dsterf, xerbla,
241* ..
242* .. Intrinsic Functions ..
243 INTRINSIC sqrt
244* ..
245* .. Executable Statements ..
246*
247* Test the input parameters.
248*
249 wantz = lsame( jobz, 'V' )
250 lower = lsame( uplo, 'L' )
251 lquery = ( lwork.EQ.-1 )
252*
253 info = 0
254 IF( .NOT.( lsame( jobz, 'N' ) ) ) THEN
255 info = -1
256 ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
257 info = -2
258 ELSE IF( n.LT.0 ) THEN
259 info = -3
260 ELSE IF( kd.LT.0 ) THEN
261 info = -4
262 ELSE IF( ldab.LT.kd+1 ) THEN
263 info = -6
264 ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
265 info = -9
266 END IF
267*
268 IF( info.EQ.0 ) THEN
269 IF( n.LE.1 ) THEN
270 lwmin = 1
271 work( 1 ) = lwmin
272 ELSE
273 ib = ilaenv2stage( 2, 'DSYTRD_SB2ST', jobz,
274 $ n, kd, -1, -1 )
275 lhtrd = ilaenv2stage( 3, 'DSYTRD_SB2ST', jobz,
276 $ n, kd, ib, -1 )
277 lwtrd = ilaenv2stage( 4, 'DSYTRD_SB2ST', jobz,
278 $ n, kd, ib, -1 )
279 lwmin = n + lhtrd + lwtrd
280 work( 1 ) = lwmin
281 ENDIF
282*
283 IF( lwork.LT.lwmin .AND. .NOT.lquery )
284 $ info = -11
285 END IF
286*
287 IF( info.NE.0 ) THEN
288 CALL xerbla( 'DSBEV_2STAGE ', -info )
289 RETURN
290 ELSE IF( lquery ) THEN
291 RETURN
292 END IF
293*
294* Quick return if possible
295*
296 IF( n.EQ.0 )
297 $ RETURN
298*
299 IF( n.EQ.1 ) THEN
300 IF( lower ) THEN
301 w( 1 ) = ab( 1, 1 )
302 ELSE
303 w( 1 ) = ab( kd+1, 1 )
304 END IF
305 IF( wantz )
306 $ z( 1, 1 ) = one
307 RETURN
308 END IF
309*
310* Get machine constants.
311*
312 safmin = dlamch( 'Safe minimum' )
313 eps = dlamch( 'Precision' )
314 smlnum = safmin / eps
315 bignum = one / smlnum
316 rmin = sqrt( smlnum )
317 rmax = sqrt( bignum )
318*
319* Scale matrix to allowable range, if necessary.
320*
321 anrm = dlansb( 'M', uplo, n, kd, ab, ldab, work )
322 iscale = 0
323 IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
324 iscale = 1
325 sigma = rmin / anrm
326 ELSE IF( anrm.GT.rmax ) THEN
327 iscale = 1
328 sigma = rmax / anrm
329 END IF
330 IF( iscale.EQ.1 ) THEN
331 IF( lower ) THEN
332 CALL dlascl( 'B', kd, kd, one, sigma, n, n, ab, ldab, info )
333 ELSE
334 CALL dlascl( 'Q', kd, kd, one, sigma, n, n, ab, ldab, info )
335 END IF
336 END IF
337*
338* Call DSYTRD_SB2ST to reduce symmetric band matrix to tridiagonal form.
339*
340 inde = 1
341 indhous = inde + n
342 indwrk = indhous + lhtrd
343 llwork = lwork - indwrk + 1
344*
345 CALL dsytrd_sb2st( "N", jobz, uplo, n, kd, ab, ldab, w,
346 $ work( inde ), work( indhous ), lhtrd,
347 $ work( indwrk ), llwork, iinfo )
348*
349* For eigenvalues only, call DSTERF. For eigenvectors, call SSTEQR.
350*
351 IF( .NOT.wantz ) THEN
352 CALL dsterf( n, w, work( inde ), info )
353 ELSE
354 CALL dsteqr( jobz, n, w, work( inde ), z, ldz, work( indwrk ),
355 $ info )
356 END IF
357*
358* If matrix was scaled, then rescale eigenvalues appropriately.
359*
360 IF( iscale.EQ.1 ) THEN
361 IF( info.EQ.0 ) THEN
362 imax = n
363 ELSE
364 imax = info - 1
365 END IF
366 CALL dscal( imax, one / sigma, w, 1 )
367 END IF
368*
369* Set WORK(1) to optimal workspace size.
370*
371 work( 1 ) = lwmin
372*
373 RETURN
374*
375* End of DSBEV_2STAGE
376*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dsytrd_sb2st(stage1, vect, uplo, n, kd, ab, ldab, d, e, hous, lhous, work, lwork, info)
DSYTRD_SB2ST reduces a real symmetric band matrix A to real symmetric tridiagonal form T
integer function ilaenv2stage(ispec, name, opts, n1, n2, n3, n4)
ILAENV2STAGE
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
double precision function dlansb(norm, uplo, n, k, ab, ldab, work)
DLANSB returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition dlansb.f:129
subroutine dlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition dlascl.f:143
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79
subroutine dsteqr(compz, n, d, e, z, ldz, work, info)
DSTEQR
Definition dsteqr.f:131
subroutine dsterf(n, d, e, info)
DSTERF
Definition dsterf.f:86
Here is the call graph for this function:
Here is the caller graph for this function: