ergo
template_lapack_stevr.h
Go to the documentation of this file.
1 /* Ergo, version 3.7, a program for linear scaling electronic structure
2  * calculations.
3  * Copyright (C) 2018 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek,
4  * and Anastasia Kruchinina.
5  *
6  * This program is free software: you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation, either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this program. If not, see <http://www.gnu.org/licenses/>.
18  *
19  * Primary academic reference:
20  * Ergo: An open-source program for linear-scaling electronic structure
21  * calculations,
22  * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia
23  * Kruchinina,
24  * SoftwareX 7, 107 (2018),
25  * <http://dx.doi.org/10.1016/j.softx.2018.03.005>
26  *
27  * For further information about Ergo, see <http://www.ergoscf.org>.
28  */
29 
30  /* This file belongs to the template_lapack part of the Ergo source
31  * code. The source files in the template_lapack directory are modified
32  * versions of files originally distributed as CLAPACK, see the
33  * Copyright/license notice in the file template_lapack/COPYING.
34  */
35 
36 
37 #ifndef TEMPLATE_LAPACK_STEVR_HEADER
38 #define TEMPLATE_LAPACK_STEVR_HEADER
39 
40 template<class Treal>
41 int template_lapack_stevr(const char *jobz, const char *range, const integer *n,
42  Treal * d__, Treal *e, const Treal *vl,
43  const Treal *vu, const integer *il,
44  const integer *iu, const Treal *abstol,
45  integer *m, Treal *w,
46  Treal *z__, const integer *ldz, integer *isuppz,
47  Treal *work,
48  integer *lwork, integer *iwork, integer *liwork,
49  integer *info)
50 {
51  /* System generated locals */
52  integer z_dim1, z_offset, i__1, i__2;
53  Treal d__1, d__2;
54 
55  /* Builtin functions */
56 
57  /* Local variables */
58  integer i__, j, jj;
59  Treal eps, vll, vuu, tmp1;
60  integer imax;
61  Treal rmin, rmax;
62  logical test;
63  Treal tnrm;
64  integer itmp1;
65  Treal sigma;
66  char order[1];
67  integer lwmin;
68  logical wantz;
69  logical alleig, indeig;
70  integer iscale, ieeeok, indibl, indifl;
71  logical valeig;
72  Treal safmin;
73  Treal bignum;
74  integer indisp;
75  integer indiwo;
76  integer liwmin;
77  logical tryrac;
78  integer nsplit;
79  Treal smlnum;
80  logical lquery;
81 
82 
83 /* -- LAPACK driver routine (version 3.2) -- */
84 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
85 /* November 2006 */
86 
87 /* .. Scalar Arguments .. */
88 /* .. */
89 /* .. Array Arguments .. */
90 /* .. */
91 
92 /* Purpose */
93 /* ======= */
94 
95 /* DSTEVR computes selected eigenvalues and, optionally, eigenvectors */
96 /* of a real symmetric tridiagonal matrix T. Eigenvalues and */
97 /* eigenvectors can be selected by specifying either a range of values */
98 /* or a range of indices for the desired eigenvalues. */
99 
100 /* Whenever possible, DSTEVR calls DSTEMR to compute the */
101 /* eigenspectrum using Relatively Robust Representations. DSTEMR */
102 /* computes eigenvalues by the dqds algorithm, while orthogonal */
103 /* eigenvectors are computed from various "good" L D L^T representations */
104 /* (also known as Relatively Robust Representations). Gram-Schmidt */
105 /* orthogonalization is avoided as far as possible. More specifically, */
106 /* the various steps of the algorithm are as follows. For the i-th */
107 /* unreduced block of T, */
108 /* (a) Compute T - sigma_i = L_i D_i L_i^T, such that L_i D_i L_i^T */
109 /* is a relatively robust representation, */
110 /* (b) Compute the eigenvalues, lambda_j, of L_i D_i L_i^T to high */
111 /* relative accuracy by the dqds algorithm, */
112 /* (c) If there is a cluster of close eigenvalues, "choose" sigma_i */
113 /* close to the cluster, and go to step (a), */
114 /* (d) Given the approximate eigenvalue lambda_j of L_i D_i L_i^T, */
115 /* compute the corresponding eigenvector by forming a */
116 /* rank-revealing twisted factorization. */
117 /* The desired accuracy of the output can be specified by the input */
118 /* parameter ABSTOL. */
119 
120 /* For more details, see "A new O(n^2) algorithm for the symmetric */
121 /* tridiagonal eigenvalue/eigenvector problem", by Inderjit Dhillon, */
122 /* Computer Science Division Technical Report No. UCB//CSD-97-971, */
123 /* UC Berkeley, May 1997. */
124 
125 
126 /* Note 1 : DSTEVR calls DSTEMR when the full spectrum is requested */
127 /* on machines which conform to the ieee-754 floating point standard. */
128 /* DSTEVR calls DSTEBZ and DSTEIN on non-ieee machines and */
129 /* when partial spectrum requests are made. */
130 
131 /* Normal execution of DSTEMR may create NaNs and infinities and */
132 /* hence may abort due to a floating point exception in environments */
133 /* which do not handle NaNs and infinities in the ieee standard default */
134 /* manner. */
135 
136 /* Arguments */
137 /* ========= */
138 
139 /* JOBZ (input) CHARACTER*1 */
140 /* = 'N': Compute eigenvalues only; */
141 /* = 'V': Compute eigenvalues and eigenvectors. */
142 
143 /* RANGE (input) CHARACTER*1 */
144 /* = 'A': all eigenvalues will be found. */
145 /* = 'V': all eigenvalues in the half-open interval (VL,VU] */
146 /* will be found. */
147 /* = 'I': the IL-th through IU-th eigenvalues will be found. */
148 /* ********* For RANGE = 'V' or 'I' and IU - IL < N - 1, DSTEBZ and */
149 /* ********* DSTEIN are called */
150 
151 /* N (input) INTEGER */
152 /* The order of the matrix. N >= 0. */
153 
154 /* D (input/output) DOUBLE PRECISION array, dimension (N) */
155 /* On entry, the n diagonal elements of the tridiagonal matrix */
156 /* A. */
157 /* On exit, D may be multiplied by a constant factor chosen */
158 /* to avoid over/underflow in computing the eigenvalues. */
159 
160 /* E (input/output) DOUBLE PRECISION array, dimension (max(1,N-1)) */
161 /* On entry, the (n-1) subdiagonal elements of the tridiagonal */
162 /* matrix A in elements 1 to N-1 of E. */
163 /* On exit, E may be multiplied by a constant factor chosen */
164 /* to avoid over/underflow in computing the eigenvalues. */
165 
166 /* VL (input) DOUBLE PRECISION */
167 /* VU (input) DOUBLE PRECISION */
168 /* If RANGE='V', the lower and upper bounds of the interval to */
169 /* be searched for eigenvalues. VL < VU. */
170 /* Not referenced if RANGE = 'A' or 'I'. */
171 
172 /* IL (input) INTEGER */
173 /* IU (input) INTEGER */
174 /* If RANGE='I', the indices (in ascending order) of the */
175 /* smallest and largest eigenvalues to be returned. */
176 /* 1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0. */
177 /* Not referenced if RANGE = 'A' or 'V'. */
178 
179 /* ABSTOL (input) DOUBLE PRECISION */
180 /* The absolute error tolerance for the eigenvalues. */
181 /* An approximate eigenvalue is accepted as converged */
182 /* when it is determined to lie in an interval [a,b] */
183 /* of width less than or equal to */
184 
185 /* ABSTOL + EPS * max( |a|,|b| ) , */
186 
187 /* where EPS is the machine precision. If ABSTOL is less than */
188 /* or equal to zero, then EPS*|T| will be used in its place, */
189 /* where |T| is the 1-norm of the tridiagonal matrix obtained */
190 /* by reducing A to tridiagonal form. */
191 
192 /* See "Computing Small Singular Values of Bidiagonal Matrices */
193 /* with Guaranteed High Relative Accuracy," by Demmel and */
194 /* Kahan, LAPACK Working Note #3. */
195 
196 /* If high relative accuracy is important, set ABSTOL to */
197 /* DLAMCH( 'Safe minimum' ). Doing so will guarantee that */
198 /* eigenvalues are computed to high relative accuracy when */
199 /* possible in future releases. The current code does not */
200 /* make any guarantees about high relative accuracy, but */
201 /* future releases will. See J. Barlow and J. Demmel, */
202 /* "Computing Accurate Eigensystems of Scaled Diagonally */
203 /* Dominant Matrices", LAPACK Working Note #7, for a discussion */
204 /* of which matrices define their eigenvalues to high relative */
205 /* accuracy. */
206 
207 /* M (output) INTEGER */
208 /* The total number of eigenvalues found. 0 <= M <= N. */
209 /* If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1. */
210 
211 /* W (output) DOUBLE PRECISION array, dimension (N) */
212 /* The first M elements contain the selected eigenvalues in */
213 /* ascending order. */
214 
215 /* Z (output) DOUBLE PRECISION array, dimension (LDZ, max(1,M) ) */
216 /* If JOBZ = 'V', then if INFO = 0, the first M columns of Z */
217 /* contain the orthonormal eigenvectors of the matrix A */
218 /* corresponding to the selected eigenvalues, with the i-th */
219 /* column of Z holding the eigenvector associated with W(i). */
220 /* Note: the user must ensure that at least max(1,M) columns are */
221 /* supplied in the array Z; if RANGE = 'V', the exact value of M */
222 /* is not known in advance and an upper bound must be used. */
223 
224 /* LDZ (input) INTEGER */
225 /* The leading dimension of the array Z. LDZ >= 1, and if */
226 /* JOBZ = 'V', LDZ >= max(1,N). */
227 
228 /* ISUPPZ (output) INTEGER array, dimension ( 2*max(1,M) ) */
229 /* The support of the eigenvectors in Z, i.e., the indices */
230 /* indicating the nonzero elements in Z. The i-th eigenvector */
231 /* is nonzero only in elements ISUPPZ( 2*i-1 ) through */
232 /* ISUPPZ( 2*i ). */
233 /* ********* Implemented only for RANGE = 'A' or 'I' and IU - IL = N - 1 */
234 
235 /* WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK)) */
236 /* On exit, if INFO = 0, WORK(1) returns the optimal (and */
237 /* minimal) LWORK. */
238 
239 /* LWORK (input) INTEGER */
240 /* The dimension of the array WORK. LWORK >= max(1,20*N). */
241 
242 /* If LWORK = -1, then a workspace query is assumed; the routine */
243 /* only calculates the optimal sizes of the WORK and IWORK */
244 /* arrays, returns these values as the first entries of the WORK */
245 /* and IWORK arrays, and no error message related to LWORK or */
246 /* LIWORK is issued by XERBLA. */
247 
248 /* IWORK (workspace/output) INTEGER array, dimension (MAX(1,LIWORK)) */
249 /* On exit, if INFO = 0, IWORK(1) returns the optimal (and */
250 /* minimal) LIWORK. */
251 
252 /* LIWORK (input) INTEGER */
253 /* The dimension of the array IWORK. LIWORK >= max(1,10*N). */
254 
255 /* If LIWORK = -1, then a workspace query is assumed; the */
256 /* routine only calculates the optimal sizes of the WORK and */
257 /* IWORK arrays, returns these values as the first entries of */
258 /* the WORK and IWORK arrays, and no error message related to */
259 /* LWORK or LIWORK is issued by XERBLA. */
260 
261 /* INFO (output) INTEGER */
262 /* = 0: successful exit */
263 /* < 0: if INFO = -i, the i-th argument had an illegal value */
264 /* > 0: Internal error */
265 
266 /* Further Details */
267 /* =============== */
268 
269 /* Based on contributions by */
270 /* Inderjit Dhillon, IBM Almaden, USA */
271 /* Osni Marques, LBNL/NERSC, USA */
272 /* Ken Stanley, Computer Science Division, University of */
273 /* California at Berkeley, USA */
274 
275 /* ===================================================================== */
276 
277 /* .. Parameters .. */
278 /* .. */
279 /* .. Local Scalars .. */
280 /* .. */
281 /* .. External Functions .. */
282 /* .. */
283 /* .. External Subroutines .. */
284 /* .. */
285 /* .. Intrinsic Functions .. */
286 /* .. */
287 /* .. Executable Statements .. */
288 
289 
290 /* Test the input parameters. */
291 
292  /* Parameter adjustments */
293  /* Table of constant values */
294 
295  integer c__10 = 10;
296  integer c__1 = 1;
297  integer c__2 = 2;
298  integer c__3 = 3;
299  integer c__4 = 4;
300 
301  --d__;
302  --e;
303  --w;
304  z_dim1 = *ldz;
305  z_offset = 1 + z_dim1;
306  z__ -= z_offset;
307  --isuppz;
308  --work;
309  --iwork;
310 
311  /* Function Body */
312  ieeeok = template_lapack_ilaenv(&c__10, "DSTEVR", "N", &c__1, &c__2, &c__3, &c__4, (ftnlen)6, (ftnlen)1);
313 
314  wantz = template_blas_lsame(jobz, "V");
315  alleig = template_blas_lsame(range, "A");
316  valeig = template_blas_lsame(range, "V");
317  indeig = template_blas_lsame(range, "I");
318 
319  lquery = *lwork == -1 || *liwork == -1;
320 /* Computing MAX */
321  i__1 = 1, i__2 = *n * 20;
322  lwmin = maxMACRO(i__1,i__2);
323 /* Computing MAX */
324  i__1 = 1, i__2 = *n * 10;
325  liwmin = maxMACRO(i__1,i__2);
326 
327 
328  *info = 0;
329  if (! (wantz || template_blas_lsame(jobz, "N"))) {
330  *info = -1;
331  } else if (! (alleig || valeig || indeig)) {
332  *info = -2;
333  } else if (*n < 0) {
334  *info = -3;
335  } else {
336  if (valeig) {
337  if (*n > 0 && *vu <= *vl) {
338  *info = -7;
339  }
340  } else if (indeig) {
341  if (*il < 1 || *il > maxMACRO(1,*n)) {
342  *info = -8;
343  } else if (*iu < minMACRO(*n,*il) || *iu > *n) {
344  *info = -9;
345  }
346  }
347  }
348  if (*info == 0) {
349  if (*ldz < 1 || ( wantz && *ldz < *n ) ) {
350  *info = -14;
351  }
352  }
353 
354  if (*info == 0) {
355  work[1] = (Treal) lwmin;
356  iwork[1] = liwmin;
357 
358  if (*lwork < lwmin && ! lquery) {
359  *info = -17;
360  } else if (*liwork < liwmin && ! lquery) {
361  *info = -19;
362  }
363  }
364 
365  if (*info != 0) {
366  i__1 = -(*info);
367  template_blas_erbla("STEVR", &i__1);
368  return 0;
369  } else if (lquery) {
370  return 0;
371  }
372 
373 /* Quick return if possible */
374 
375  *m = 0;
376  if (*n == 0) {
377  return 0;
378  }
379 
380  if (*n == 1) {
381  if (alleig || indeig) {
382  *m = 1;
383  w[1] = d__[1];
384  } else {
385  if (*vl < d__[1] && *vu >= d__[1]) {
386  *m = 1;
387  w[1] = d__[1];
388  }
389  }
390  if (wantz) {
391  z__[z_dim1 + 1] = 1.;
392  }
393  return 0;
394  }
395 
396 /* Get machine constants. */
397 
398  safmin = template_lapack_lamch("Safe minimum", (Treal)0);
399  eps = template_lapack_lamch("Precision", (Treal)0);
400  smlnum = safmin / eps;
401  bignum = 1. / smlnum;
402  rmin = template_blas_sqrt(smlnum);
403 /* Computing MIN */
404  d__1 = template_blas_sqrt(bignum), d__2 = 1. / template_blas_sqrt(template_blas_sqrt(safmin));
405  rmax = minMACRO(d__1,d__2);
406 
407 
408 /* Scale matrix to allowable range, if necessary. */
409 
410  iscale = 0;
411  vll = *vl;
412  vuu = *vu;
413 
414  tnrm = template_lapack_lanst("M", n, &d__[1], &e[1]);
415  if (tnrm > 0. && tnrm < rmin) {
416  iscale = 1;
417  sigma = rmin / tnrm;
418  } else if (tnrm > rmax) {
419  iscale = 1;
420  sigma = rmax / tnrm;
421  }
422  if (iscale == 1) {
423  template_blas_scal(n, &sigma, &d__[1], &c__1);
424  i__1 = *n - 1;
425  template_blas_scal(&i__1, &sigma, &e[1], &c__1);
426  if (valeig) {
427  vll = *vl * sigma;
428  vuu = *vu * sigma;
429  }
430  }
431 /* Initialize indices into workspaces. Note: These indices are used only */
432 /* if DSTERF or DSTEMR fail. */
433 /* IWORK(INDIBL:INDIBL+M-1) corresponds to IBLOCK in DSTEBZ and */
434 /* stores the block indices of each of the M<=N eigenvalues. */
435  indibl = 1;
436 /* IWORK(INDISP:INDISP+NSPLIT-1) corresponds to ISPLIT in DSTEBZ and */
437 /* stores the starting and finishing indices of each block. */
438  indisp = indibl + *n;
439 /* IWORK(INDIFL:INDIFL+N-1) stores the indices of eigenvectors */
440 /* that corresponding to eigenvectors that fail to converge in */
441 /* DSTEIN. This information is discarded; if any fail, the driver */
442 /* returns INFO > 0. */
443  indifl = indisp + *n;
444 /* INDIWO is the offset of the remaining integer workspace. */
445  indiwo = indisp + *n;
446 
447 /* If all eigenvalues are desired, then */
448 /* call DSTERF or DSTEMR. If this fails for some eigenvalue, then */
449 /* try DSTEBZ. */
450 
451 
452  test = FALSE_;
453  if (indeig) {
454  if (*il == 1 && *iu == *n) {
455  test = TRUE_;
456  }
457  }
458  if ((alleig || test) && ieeeok == 1) {
459  i__1 = *n - 1;
460  template_blas_copy(&i__1, &e[1], &c__1, &work[1], &c__1);
461  if (! wantz) {
462  template_blas_copy(n, &d__[1], &c__1, &w[1], &c__1);
463  template_lapack_sterf(n, &w[1], &work[1], info);
464  } else {
465  template_blas_copy(n, &d__[1], &c__1, &work[*n + 1], &c__1);
466  if (*abstol <= *n * 2. * eps) {
467  tryrac = TRUE_;
468  } else {
469  tryrac = FALSE_;
470  }
471  i__1 = *lwork - (*n << 1);
472  template_lapack_stemr(jobz, "A", n, &work[*n + 1], &work[1], vl, vu, il, iu, m,
473  &w[1], &z__[z_offset], ldz, n, &isuppz[1], &tryrac, &work[
474  (*n << 1) + 1], &i__1, &iwork[1], liwork, info);
475 
476  }
477  if (*info == 0) {
478  *m = *n;
479  goto L10;
480  }
481  *info = 0;
482  }
483 
484 /* Otherwise, call DSTEBZ and, if eigenvectors are desired, DSTEIN. */
485 
486  if (wantz) {
487  *(unsigned char *)order = 'B';
488  } else {
489  *(unsigned char *)order = 'E';
490  }
491  template_lapack_stebz(range, order, n, &vll, &vuu, il, iu, abstol, &d__[1], &e[1], m, &
492  nsplit, &w[1], &iwork[indibl], &iwork[indisp], &work[1], &iwork[
493  indiwo], info);
494 
495  if (wantz) {
496  template_lapack_stein(n, &d__[1], &e[1], m, &w[1], &iwork[indibl], &iwork[indisp], &
497  z__[z_offset], ldz, &work[1], &iwork[indiwo], &iwork[indifl],
498  info);
499  }
500 
501 /* If matrix was scaled, then rescale eigenvalues appropriately. */
502 
503 L10:
504  if (iscale == 1) {
505  if (*info == 0) {
506  imax = *m;
507  } else {
508  imax = *info - 1;
509  }
510  d__1 = 1. / sigma;
511  template_blas_scal(&imax, &d__1, &w[1], &c__1);
512  }
513 
514 /* If eigenvalues are not in order, then sort them, along with */
515 /* eigenvectors. */
516 
517  if (wantz) {
518  i__1 = *m - 1;
519  for (j = 1; j <= i__1; ++j) {
520  i__ = 0;
521  tmp1 = w[j];
522  i__2 = *m;
523  for (jj = j + 1; jj <= i__2; ++jj) {
524  if (w[jj] < tmp1) {
525  i__ = jj;
526  tmp1 = w[jj];
527  }
528 /* L20: */
529  }
530 
531  if (i__ != 0) {
532  itmp1 = iwork[i__];
533  w[i__] = w[j];
534  iwork[i__] = iwork[j];
535  w[j] = tmp1;
536  iwork[j] = itmp1;
537  template_blas_swap(n, &z__[i__ * z_dim1 + 1], &c__1, &z__[j * z_dim1 + 1],
538  &c__1);
539  }
540 /* L30: */
541  }
542  }
543 
544 /* Causes problems with tests 19 & 20: */
545 /* IF (wantz .and. INDEIG ) Z( 1,1) = Z(1,1) / 1.002 + .002 */
546 
547 
548  work[1] = (Treal) lwmin;
549  iwork[1] = liwmin;
550  return 0;
551 
552 /* End of DSTEVR */
553 
554 } /* dstevr_ */
555 
556 #endif
557 
int template_blas_scal(const integer *n, const Treal *da, Treal *dx, const integer *incx)
Definition: template_blas_scal.h:43
int integer
Definition: template_blas_common.h:40
integer template_lapack_ilaenv(const integer *ispec, const char *name__, const char *opts, const integer *n1, const integer *n2, const integer *n3, const integer *n4, ftnlen name_len, ftnlen opts_len)
Definition: template_lapack_common.cc:281
#define maxMACRO(a, b)
Definition: template_blas_common.h:45
#define minMACRO(a, b)
Definition: template_blas_common.h:46
int template_lapack_stebz(const char *range, const char *order, const integer *n, const Treal *vl, const Treal *vu, const integer *il, const integer *iu, const Treal *abstol, const Treal *d__, const Treal *e, integer *m, integer *nsplit, Treal *w, integer *iblock, integer *isplit, Treal *work, integer *iwork, integer *info)
Definition: template_lapack_stebz.h:42
int template_blas_erbla(const char *srname, integer *info)
Definition: template_blas_common.cc:146
int template_blas_swap(const integer *n, Treal *dx, const integer *incx, Treal *dy, const integer *incy)
Definition: template_blas_swap.h:42
int template_lapack_stein(const integer *n, const Treal *d__, const Treal *e, const integer *m, const Treal *w, const integer *iblock, const integer *isplit, Treal *z__, const integer *ldz, Treal *work, integer *iwork, integer *ifail, integer *info)
Definition: template_lapack_stein.h:42
Treal template_lapack_lamch(const char *cmach, Treal dummyReal)
Definition: template_lapack_lamch.h:202
bool logical
Definition: template_blas_common.h:41
#define TRUE_
Definition: template_lapack_common.h:42
int template_blas_copy(const integer *n, const Treal *dx, const integer *incx, Treal *dy, const integer *incy)
Definition: template_blas_copy.h:42
#define FALSE_
Definition: template_lapack_common.h:43
int ftnlen
Definition: template_blas_common.h:42
int template_lapack_sterf(const integer *n, Treal *d__, Treal *e, integer *info)
Definition: template_lapack_sterf.h:43
int template_lapack_stemr(const char *jobz, const char *range, const integer *n, Treal *d__, Treal *e, const Treal *vl, const Treal *vu, const integer *il, const integer *iu, integer *m, Treal *w, Treal *z__, const integer *ldz, const integer *nzc, integer *isuppz, logical *tryrac, Treal *work, integer *lwork, integer *iwork, integer *liwork, integer *info)
Definition: template_lapack_stemr.h:41
int template_lapack_stevr(const char *jobz, const char *range, const integer *n, Treal *d__, Treal *e, const Treal *vl, const Treal *vu, const integer *il, const integer *iu, const Treal *abstol, integer *m, Treal *w, Treal *z__, const integer *ldz, integer *isuppz, Treal *work, integer *lwork, integer *iwork, integer *liwork, integer *info)
Definition: template_lapack_stevr.h:41
Treal template_lapack_lanst(const char *norm, const integer *n, const Treal *d__, const Treal *e)
Definition: template_lapack_lanst.h:42
Treal template_blas_sqrt(Treal x)
logical template_blas_lsame(const char *ca, const char *cb)
Definition: template_blas_common.cc:46