Doxygen Source Code Documentation
Main Page Alphabetical List Data Structures File List Data Fields Globals Search
eis_reduc2.c
Go to the documentation of this file.00001
00002
00003
00004
00005
00006 #include "f2c.h"
00007
00008 int reduc2_(integer *nm, integer *n, doublereal *a,
00009 doublereal *b, doublereal *dl, integer *ierr)
00010 {
00011
00012 integer a_dim1, a_offset, b_dim1, b_offset, i__1, i__2, i__3;
00013
00014
00015 double sqrt(doublereal);
00016
00017
00018 static integer i__, j, k;
00019 static doublereal x, y;
00020 static integer i1, j1, nn;
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078 --dl;
00079 b_dim1 = *nm;
00080 b_offset = b_dim1 + 1;
00081 b -= b_offset;
00082 a_dim1 = *nm;
00083 a_offset = a_dim1 + 1;
00084 a -= a_offset;
00085
00086
00087 *ierr = 0;
00088 nn = abs(*n);
00089 if (*n < 0) {
00090 goto L100;
00091 }
00092
00093 i__1 = *n;
00094 for (i__ = 1; i__ <= i__1; ++i__) {
00095 i1 = i__ - 1;
00096
00097 i__2 = *n;
00098 for (j = i__; j <= i__2; ++j) {
00099 x = b[i__ + j * b_dim1];
00100 if (i__ == 1) {
00101 goto L40;
00102 }
00103
00104 i__3 = i1;
00105 for (k = 1; k <= i__3; ++k) {
00106
00107 x -= b[i__ + k * b_dim1] * b[j + k * b_dim1];
00108 }
00109
00110 L40:
00111 if (j != i__) {
00112 goto L60;
00113 }
00114 if (x <= 0.) {
00115 goto L1000;
00116 }
00117 y = sqrt(x);
00118 dl[i__] = y;
00119 goto L80;
00120 L60:
00121 b[j + i__ * b_dim1] = x / y;
00122 L80:
00123 ;
00124 }
00125 }
00126
00127
00128 L100:
00129 i__2 = nn;
00130 for (i__ = 1; i__ <= i__2; ++i__) {
00131 i1 = i__ + 1;
00132
00133 i__1 = i__;
00134 for (j = 1; j <= i__1; ++j) {
00135 x = a[j + i__ * a_dim1] * dl[j];
00136 if (j == i__) {
00137 goto L140;
00138 }
00139 j1 = j + 1;
00140
00141 i__3 = i__;
00142 for (k = j1; k <= i__3; ++k) {
00143
00144 x += a[k + i__ * a_dim1] * b[k + j * b_dim1];
00145 }
00146
00147 L140:
00148 if (i__ == nn) {
00149 goto L180;
00150 }
00151
00152 i__3 = nn;
00153 for (k = i1; k <= i__3; ++k) {
00154
00155 x += a[i__ + k * a_dim1] * b[k + j * b_dim1];
00156 }
00157
00158 L180:
00159 a[i__ + j * a_dim1] = x;
00160
00161 }
00162 }
00163
00164 i__1 = nn;
00165 for (i__ = 1; i__ <= i__1; ++i__) {
00166 i1 = i__ + 1;
00167 y = dl[i__];
00168
00169 i__2 = i__;
00170 for (j = 1; j <= i__2; ++j) {
00171 x = y * a[i__ + j * a_dim1];
00172 if (i__ == nn) {
00173 goto L280;
00174 }
00175
00176 i__3 = nn;
00177 for (k = i1; k <= i__3; ++k) {
00178
00179 x += a[k + j * a_dim1] * b[k + i__ * b_dim1];
00180 }
00181
00182 L280:
00183 a[i__ + j * a_dim1] = x;
00184
00185 }
00186 }
00187
00188 goto L1001;
00189
00190 L1000:
00191 *ierr = *n * 7 + 1;
00192 L1001:
00193 return 0;
00194 }
00195