00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015
00016
00017
00018 static integer c__3 = 3;
00019 static integer c__1 = 1;
00020 static integer c__5 = 5;
00021 static integer c__20 = 20;
00022
00023 int dchkgl_(integer *nin, integer *nout)
00024 {
00025
00026 static char fmt_9999[] = "(1x,\002.. test output of DGGBAL .. \002)";
00027 static char fmt_9998[] = "(1x,\002value of largest test error "
00028 " = \002,d12.3)";
00029 static char fmt_9997[] = "(1x,\002example number where info is not zero "
00030 " = \002,i4)";
00031 static char fmt_9996[] = "(1x,\002example number where ILO or IHI wrong "
00032 " = \002,i4)";
00033 static char fmt_9995[] = "(1x,\002example number having largest error "
00034 " = \002,i4)";
00035 static char fmt_9994[] = "(1x,\002number of examples where info is not 0"
00036 " = \002,i4)";
00037 static char fmt_9993[] = "(1x,\002total number of examples tested "
00038 " = \002,i4)";
00039
00040
00041 integer i__1, i__2;
00042 doublereal d__1, d__2, d__3;
00043
00044
00045 integer s_rsle(cilist *), do_lio(integer *, integer *, char *, ftnlen),
00046 e_rsle(void), s_wsfe(cilist *), e_wsfe(void), do_fio(integer *,
00047 char *, ftnlen);
00048
00049
00050 doublereal a[400] , b[400] ;
00051 integer i__, j, n;
00052 doublereal ain[400] , bin[400] ;
00053 integer ihi, ilo;
00054 doublereal eps;
00055 integer knt, info, lmax[5];
00056 doublereal rmax, vmax, work[120];
00057 integer ihiin, ninfo, iloin;
00058 doublereal anorm, bnorm;
00059 extern int dggbal_(char *, integer *, doublereal *,
00060 integer *, doublereal *, integer *, integer *, integer *,
00061 doublereal *, doublereal *, doublereal *, integer *);
00062 extern doublereal dlamch_(char *), dlange_(char *, integer *,
00063 integer *, doublereal *, integer *, doublereal *);
00064 doublereal lscale[20], rscale[20], lsclin[20], rsclin[20];
00065
00066
00067 static cilist io___6 = { 0, 0, 0, 0, 0 };
00068 static cilist io___9 = { 0, 0, 0, 0, 0 };
00069 static cilist io___12 = { 0, 0, 0, 0, 0 };
00070 static cilist io___14 = { 0, 0, 0, 0, 0 };
00071 static cilist io___17 = { 0, 0, 0, 0, 0 };
00072 static cilist io___19 = { 0, 0, 0, 0, 0 };
00073 static cilist io___21 = { 0, 0, 0, 0, 0 };
00074 static cilist io___23 = { 0, 0, 0, 0, 0 };
00075 static cilist io___34 = { 0, 0, 0, fmt_9999, 0 };
00076 static cilist io___35 = { 0, 0, 0, fmt_9998, 0 };
00077 static cilist io___36 = { 0, 0, 0, fmt_9997, 0 };
00078 static cilist io___37 = { 0, 0, 0, fmt_9996, 0 };
00079 static cilist io___38 = { 0, 0, 0, fmt_9995, 0 };
00080 static cilist io___39 = { 0, 0, 0, fmt_9994, 0 };
00081 static cilist io___40 = { 0, 0, 0, fmt_9993, 0 };
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122 lmax[0] = 0;
00123 lmax[1] = 0;
00124 lmax[2] = 0;
00125 ninfo = 0;
00126 knt = 0;
00127 rmax = 0.;
00128
00129 eps = dlamch_("Precision");
00130
00131 L10:
00132
00133 io___6.ciunit = *nin;
00134 s_rsle(&io___6);
00135 do_lio(&c__3, &c__1, (char *)&n, (ftnlen)sizeof(integer));
00136 e_rsle();
00137 if (n == 0) {
00138 goto L90;
00139 }
00140 i__1 = n;
00141 for (i__ = 1; i__ <= i__1; ++i__) {
00142 io___9.ciunit = *nin;
00143 s_rsle(&io___9);
00144 i__2 = n;
00145 for (j = 1; j <= i__2; ++j) {
00146 do_lio(&c__5, &c__1, (char *)&a[i__ + j * 20 - 21], (ftnlen)
00147 sizeof(doublereal));
00148 }
00149 e_rsle();
00150
00151 }
00152
00153 i__1 = n;
00154 for (i__ = 1; i__ <= i__1; ++i__) {
00155 io___12.ciunit = *nin;
00156 s_rsle(&io___12);
00157 i__2 = n;
00158 for (j = 1; j <= i__2; ++j) {
00159 do_lio(&c__5, &c__1, (char *)&b[i__ + j * 20 - 21], (ftnlen)
00160 sizeof(doublereal));
00161 }
00162 e_rsle();
00163
00164 }
00165
00166 io___14.ciunit = *nin;
00167 s_rsle(&io___14);
00168 do_lio(&c__3, &c__1, (char *)&iloin, (ftnlen)sizeof(integer));
00169 do_lio(&c__3, &c__1, (char *)&ihiin, (ftnlen)sizeof(integer));
00170 e_rsle();
00171 i__1 = n;
00172 for (i__ = 1; i__ <= i__1; ++i__) {
00173 io___17.ciunit = *nin;
00174 s_rsle(&io___17);
00175 i__2 = n;
00176 for (j = 1; j <= i__2; ++j) {
00177 do_lio(&c__5, &c__1, (char *)&ain[i__ + j * 20 - 21], (ftnlen)
00178 sizeof(doublereal));
00179 }
00180 e_rsle();
00181
00182 }
00183 i__1 = n;
00184 for (i__ = 1; i__ <= i__1; ++i__) {
00185 io___19.ciunit = *nin;
00186 s_rsle(&io___19);
00187 i__2 = n;
00188 for (j = 1; j <= i__2; ++j) {
00189 do_lio(&c__5, &c__1, (char *)&bin[i__ + j * 20 - 21], (ftnlen)
00190 sizeof(doublereal));
00191 }
00192 e_rsle();
00193
00194 }
00195
00196 io___21.ciunit = *nin;
00197 s_rsle(&io___21);
00198 i__1 = n;
00199 for (i__ = 1; i__ <= i__1; ++i__) {
00200 do_lio(&c__5, &c__1, (char *)&lsclin[i__ - 1], (ftnlen)sizeof(
00201 doublereal));
00202 }
00203 e_rsle();
00204 io___23.ciunit = *nin;
00205 s_rsle(&io___23);
00206 i__1 = n;
00207 for (i__ = 1; i__ <= i__1; ++i__) {
00208 do_lio(&c__5, &c__1, (char *)&rsclin[i__ - 1], (ftnlen)sizeof(
00209 doublereal));
00210 }
00211 e_rsle();
00212
00213 anorm = dlange_("M", &n, &n, a, &c__20, work);
00214 bnorm = dlange_("M", &n, &n, b, &c__20, work);
00215
00216 ++knt;
00217
00218 dggbal_("B", &n, a, &c__20, b, &c__20, &ilo, &ihi, lscale, rscale, work, &
00219 info);
00220
00221 if (info != 0) {
00222 ++ninfo;
00223 lmax[0] = knt;
00224 }
00225
00226 if (ilo != iloin || ihi != ihiin) {
00227 ++ninfo;
00228 lmax[1] = knt;
00229 }
00230
00231 vmax = 0.;
00232 i__1 = n;
00233 for (i__ = 1; i__ <= i__1; ++i__) {
00234 i__2 = n;
00235 for (j = 1; j <= i__2; ++j) {
00236
00237 d__2 = vmax, d__3 = (d__1 = a[i__ + j * 20 - 21] - ain[i__ + j *
00238 20 - 21], abs(d__1));
00239 vmax = max(d__2,d__3);
00240
00241 d__2 = vmax, d__3 = (d__1 = b[i__ + j * 20 - 21] - bin[i__ + j *
00242 20 - 21], abs(d__1));
00243 vmax = max(d__2,d__3);
00244
00245 }
00246
00247 }
00248
00249 i__1 = n;
00250 for (i__ = 1; i__ <= i__1; ++i__) {
00251
00252 d__2 = vmax, d__3 = (d__1 = lscale[i__ - 1] - lsclin[i__ - 1], abs(
00253 d__1));
00254 vmax = max(d__2,d__3);
00255
00256 d__2 = vmax, d__3 = (d__1 = rscale[i__ - 1] - rsclin[i__ - 1], abs(
00257 d__1));
00258 vmax = max(d__2,d__3);
00259
00260 }
00261
00262 vmax /= eps * max(anorm,bnorm);
00263
00264 if (vmax > rmax) {
00265 lmax[2] = knt;
00266 rmax = vmax;
00267 }
00268
00269 goto L10;
00270
00271 L90:
00272
00273 io___34.ciunit = *nout;
00274 s_wsfe(&io___34);
00275 e_wsfe();
00276
00277 io___35.ciunit = *nout;
00278 s_wsfe(&io___35);
00279 do_fio(&c__1, (char *)&rmax, (ftnlen)sizeof(doublereal));
00280 e_wsfe();
00281 io___36.ciunit = *nout;
00282 s_wsfe(&io___36);
00283 do_fio(&c__1, (char *)&lmax[0], (ftnlen)sizeof(integer));
00284 e_wsfe();
00285 io___37.ciunit = *nout;
00286 s_wsfe(&io___37);
00287 do_fio(&c__1, (char *)&lmax[1], (ftnlen)sizeof(integer));
00288 e_wsfe();
00289 io___38.ciunit = *nout;
00290 s_wsfe(&io___38);
00291 do_fio(&c__1, (char *)&lmax[2], (ftnlen)sizeof(integer));
00292 e_wsfe();
00293 io___39.ciunit = *nout;
00294 s_wsfe(&io___39);
00295 do_fio(&c__1, (char *)&ninfo, (ftnlen)sizeof(integer));
00296 e_wsfe();
00297 io___40.ciunit = *nout;
00298 s_wsfe(&io___40);
00299 do_fio(&c__1, (char *)&knt, (ftnlen)sizeof(integer));
00300 e_wsfe();
00301
00302 return 0;
00303
00304
00305
00306 }