Doxygen Source Code Documentation
eis_cortb.c File Reference
#include "f2c.h"
Go to the source code of this file.
Functions | |
int | cortb_ (integer *nm, integer *low, integer *igh, doublereal *ar, doublereal *ai, doublereal *ortr, doublereal *orti, integer *m, doublereal *zr, doublereal *zi) |
Function Documentation
|
Definition at line 8 of file eis_cortb.c. References mp.
00011 { 00012 /* System generated locals */ 00013 integer ar_dim1, ar_offset, ai_dim1, ai_offset, zr_dim1, zr_offset, 00014 zi_dim1, zi_offset, i__1, i__2, i__3; 00015 00016 /* Local variables */ 00017 static doublereal h__; 00018 static integer i__, j, la; 00019 static doublereal gi, gr; 00020 static integer mm, mp, kp1, mp1; 00021 00022 00023 00024 /* THIS SUBROUTINE IS A TRANSLATION OF A COMPLEX ANALOGUE OF */ 00025 /* THE ALGOL PROCEDURE ORTBAK, NUM. MATH. 12, 349-368(1968) */ 00026 /* BY MARTIN AND WILKINSON. */ 00027 /* HANDBOOK FOR AUTO. COMP., VOL.II-LINEAR ALGEBRA, 339-358(1971). */ 00028 00029 /* THIS SUBROUTINE FORMS THE EIGENVECTORS OF A COMPLEX GENERAL */ 00030 /* MATRIX BY BACK TRANSFORMING THOSE OF THE CORRESPONDING */ 00031 /* UPPER HESSENBERG MATRIX DETERMINED BY CORTH. */ 00032 00033 /* ON INPUT */ 00034 00035 /* NM MUST BE SET TO THE ROW DIMENSION OF TWO-DIMENSIONAL */ 00036 /* ARRAY PARAMETERS AS DECLARED IN THE CALLING PROGRAM */ 00037 /* DIMENSION STATEMENT. */ 00038 00039 /* LOW AND IGH ARE INTEGERS DETERMINED BY THE BALANCING */ 00040 /* SUBROUTINE CBAL. IF CBAL HAS NOT BEEN USED, */ 00041 /* SET LOW=1 AND IGH EQUAL TO THE ORDER OF THE MATRIX. */ 00042 00043 /* AR AND AI CONTAIN INFORMATION ABOUT THE UNITARY */ 00044 /* TRANSFORMATIONS USED IN THE REDUCTION BY CORTH */ 00045 /* IN THEIR STRICT LOWER TRIANGLES. */ 00046 00047 /* ORTR AND ORTI CONTAIN FURTHER INFORMATION ABOUT THE */ 00048 /* TRANSFORMATIONS USED IN THE REDUCTION BY CORTH. */ 00049 /* ONLY ELEMENTS LOW THROUGH IGH ARE USED. */ 00050 00051 /* M IS THE NUMBER OF COLUMNS OF ZR AND ZI TO BE BACK TRANSFORMED. 00052 */ 00053 00054 /* ZR AND ZI CONTAIN THE REAL AND IMAGINARY PARTS, */ 00055 /* RESPECTIVELY, OF THE EIGENVECTORS TO BE */ 00056 /* BACK TRANSFORMED IN THEIR FIRST M COLUMNS. */ 00057 00058 /* ON OUTPUT */ 00059 00060 /* ZR AND ZI CONTAIN THE REAL AND IMAGINARY PARTS, */ 00061 /* RESPECTIVELY, OF THE TRANSFORMED EIGENVECTORS */ 00062 /* IN THEIR FIRST M COLUMNS. */ 00063 00064 /* ORTR AND ORTI HAVE BEEN ALTERED. */ 00065 00066 /* NOTE THAT CORTB PRESERVES VECTOR EUCLIDEAN NORMS. */ 00067 00068 /* QUESTIONS AND COMMENTS SHOULD BE DIRECTED TO BURTON S. GARBOW, */ 00069 /* MATHEMATICS AND COMPUTER SCIENCE DIV, ARGONNE NATIONAL LABORATORY 00070 */ 00071 00072 /* THIS VERSION DATED AUGUST 1983. */ 00073 00074 /* ------------------------------------------------------------------ 00075 */ 00076 00077 /* Parameter adjustments */ 00078 --orti; 00079 --ortr; 00080 ai_dim1 = *nm; 00081 ai_offset = ai_dim1 + 1; 00082 ai -= ai_offset; 00083 ar_dim1 = *nm; 00084 ar_offset = ar_dim1 + 1; 00085 ar -= ar_offset; 00086 zi_dim1 = *nm; 00087 zi_offset = zi_dim1 + 1; 00088 zi -= zi_offset; 00089 zr_dim1 = *nm; 00090 zr_offset = zr_dim1 + 1; 00091 zr -= zr_offset; 00092 00093 /* Function Body */ 00094 if (*m == 0) { 00095 goto L200; 00096 } 00097 la = *igh - 1; 00098 kp1 = *low + 1; 00099 if (la < kp1) { 00100 goto L200; 00101 } 00102 /* .......... FOR MP=IGH-1 STEP -1 UNTIL LOW+1 DO -- .......... */ 00103 i__1 = la; 00104 for (mm = kp1; mm <= i__1; ++mm) { 00105 mp = *low + *igh - mm; 00106 if (ar[mp + (mp - 1) * ar_dim1] == 0. && ai[mp + (mp - 1) * ai_dim1] 00107 == 0.) { 00108 goto L140; 00109 } 00110 /* .......... H BELOW IS NEGATIVE OF H FORMED IN CORTH .......... 00111 */ 00112 h__ = ar[mp + (mp - 1) * ar_dim1] * ortr[mp] + ai[mp + (mp - 1) * 00113 ai_dim1] * orti[mp]; 00114 mp1 = mp + 1; 00115 00116 i__2 = *igh; 00117 for (i__ = mp1; i__ <= i__2; ++i__) { 00118 ortr[i__] = ar[i__ + (mp - 1) * ar_dim1]; 00119 orti[i__] = ai[i__ + (mp - 1) * ai_dim1]; 00120 /* L100: */ 00121 } 00122 00123 i__2 = *m; 00124 for (j = 1; j <= i__2; ++j) { 00125 gr = 0.; 00126 gi = 0.; 00127 00128 i__3 = *igh; 00129 for (i__ = mp; i__ <= i__3; ++i__) { 00130 gr = gr + ortr[i__] * zr[i__ + j * zr_dim1] + orti[i__] * zi[ 00131 i__ + j * zi_dim1]; 00132 gi = gi + ortr[i__] * zi[i__ + j * zi_dim1] - orti[i__] * zr[ 00133 i__ + j * zr_dim1]; 00134 /* L110: */ 00135 } 00136 00137 gr /= h__; 00138 gi /= h__; 00139 00140 i__3 = *igh; 00141 for (i__ = mp; i__ <= i__3; ++i__) { 00142 zr[i__ + j * zr_dim1] = zr[i__ + j * zr_dim1] + gr * ortr[i__] 00143 - gi * orti[i__]; 00144 zi[i__ + j * zi_dim1] = zi[i__ + j * zi_dim1] + gr * orti[i__] 00145 + gi * ortr[i__]; 00146 /* L120: */ 00147 } 00148 00149 /* L130: */ 00150 } 00151 00152 L140: 00153 ; 00154 } 00155 00156 L200: 00157 return 0; 00158 } /* cortb_ */ |