1 /* ./src_f77/dstev.f -- translated by f2c (version 20030320).
2    You must link the resulting object file with the libraries:
3 	-lf2c -lm   (in that order)
4 */
5 
6 #include <punc/vf2c.h>
7 
8 /* Table of constant values */
9 
10 static integer c__1 = 1;
11 
dstev_(char * jobz,integer * n,doublereal * d__,doublereal * e,doublereal * z__,integer * ldz,doublereal * work,integer * info,ftnlen jobz_len)12 /* Subroutine */ int dstev_(char *jobz, integer *n, doublereal *d__,
13 	doublereal *e, doublereal *z__, integer *ldz, doublereal *work,
14 	integer *info, ftnlen jobz_len)
15 {
16     /* System generated locals */
17     integer z_dim1, z_offset, i__1;
18     doublereal d__1;
19 
20     /* Builtin functions */
21     double sqrt(doublereal);
22 
23     /* Local variables */
24     static doublereal eps;
25     static integer imax;
26     static doublereal rmin, rmax, tnrm;
27     extern /* Subroutine */ int dscal_(integer *, doublereal *, doublereal *,
28 	    integer *);
29     static doublereal sigma;
30     extern logical lsame_(char *, char *, ftnlen, ftnlen);
31     static logical wantz;
32     extern doublereal dlamch_(char *, ftnlen);
33     static integer iscale;
34     static doublereal safmin;
35     extern /* Subroutine */ int xerbla_(char *, integer *, ftnlen);
36     static doublereal bignum;
37     extern doublereal dlanst_(char *, integer *, doublereal *, doublereal *,
38 	    ftnlen);
39     extern /* Subroutine */ int dsterf_(integer *, doublereal *, doublereal *,
40 	     integer *), dsteqr_(char *, integer *, doublereal *, doublereal *
41 	    , doublereal *, integer *, doublereal *, integer *, ftnlen);
42     static doublereal smlnum;
43 
44 
45 /*  -- LAPACK driver routine (version 3.0) -- */
46 /*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., */
47 /*     Courant Institute, Argonne National Lab, and Rice University */
48 /*     September 30, 1994 */
49 
50 /*     .. Scalar Arguments .. */
51 /*     .. */
52 /*     .. Array Arguments .. */
53 /*     .. */
54 
55 /*  Purpose */
56 /*  ======= */
57 
58 /*  DSTEV computes all eigenvalues and, optionally, eigenvectors of a */
59 /*  real symmetric tridiagonal matrix A. */
60 
61 /*  Arguments */
62 /*  ========= */
63 
64 /*  JOBZ    (input) CHARACTER*1 */
65 /*          = 'N':  Compute eigenvalues only; */
66 /*          = 'V':  Compute eigenvalues and eigenvectors. */
67 
68 /*  N       (input) INTEGER */
69 /*          The order of the matrix.  N >= 0. */
70 
71 /*  D       (input/output) DOUBLE PRECISION array, dimension (N) */
72 /*          On entry, the n diagonal elements of the tridiagonal matrix */
73 /*          A. */
74 /*          On exit, if INFO = 0, the eigenvalues in ascending order. */
75 
76 /*  E       (input/output) DOUBLE PRECISION array, dimension (N) */
77 /*          On entry, the (n-1) subdiagonal elements of the tridiagonal */
78 /*          matrix A, stored in elements 1 to N-1 of E; E(N) need not */
79 /*          be set, but is used by the routine. */
80 /*          On exit, the contents of E are destroyed. */
81 
82 /*  Z       (output) DOUBLE PRECISION array, dimension (LDZ, N) */
83 /*          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal */
84 /*          eigenvectors of the matrix A, with the i-th column of Z */
85 /*          holding the eigenvector associated with D(i). */
86 /*          If JOBZ = 'N', then Z is not referenced. */
87 
88 /*  LDZ     (input) INTEGER */
89 /*          The leading dimension of the array Z.  LDZ >= 1, and if */
90 /*          JOBZ = 'V', LDZ >= max(1,N). */
91 
92 /*  WORK    (workspace) DOUBLE PRECISION array, dimension (max(1,2*N-2)) */
93 /*          If JOBZ = 'N', WORK is not referenced. */
94 
95 /*  INFO    (output) INTEGER */
96 /*          = 0:  successful exit */
97 /*          < 0:  if INFO = -i, the i-th argument had an illegal value */
98 /*          > 0:  if INFO = i, the algorithm failed to converge; i */
99 /*                off-diagonal elements of E did not converge to zero. */
100 
101 /*  ===================================================================== */
102 
103 /*     .. Parameters .. */
104 /*     .. */
105 /*     .. Local Scalars .. */
106 /*     .. */
107 /*     .. External Functions .. */
108 /*     .. */
109 /*     .. External Subroutines .. */
110 /*     .. */
111 /*     .. Intrinsic Functions .. */
112 /*     .. */
113 /*     .. Executable Statements .. */
114 
115 /*     Test the input parameters. */
116 
117     /* Parameter adjustments */
118     --d__;
119     --e;
120     z_dim1 = *ldz;
121     z_offset = 1 + z_dim1;
122     z__ -= z_offset;
123     --work;
124 
125     /* Function Body */
126     wantz = lsame_(jobz, "V", (ftnlen)1, (ftnlen)1);
127 
128     *info = 0;
129     if (! (wantz || lsame_(jobz, "N", (ftnlen)1, (ftnlen)1))) {
130 	*info = -1;
131     } else if (*n < 0) {
132 	*info = -2;
133     } else if (*ldz < 1 || wantz && *ldz < *n) {
134 	*info = -6;
135     }
136 
137     if (*info != 0) {
138 	i__1 = -(*info);
139 	xerbla_("DSTEV ", &i__1, (ftnlen)6);
140 	return 0;
141     }
142 
143 /*     Quick return if possible */
144 
145     if (*n == 0) {
146 	return 0;
147     }
148 
149     if (*n == 1) {
150 	if (wantz) {
151 	    z__[z_dim1 + 1] = 1.;
152 	}
153 	return 0;
154     }
155 
156 /*     Get machine constants. */
157 
158     safmin = dlamch_("Safe minimum", (ftnlen)12);
159     eps = dlamch_("Precision", (ftnlen)9);
160     smlnum = safmin / eps;
161     bignum = 1. / smlnum;
162     rmin = sqrt(smlnum);
163     rmax = sqrt(bignum);
164 
165 /*     Scale matrix to allowable range, if necessary. */
166 
167     iscale = 0;
168     tnrm = dlanst_("M", n, &d__[1], &e[1], (ftnlen)1);
169     if (tnrm > 0. && tnrm < rmin) {
170 	iscale = 1;
171 	sigma = rmin / tnrm;
172     } else if (tnrm > rmax) {
173 	iscale = 1;
174 	sigma = rmax / tnrm;
175     }
176     if (iscale == 1) {
177 	dscal_(n, &sigma, &d__[1], &c__1);
178 	i__1 = *n - 1;
179 	dscal_(&i__1, &sigma, &e[1], &c__1);
180     }
181 
182 /*     For eigenvalues only, call DSTERF.  For eigenvalues and */
183 /*     eigenvectors, call DSTEQR. */
184 
185     if (! wantz) {
186 	dsterf_(n, &d__[1], &e[1], info);
187     } else {
188 	dsteqr_("I", n, &d__[1], &e[1], &z__[z_offset], ldz, &work[1], info, (
189 		ftnlen)1);
190     }
191 
192 /*     If matrix was scaled, then rescale eigenvalues appropriately. */
193 
194     if (iscale == 1) {
195 	if (*info == 0) {
196 	    imax = *n;
197 	} else {
198 	    imax = *info - 1;
199 	}
200 	d__1 = 1. / sigma;
201 	dscal_(&imax, &d__1, &d__[1], &c__1);
202     }
203 
204     return 0;
205 
206 /*     End of DSTEV */
207 
208 } /* dstev_ */
209 
210