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 complex c_b1 = {0.f,0.f};
00019 static complex c_b2 = {1.f,0.f};
00020 static integer c__4 = 4;
00021 static integer c__1 = 1;
00022 static integer c__0 = 0;
00023
00024 int cdrvbd_(integer *nsizes, integer *mm, integer *nn,
00025 integer *ntypes, logical *dotype, integer *iseed, real *thresh,
00026 complex *a, integer *lda, complex *u, integer *ldu, complex *vt,
00027 integer *ldvt, complex *asav, complex *usav, complex *vtsav, real *s,
00028 real *ssav, real *e, complex *work, integer *lwork, real *rwork,
00029 integer *iwork, integer *nounit, integer *info)
00030 {
00031
00032
00033 static char cjob[1*4] = "N" "O" "S" "A";
00034
00035
00036 static char fmt_9996[] = "(\002 CDRVBD: \002,a,\002 returned INFO=\002,i"
00037 "6,\002.\002,/9x,\002M=\002,i6,\002, N=\002,i6,\002, JTYPE=\002,i"
00038 "6,\002, ISEED=(\002,3(i5,\002,\002),i5,\002)\002)";
00039 static char fmt_9995[] = "(\002 CDRVBD: \002,a,\002 returned INFO=\002,i"
00040 "6,\002.\002,/9x,\002M=\002,i6,\002, N=\002,i6,\002, JTYPE=\002,i"
00041 "6,\002, LSWORK=\002,i6,/9x,\002ISEED=(\002,3(i5,\002,\002),i5"
00042 ",\002)\002)";
00043 static char fmt_9999[] = "(\002 SVD -- Complex Singular Value Decomposit"
00044 "ion Driver \002,/\002 Matrix types (see CDRVBD for details):\002"
00045 ",//\002 1 = Zero matrix\002,/\002 2 = Identity matrix\002,/\002 "
00046 "3 = Evenly spaced singular values near 1\002,/\002 4 = Evenly sp"
00047 "aced singular values near underflow\002,/\002 5 = Evenly spaced "
00048 "singular values near overflow\002,//\002 Tests performed: ( A is"
00049 " dense, U and V are unitary,\002,/19x,\002 S is an array, and Up"
00050 "artial, VTpartial, and\002,/19x,\002 Spartial are partially comp"
00051 "uted U, VT and S),\002,/)";
00052 static char fmt_9998[] = "(\002 Tests performed with Test Threshold ="
00053 " \002,f8.2,/\002 CGESVD: \002,/\002 1 = | A - U diag(S) VT | / ("
00054 " |A| max(M,N) ulp ) \002,/\002 2 = | I - U**T U | / ( M ulp )"
00055 " \002,/\002 3 = | I - VT VT**T | / ( N ulp ) \002,/\002 4 = 0 if"
00056 " S contains min(M,N) nonnegative values in\002,\002 decreasing o"
00057 "rder, else 1/ulp\002,/\002 5 = | U - Upartial | / ( M ulp )\002,/"
00058 "\002 6 = | VT - VTpartial | / ( N ulp )\002,/\002 7 = | S - Spar"
00059 "tial | / ( min(M,N) ulp |S| )\002,/\002 CGESDD: \002,/\002 8 = |"
00060 " A - U diag(S) VT | / ( |A| max(M,N) ulp ) \002,/\002 9 = | I - "
00061 "U**T U | / ( M ulp ) \002,/\00210 = | I - VT VT**T | / ( N ulp ) "
00062 "\002,/\00211 = 0 if S contains min(M,N) nonnegative values in"
00063 "\002,\002 decreasing order, else 1/ulp\002,/\00212 = | U - Upart"
00064 "ial | / ( M ulp )\002,/\00213 = | VT - VTpartial | / ( N ulp "
00065 ")\002,/\00214 = | S - Spartial | / ( min(M,N) ulp |S| )\002,//)";
00066 static char fmt_9997[] = "(\002 M=\002,i5,\002, N=\002,i5,\002, type "
00067 "\002,i1,\002, IWS=\002,i1,\002, seed=\002,4(i4,\002,\002),\002 t"
00068 "est(\002,i1,\002)=\002,g11.4)";
00069
00070
00071 integer a_dim1, a_offset, asav_dim1, asav_offset, u_dim1, u_offset,
00072 usav_dim1, usav_offset, vt_dim1, vt_offset, vtsav_dim1,
00073 vtsav_offset, i__1, i__2, i__3, i__4, i__5, i__6, i__7, i__8,
00074 i__9, i__10, i__11, i__12, i__13, i__14;
00075 real r__1, r__2, r__3;
00076
00077
00078 integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00079
00080
00081 integer i__, j, m, n;
00082 real dif, div;
00083 integer ijq, iju;
00084 real ulp;
00085 char jobq[1], jobu[1];
00086 integer mmax, nmax;
00087 real unfl, ovfl;
00088 integer ijvt;
00089 extern int cbdt01_(integer *, integer *, integer *,
00090 complex *, integer *, complex *, integer *, real *, real *,
00091 complex *, integer *, complex *, real *, real *);
00092 logical badmm, badnn;
00093 integer nfail, iinfo;
00094 extern int cunt01_(char *, integer *, integer *, complex
00095 *, integer *, complex *, integer *, real *, real *);
00096 real anorm;
00097 extern int cunt03_(char *, integer *, integer *, integer
00098 *, integer *, complex *, integer *, complex *, integer *, complex
00099 *, integer *, real *, real *, integer *);
00100 integer mnmin, mnmax;
00101 char jobvt[1];
00102 integer iwspc, jsize, nerrs, jtype, ntest, iwtmp;
00103 extern int cgesdd_(char *, integer *, integer *, complex
00104 *, integer *, real *, complex *, integer *, complex *, integer *,
00105 complex *, integer *, real *, integer *, integer *);
00106 extern doublereal slamch_(char *);
00107 extern int cgesvd_(char *, char *, integer *, integer *,
00108 complex *, integer *, real *, complex *, integer *, complex *,
00109 integer *, complex *, integer *, real *, integer *), clacpy_(char *, integer *, integer *, complex *, integer
00110 *, complex *, integer *), claset_(char *, integer *,
00111 integer *, complex *, complex *, complex *, integer *);
00112 integer ioldsd[4];
00113 extern int xerbla_(char *, integer *), alasvm_(
00114 char *, integer *, integer *, integer *, integer *),
00115 clatms_(integer *, integer *, char *, integer *, char *, real *,
00116 integer *, real *, real *, integer *, integer *, char *, complex *
00117 , integer *, complex *, integer *);
00118 integer ntestf, minwrk;
00119 real ulpinv, result[14];
00120 integer lswork, mtypes, ntestt;
00121
00122
00123 static cilist io___27 = { 0, 0, 0, fmt_9996, 0 };
00124 static cilist io___32 = { 0, 0, 0, fmt_9995, 0 };
00125 static cilist io___39 = { 0, 0, 0, fmt_9995, 0 };
00126 static cilist io___43 = { 0, 0, 0, fmt_9999, 0 };
00127 static cilist io___44 = { 0, 0, 0, fmt_9998, 0 };
00128 static cilist io___45 = { 0, 0, 0, fmt_9997, 0 };
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199
00200
00201
00202
00203
00204
00205
00206
00207
00208
00209
00210
00211
00212
00213
00214
00215
00216
00217
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232
00233
00234
00235
00236
00237
00238
00239
00240
00241
00242
00243
00244
00245
00246
00247
00248
00249
00250
00251
00252
00253
00254
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265
00266
00267
00268
00269
00270
00271
00272
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337
00338
00339
00340
00341
00342
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365
00366
00367
00368
00369 --mm;
00370 --nn;
00371 --dotype;
00372 --iseed;
00373 asav_dim1 = *lda;
00374 asav_offset = 1 + asav_dim1;
00375 asav -= asav_offset;
00376 a_dim1 = *lda;
00377 a_offset = 1 + a_dim1;
00378 a -= a_offset;
00379 usav_dim1 = *ldu;
00380 usav_offset = 1 + usav_dim1;
00381 usav -= usav_offset;
00382 u_dim1 = *ldu;
00383 u_offset = 1 + u_dim1;
00384 u -= u_offset;
00385 vtsav_dim1 = *ldvt;
00386 vtsav_offset = 1 + vtsav_dim1;
00387 vtsav -= vtsav_offset;
00388 vt_dim1 = *ldvt;
00389 vt_offset = 1 + vt_dim1;
00390 vt -= vt_offset;
00391 --s;
00392 --ssav;
00393 --e;
00394 --work;
00395 --rwork;
00396 --iwork;
00397
00398
00399
00400
00401
00402
00403
00404 *info = 0;
00405
00406
00407
00408 nerrs = 0;
00409 ntestt = 0;
00410 ntestf = 0;
00411 badmm = FALSE_;
00412 badnn = FALSE_;
00413 mmax = 1;
00414 nmax = 1;
00415 mnmax = 1;
00416 minwrk = 1;
00417 i__1 = *nsizes;
00418 for (j = 1; j <= i__1; ++j) {
00419
00420 i__2 = mmax, i__3 = mm[j];
00421 mmax = max(i__2,i__3);
00422 if (mm[j] < 0) {
00423 badmm = TRUE_;
00424 }
00425
00426 i__2 = nmax, i__3 = nn[j];
00427 nmax = max(i__2,i__3);
00428 if (nn[j] < 0) {
00429 badnn = TRUE_;
00430 }
00431
00432
00433 i__4 = mm[j], i__5 = nn[j];
00434 i__2 = mnmax, i__3 = min(i__4,i__5);
00435 mnmax = max(i__2,i__3);
00436
00437
00438
00439 i__6 = mm[j], i__7 = nn[j];
00440
00441 i__9 = mm[j], i__10 = nn[j];
00442
00443 i__8 = max(i__9,i__10);
00444
00445 i__11 = mm[j], i__12 = nn[j];
00446
00447 i__13 = mm[j], i__14 = nn[j];
00448 i__4 = min(i__6,i__7) * 3 + i__8 * i__8, i__5 = min(i__11,i__12) * 5,
00449 i__4 = max(i__4,i__5), i__5 = max(i__13,i__14) * 3;
00450 i__2 = minwrk, i__3 = max(i__4,i__5);
00451 minwrk = max(i__2,i__3);
00452
00453 }
00454
00455
00456
00457 if (*nsizes < 0) {
00458 *info = -1;
00459 } else if (badmm) {
00460 *info = -2;
00461 } else if (badnn) {
00462 *info = -3;
00463 } else if (*ntypes < 0) {
00464 *info = -4;
00465 } else if (*lda < max(1,mmax)) {
00466 *info = -10;
00467 } else if (*ldu < max(1,mmax)) {
00468 *info = -12;
00469 } else if (*ldvt < max(1,nmax)) {
00470 *info = -14;
00471 } else if (minwrk > *lwork) {
00472 *info = -21;
00473 }
00474
00475 if (*info != 0) {
00476 i__1 = -(*info);
00477 xerbla_("CDRVBD", &i__1);
00478 return 0;
00479 }
00480
00481
00482
00483 if (*nsizes == 0 || *ntypes == 0) {
00484 return 0;
00485 }
00486
00487
00488
00489 unfl = slamch_("S");
00490 ovfl = 1.f / unfl;
00491 ulp = slamch_("E");
00492 ulpinv = 1.f / ulp;
00493
00494
00495
00496 nerrs = 0;
00497
00498 i__1 = *nsizes;
00499 for (jsize = 1; jsize <= i__1; ++jsize) {
00500 m = mm[jsize];
00501 n = nn[jsize];
00502 mnmin = min(m,n);
00503
00504 if (*nsizes != 1) {
00505 mtypes = min(5,*ntypes);
00506 } else {
00507 mtypes = min(6,*ntypes);
00508 }
00509
00510 i__2 = mtypes;
00511 for (jtype = 1; jtype <= i__2; ++jtype) {
00512 if (! dotype[jtype]) {
00513 goto L170;
00514 }
00515 ntest = 0;
00516
00517 for (j = 1; j <= 4; ++j) {
00518 ioldsd[j - 1] = iseed[j];
00519
00520 }
00521
00522
00523
00524 if (mtypes > 5) {
00525 goto L50;
00526 }
00527
00528 if (jtype == 1) {
00529
00530
00531
00532 claset_("Full", &m, &n, &c_b1, &c_b1, &a[a_offset], lda);
00533 i__3 = min(m,n);
00534 for (i__ = 1; i__ <= i__3; ++i__) {
00535 s[i__] = 0.f;
00536
00537 }
00538
00539 } else if (jtype == 2) {
00540
00541
00542
00543 claset_("Full", &m, &n, &c_b1, &c_b2, &a[a_offset], lda);
00544 i__3 = min(m,n);
00545 for (i__ = 1; i__ <= i__3; ++i__) {
00546 s[i__] = 1.f;
00547
00548 }
00549
00550 } else {
00551
00552
00553
00554 if (jtype == 3) {
00555 anorm = 1.f;
00556 }
00557 if (jtype == 4) {
00558 anorm = unfl / ulp;
00559 }
00560 if (jtype == 5) {
00561 anorm = ovfl * ulp;
00562 }
00563 r__1 = (real) mnmin;
00564 i__3 = m - 1;
00565 i__4 = n - 1;
00566 clatms_(&m, &n, "U", &iseed[1], "N", &s[1], &c__4, &r__1, &
00567 anorm, &i__3, &i__4, "N", &a[a_offset], lda, &work[1],
00568 &iinfo);
00569 if (iinfo != 0) {
00570 io___27.ciunit = *nounit;
00571 s_wsfe(&io___27);
00572 do_fio(&c__1, "Generator", (ftnlen)9);
00573 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00574 do_fio(&c__1, (char *)&m, (ftnlen)sizeof(integer));
00575 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00576 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00577 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
00578 ;
00579 e_wsfe();
00580 *info = abs(iinfo);
00581 return 0;
00582 }
00583 }
00584
00585 L50:
00586 clacpy_("F", &m, &n, &a[a_offset], lda, &asav[asav_offset], lda);
00587
00588
00589
00590 for (iwspc = 1; iwspc <= 4; ++iwspc) {
00591
00592
00593
00594 iwtmp = (min(m,n) << 1) + max(m,n);
00595 lswork = iwtmp + (iwspc - 1) * (*lwork - iwtmp) / 3;
00596 lswork = min(lswork,*lwork);
00597 lswork = max(lswork,1);
00598 if (iwspc == 4) {
00599 lswork = *lwork;
00600 }
00601
00602 for (j = 1; j <= 14; ++j) {
00603 result[j - 1] = -1.f;
00604
00605 }
00606
00607
00608
00609 if (iwspc > 1) {
00610 clacpy_("F", &m, &n, &asav[asav_offset], lda, &a[a_offset]
00611 , lda);
00612 }
00613 cgesvd_("A", "A", &m, &n, &a[a_offset], lda, &ssav[1], &usav[
00614 usav_offset], ldu, &vtsav[vtsav_offset], ldvt, &work[
00615 1], &lswork, &rwork[1], &iinfo);
00616 if (iinfo != 0) {
00617 io___32.ciunit = *nounit;
00618 s_wsfe(&io___32);
00619 do_fio(&c__1, "GESVD", (ftnlen)5);
00620 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00621 do_fio(&c__1, (char *)&m, (ftnlen)sizeof(integer));
00622 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00623 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00624 do_fio(&c__1, (char *)&lswork, (ftnlen)sizeof(integer));
00625 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
00626 ;
00627 e_wsfe();
00628 *info = abs(iinfo);
00629 return 0;
00630 }
00631
00632
00633
00634 cbdt01_(&m, &n, &c__0, &asav[asav_offset], lda, &usav[
00635 usav_offset], ldu, &ssav[1], &e[1], &vtsav[
00636 vtsav_offset], ldvt, &work[1], &rwork[1], result);
00637 if (m != 0 && n != 0) {
00638 cunt01_("Columns", &mnmin, &m, &usav[usav_offset], ldu, &
00639 work[1], lwork, &rwork[1], &result[1]);
00640 cunt01_("Rows", &mnmin, &n, &vtsav[vtsav_offset], ldvt, &
00641 work[1], lwork, &rwork[1], &result[2]);
00642 }
00643 result[3] = 0.f;
00644 i__3 = mnmin - 1;
00645 for (i__ = 1; i__ <= i__3; ++i__) {
00646 if (ssav[i__] < ssav[i__ + 1]) {
00647 result[3] = ulpinv;
00648 }
00649 if (ssav[i__] < 0.f) {
00650 result[3] = ulpinv;
00651 }
00652
00653 }
00654 if (mnmin >= 1) {
00655 if (ssav[mnmin] < 0.f) {
00656 result[3] = ulpinv;
00657 }
00658 }
00659
00660
00661
00662 result[4] = 0.f;
00663 result[5] = 0.f;
00664 result[6] = 0.f;
00665 for (iju = 0; iju <= 3; ++iju) {
00666 for (ijvt = 0; ijvt <= 3; ++ijvt) {
00667 if (iju == 3 && ijvt == 3 || iju == 1 && ijvt == 1) {
00668 goto L90;
00669 }
00670 *(unsigned char *)jobu = *(unsigned char *)&cjob[iju];
00671 *(unsigned char *)jobvt = *(unsigned char *)&cjob[
00672 ijvt];
00673 clacpy_("F", &m, &n, &asav[asav_offset], lda, &a[
00674 a_offset], lda);
00675 cgesvd_(jobu, jobvt, &m, &n, &a[a_offset], lda, &s[1],
00676 &u[u_offset], ldu, &vt[vt_offset], ldvt, &
00677 work[1], &lswork, &rwork[1], &iinfo);
00678
00679
00680
00681 dif = 0.f;
00682 if (m > 0 && n > 0) {
00683 if (iju == 1) {
00684 cunt03_("C", &m, &mnmin, &m, &mnmin, &usav[
00685 usav_offset], ldu, &a[a_offset], lda,
00686 &work[1], lwork, &rwork[1], &dif, &
00687 iinfo);
00688 } else if (iju == 2) {
00689 cunt03_("C", &m, &mnmin, &m, &mnmin, &usav[
00690 usav_offset], ldu, &u[u_offset], ldu,
00691 &work[1], lwork, &rwork[1], &dif, &
00692 iinfo);
00693 } else if (iju == 3) {
00694 cunt03_("C", &m, &m, &m, &mnmin, &usav[
00695 usav_offset], ldu, &u[u_offset], ldu,
00696 &work[1], lwork, &rwork[1], &dif, &
00697 iinfo);
00698 }
00699 }
00700 result[4] = dmax(result[4],dif);
00701
00702
00703
00704 dif = 0.f;
00705 if (m > 0 && n > 0) {
00706 if (ijvt == 1) {
00707 cunt03_("R", &n, &mnmin, &n, &mnmin, &vtsav[
00708 vtsav_offset], ldvt, &a[a_offset],
00709 lda, &work[1], lwork, &rwork[1], &dif,
00710 &iinfo);
00711 } else if (ijvt == 2) {
00712 cunt03_("R", &n, &mnmin, &n, &mnmin, &vtsav[
00713 vtsav_offset], ldvt, &vt[vt_offset],
00714 ldvt, &work[1], lwork, &rwork[1], &
00715 dif, &iinfo);
00716 } else if (ijvt == 3) {
00717 cunt03_("R", &n, &n, &n, &mnmin, &vtsav[
00718 vtsav_offset], ldvt, &vt[vt_offset],
00719 ldvt, &work[1], lwork, &rwork[1], &
00720 dif, &iinfo);
00721 }
00722 }
00723 result[5] = dmax(result[5],dif);
00724
00725
00726
00727 dif = 0.f;
00728
00729 r__1 = (real) mnmin * ulp * s[1], r__2 = slamch_(
00730 "Safe minimum");
00731 div = dmax(r__1,r__2);
00732 i__3 = mnmin - 1;
00733 for (i__ = 1; i__ <= i__3; ++i__) {
00734 if (ssav[i__] < ssav[i__ + 1]) {
00735 dif = ulpinv;
00736 }
00737 if (ssav[i__] < 0.f) {
00738 dif = ulpinv;
00739 }
00740
00741 r__2 = dif, r__3 = (r__1 = ssav[i__] - s[i__],
00742 dabs(r__1)) / div;
00743 dif = dmax(r__2,r__3);
00744
00745 }
00746 result[6] = dmax(result[6],dif);
00747 L90:
00748 ;
00749 }
00750
00751 }
00752
00753
00754
00755 iwtmp = (mnmin << 1) * mnmin + (mnmin << 1) + max(m,n);
00756 lswork = iwtmp + (iwspc - 1) * (*lwork - iwtmp) / 3;
00757 lswork = min(lswork,*lwork);
00758 lswork = max(lswork,1);
00759 if (iwspc == 4) {
00760 lswork = *lwork;
00761 }
00762
00763
00764
00765 clacpy_("F", &m, &n, &asav[asav_offset], lda, &a[a_offset],
00766 lda);
00767 cgesdd_("A", &m, &n, &a[a_offset], lda, &ssav[1], &usav[
00768 usav_offset], ldu, &vtsav[vtsav_offset], ldvt, &work[
00769 1], &lswork, &rwork[1], &iwork[1], &iinfo);
00770 if (iinfo != 0) {
00771 io___39.ciunit = *nounit;
00772 s_wsfe(&io___39);
00773 do_fio(&c__1, "GESDD", (ftnlen)5);
00774 do_fio(&c__1, (char *)&iinfo, (ftnlen)sizeof(integer));
00775 do_fio(&c__1, (char *)&m, (ftnlen)sizeof(integer));
00776 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00777 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00778 do_fio(&c__1, (char *)&lswork, (ftnlen)sizeof(integer));
00779 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer))
00780 ;
00781 e_wsfe();
00782 *info = abs(iinfo);
00783 return 0;
00784 }
00785
00786
00787
00788 cbdt01_(&m, &n, &c__0, &asav[asav_offset], lda, &usav[
00789 usav_offset], ldu, &ssav[1], &e[1], &vtsav[
00790 vtsav_offset], ldvt, &work[1], &rwork[1], &result[7]);
00791 if (m != 0 && n != 0) {
00792 cunt01_("Columns", &mnmin, &m, &usav[usav_offset], ldu, &
00793 work[1], lwork, &rwork[1], &result[8]);
00794 cunt01_("Rows", &mnmin, &n, &vtsav[vtsav_offset], ldvt, &
00795 work[1], lwork, &rwork[1], &result[9]);
00796 }
00797 result[10] = 0.f;
00798 i__3 = mnmin - 1;
00799 for (i__ = 1; i__ <= i__3; ++i__) {
00800 if (ssav[i__] < ssav[i__ + 1]) {
00801 result[10] = ulpinv;
00802 }
00803 if (ssav[i__] < 0.f) {
00804 result[10] = ulpinv;
00805 }
00806
00807 }
00808 if (mnmin >= 1) {
00809 if (ssav[mnmin] < 0.f) {
00810 result[10] = ulpinv;
00811 }
00812 }
00813
00814
00815
00816 result[11] = 0.f;
00817 result[12] = 0.f;
00818 result[13] = 0.f;
00819 for (ijq = 0; ijq <= 2; ++ijq) {
00820 *(unsigned char *)jobq = *(unsigned char *)&cjob[ijq];
00821 clacpy_("F", &m, &n, &asav[asav_offset], lda, &a[a_offset]
00822 , lda);
00823 cgesdd_(jobq, &m, &n, &a[a_offset], lda, &s[1], &u[
00824 u_offset], ldu, &vt[vt_offset], ldvt, &work[1], &
00825 lswork, &rwork[1], &iwork[1], &iinfo);
00826
00827
00828
00829 dif = 0.f;
00830 if (m > 0 && n > 0) {
00831 if (ijq == 1) {
00832 if (m >= n) {
00833 cunt03_("C", &m, &mnmin, &m, &mnmin, &usav[
00834 usav_offset], ldu, &a[a_offset], lda,
00835 &work[1], lwork, &rwork[1], &dif, &
00836 iinfo);
00837 } else {
00838 cunt03_("C", &m, &mnmin, &m, &mnmin, &usav[
00839 usav_offset], ldu, &u[u_offset], ldu,
00840 &work[1], lwork, &rwork[1], &dif, &
00841 iinfo);
00842 }
00843 } else if (ijq == 2) {
00844 cunt03_("C", &m, &mnmin, &m, &mnmin, &usav[
00845 usav_offset], ldu, &u[u_offset], ldu, &
00846 work[1], lwork, &rwork[1], &dif, &iinfo);
00847 }
00848 }
00849 result[11] = dmax(result[11],dif);
00850
00851
00852
00853 dif = 0.f;
00854 if (m > 0 && n > 0) {
00855 if (ijq == 1) {
00856 if (m >= n) {
00857 cunt03_("R", &n, &mnmin, &n, &mnmin, &vtsav[
00858 vtsav_offset], ldvt, &vt[vt_offset],
00859 ldvt, &work[1], lwork, &rwork[1], &
00860 dif, &iinfo);
00861 } else {
00862 cunt03_("R", &n, &mnmin, &n, &mnmin, &vtsav[
00863 vtsav_offset], ldvt, &a[a_offset],
00864 lda, &work[1], lwork, &rwork[1], &dif,
00865 &iinfo);
00866 }
00867 } else if (ijq == 2) {
00868 cunt03_("R", &n, &mnmin, &n, &mnmin, &vtsav[
00869 vtsav_offset], ldvt, &vt[vt_offset], ldvt,
00870 &work[1], lwork, &rwork[1], &dif, &iinfo);
00871 }
00872 }
00873 result[12] = dmax(result[12],dif);
00874
00875
00876
00877 dif = 0.f;
00878
00879 r__1 = (real) mnmin * ulp * s[1], r__2 = slamch_("Safe m"
00880 "inimum");
00881 div = dmax(r__1,r__2);
00882 i__3 = mnmin - 1;
00883 for (i__ = 1; i__ <= i__3; ++i__) {
00884 if (ssav[i__] < ssav[i__ + 1]) {
00885 dif = ulpinv;
00886 }
00887 if (ssav[i__] < 0.f) {
00888 dif = ulpinv;
00889 }
00890
00891 r__2 = dif, r__3 = (r__1 = ssav[i__] - s[i__], dabs(
00892 r__1)) / div;
00893 dif = dmax(r__2,r__3);
00894
00895 }
00896 result[13] = dmax(result[13],dif);
00897
00898 }
00899
00900
00901
00902 ntest = 0;
00903 nfail = 0;
00904 for (j = 1; j <= 14; ++j) {
00905 if (result[j - 1] >= 0.f) {
00906 ++ntest;
00907 }
00908 if (result[j - 1] >= *thresh) {
00909 ++nfail;
00910 }
00911
00912 }
00913
00914 if (nfail > 0) {
00915 ++ntestf;
00916 }
00917 if (ntestf == 1) {
00918 io___43.ciunit = *nounit;
00919 s_wsfe(&io___43);
00920 e_wsfe();
00921 io___44.ciunit = *nounit;
00922 s_wsfe(&io___44);
00923 do_fio(&c__1, (char *)&(*thresh), (ftnlen)sizeof(real));
00924 e_wsfe();
00925 ntestf = 2;
00926 }
00927
00928 for (j = 1; j <= 14; ++j) {
00929 if (result[j - 1] >= *thresh) {
00930 io___45.ciunit = *nounit;
00931 s_wsfe(&io___45);
00932 do_fio(&c__1, (char *)&m, (ftnlen)sizeof(integer));
00933 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00934 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer))
00935 ;
00936 do_fio(&c__1, (char *)&iwspc, (ftnlen)sizeof(integer))
00937 ;
00938 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
00939 integer));
00940 do_fio(&c__1, (char *)&j, (ftnlen)sizeof(integer));
00941 do_fio(&c__1, (char *)&result[j - 1], (ftnlen)sizeof(
00942 real));
00943 e_wsfe();
00944 }
00945
00946 }
00947
00948 nerrs += nfail;
00949 ntestt += ntest;
00950
00951
00952 }
00953
00954 L170:
00955 ;
00956 }
00957
00958 }
00959
00960
00961
00962 alasvm_("CBD", nounit, &nerrs, &ntestt, &c__0);
00963
00964
00965 return 0;
00966
00967
00968
00969 }