Doxygen Source Code Documentation
eis_imtqlv.c File Reference
#include "f2c.h"
Go to the source code of this file.
Functions | |
int | imtqlv_ (integer *n, doublereal *d__, doublereal *e, doublereal *e2, doublereal *w, integer *ind, integer *ierr, doublereal *rv1) |
Variables | |
doublereal | c_b11 = 1. |
Function Documentation
|
Definition at line 12 of file eis_imtqlv.c. References abs, c_b11, d_sign(), ind, l, p, and pythag_(). Referenced by rsm_().
00015 { 00016 /* System generated locals */ 00017 integer i__1, i__2; 00018 doublereal d__1, d__2; 00019 00020 /* Builtin functions */ 00021 double d_sign(doublereal *, doublereal *); 00022 00023 /* Local variables */ 00024 static doublereal b, c__, f, g; 00025 static integer i__, j, k, l, m; 00026 static doublereal p, r__, s; 00027 static integer ii; 00028 extern doublereal pythag_(doublereal *, doublereal *); 00029 static integer tag, mml; 00030 static doublereal tst1, tst2; 00031 00032 00033 00034 /* THIS SUBROUTINE IS A VARIANT OF IMTQL1 WHICH IS A TRANSLATION OF 00035 */ 00036 /* ALGOL PROCEDURE IMTQL1, NUM. MATH. 12, 377-383(1968) BY MARTIN AND 00037 */ 00038 /* WILKINSON, AS MODIFIED IN NUM. MATH. 15, 450(1970) BY DUBRULLE. */ 00039 /* HANDBOOK FOR AUTO. COMP., VOL.II-LINEAR ALGEBRA, 241-248(1971). */ 00040 00041 /* THIS SUBROUTINE FINDS THE EIGENVALUES OF A SYMMETRIC TRIDIAGONAL */ 00042 /* MATRIX BY THE IMPLICIT QL METHOD AND ASSOCIATES WITH THEM */ 00043 /* THEIR CORRESPONDING SUBMATRIX INDICES. */ 00044 00045 /* ON INPUT */ 00046 00047 /* N IS THE ORDER OF THE MATRIX. */ 00048 00049 /* D CONTAINS THE DIAGONAL ELEMENTS OF THE INPUT MATRIX. */ 00050 00051 /* E CONTAINS THE SUBDIAGONAL ELEMENTS OF THE INPUT MATRIX */ 00052 /* IN ITS LAST N-1 POSITIONS. E(1) IS ARBITRARY. */ 00053 00054 /* E2 CONTAINS THE SQUARES OF THE CORRESPONDING ELEMENTS OF E. */ 00055 /* E2(1) IS ARBITRARY. */ 00056 00057 /* ON OUTPUT */ 00058 00059 /* D AND E ARE UNALTERED. */ 00060 00061 /* ELEMENTS OF E2, CORRESPONDING TO ELEMENTS OF E REGARDED */ 00062 /* AS NEGLIGIBLE, HAVE BEEN REPLACED BY ZERO CAUSING THE */ 00063 /* MATRIX TO SPLIT INTO A DIRECT SUM OF SUBMATRICES. */ 00064 /* E2(1) IS ALSO SET TO ZERO. */ 00065 00066 /* W CONTAINS THE EIGENVALUES IN ASCENDING ORDER. IF AN */ 00067 /* ERROR EXIT IS MADE, THE EIGENVALUES ARE CORRECT AND */ 00068 /* ORDERED FOR INDICES 1,2,...IERR-1, BUT MAY NOT BE */ 00069 /* THE SMALLEST EIGENVALUES. */ 00070 00071 /* IND CONTAINS THE SUBMATRIX INDICES ASSOCIATED WITH THE */ 00072 /* CORRESPONDING EIGENVALUES IN W -- 1 FOR EIGENVALUES */ 00073 /* BELONGING TO THE FIRST SUBMATRIX FROM THE TOP, */ 00074 /* 2 FOR THOSE BELONGING TO THE SECOND SUBMATRIX, ETC.. */ 00075 00076 /* IERR IS SET TO */ 00077 /* ZERO FOR NORMAL RETURN, */ 00078 /* J IF THE J-TH EIGENVALUE HAS NOT BEEN */ 00079 /* DETERMINED AFTER 30 ITERATIONS. */ 00080 00081 /* RV1 IS A TEMPORARY STORAGE ARRAY. */ 00082 00083 /* CALLS PYTHAG FOR DSQRT(A*A + B*B) . */ 00084 00085 /* QUESTIONS AND COMMENTS SHOULD BE DIRECTED TO BURTON S. GARBOW, */ 00086 /* MATHEMATICS AND COMPUTER SCIENCE DIV, ARGONNE NATIONAL LABORATORY 00087 */ 00088 00089 /* THIS VERSION DATED AUGUST 1983. */ 00090 00091 /* ------------------------------------------------------------------ 00092 */ 00093 00094 /* Parameter adjustments */ 00095 --rv1; 00096 --ind; 00097 --w; 00098 --e2; 00099 --e; 00100 --d__; 00101 00102 /* Function Body */ 00103 *ierr = 0; 00104 k = 0; 00105 tag = 0; 00106 00107 i__1 = *n; 00108 for (i__ = 1; i__ <= i__1; ++i__) { 00109 w[i__] = d__[i__]; 00110 if (i__ != 1) { 00111 rv1[i__ - 1] = e[i__]; 00112 } 00113 /* L100: */ 00114 } 00115 00116 e2[1] = 0.; 00117 rv1[*n] = 0.; 00118 00119 i__1 = *n; 00120 for (l = 1; l <= i__1; ++l) { 00121 j = 0; 00122 /* .......... LOOK FOR SMALL SUB-DIAGONAL ELEMENT .......... */ 00123 L105: 00124 i__2 = *n; 00125 for (m = l; m <= i__2; ++m) { 00126 if (m == *n) { 00127 goto L120; 00128 } 00129 tst1 = (d__1 = w[m], abs(d__1)) + (d__2 = w[m + 1], abs(d__2)); 00130 tst2 = tst1 + (d__1 = rv1[m], abs(d__1)); 00131 if (tst2 == tst1) { 00132 goto L120; 00133 } 00134 /* .......... GUARD AGAINST UNDERFLOWED ELEMENT OF E2 ........ 00135 .. */ 00136 if (e2[m + 1] == 0.) { 00137 goto L125; 00138 } 00139 /* L110: */ 00140 } 00141 00142 L120: 00143 if (m <= k) { 00144 goto L130; 00145 } 00146 if (m != *n) { 00147 e2[m + 1] = 0.; 00148 } 00149 L125: 00150 k = m; 00151 ++tag; 00152 L130: 00153 p = w[l]; 00154 if (m == l) { 00155 goto L215; 00156 } 00157 if (j == 30) { 00158 goto L1000; 00159 } 00160 ++j; 00161 /* .......... FORM SHIFT .......... */ 00162 g = (w[l + 1] - p) / (rv1[l] * 2.); 00163 r__ = pythag_(&g, &c_b11); 00164 g = w[m] - p + rv1[l] / (g + d_sign(&r__, &g)); 00165 s = 1.; 00166 c__ = 1.; 00167 p = 0.; 00168 mml = m - l; 00169 /* .......... FOR I=M-1 STEP -1 UNTIL L DO -- .......... */ 00170 i__2 = mml; 00171 for (ii = 1; ii <= i__2; ++ii) { 00172 i__ = m - ii; 00173 f = s * rv1[i__]; 00174 b = c__ * rv1[i__]; 00175 r__ = pythag_(&f, &g); 00176 rv1[i__ + 1] = r__; 00177 if (r__ == 0.) { 00178 goto L210; 00179 } 00180 s = f / r__; 00181 c__ = g / r__; 00182 g = w[i__ + 1] - p; 00183 r__ = (w[i__] - g) * s + c__ * 2. * b; 00184 p = s * r__; 00185 w[i__ + 1] = g + p; 00186 g = c__ * r__ - b; 00187 /* L200: */ 00188 } 00189 00190 w[l] -= p; 00191 rv1[l] = g; 00192 rv1[m] = 0.; 00193 goto L105; 00194 /* .......... RECOVER FROM UNDERFLOW .......... */ 00195 L210: 00196 w[i__ + 1] -= p; 00197 rv1[m] = 0.; 00198 goto L105; 00199 /* .......... ORDER EIGENVALUES .......... */ 00200 L215: 00201 if (l == 1) { 00202 goto L250; 00203 } 00204 /* .......... FOR I=L STEP -1 UNTIL 2 DO -- .......... */ 00205 i__2 = l; 00206 for (ii = 2; ii <= i__2; ++ii) { 00207 i__ = l + 2 - ii; 00208 if (p >= w[i__ - 1]) { 00209 goto L270; 00210 } 00211 w[i__] = w[i__ - 1]; 00212 ind[i__] = ind[i__ - 1]; 00213 /* L230: */ 00214 } 00215 00216 L250: 00217 i__ = 1; 00218 L270: 00219 w[i__] = p; 00220 ind[i__] = tag; 00221 /* L290: */ 00222 } 00223 00224 goto L1001; 00225 /* .......... SET ERROR -- NO CONVERGENCE TO AN */ 00226 /* EIGENVALUE AFTER 30 ITERATIONS .......... */ 00227 L1000: 00228 *ierr = l; 00229 L1001: 00230 return 0; 00231 } /* imtqlv_ */ |
Variable Documentation
|
Definition at line 10 of file eis_imtqlv.c. Referenced by imtqlv_(). |