strsyl.c
Go to the documentation of this file.
00001 /* strsyl.f -- translated by f2c (version 20061008).
00002    You must link the resulting object file with libf2c:
00003         on Microsoft Windows system, link with libf2c.lib;
00004         on Linux or Unix systems, link with .../path/to/libf2c.a -lm
00005         or, if you install libf2c.a in a standard place, with -lf2c -lm
00006         -- in that order, at the end of the command line, as in
00007                 cc *.o -lf2c -lm
00008         Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
00009 
00010                 http://www.netlib.org/f2c/libf2c.zip
00011 */
00012 
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015 
00016 /* Table of constant values */
00017 
00018 static integer c__1 = 1;
00019 static logical c_false = FALSE_;
00020 static integer c__2 = 2;
00021 static real c_b26 = 1.f;
00022 static real c_b30 = 0.f;
00023 static logical c_true = TRUE_;
00024 
00025 /* Subroutine */ int strsyl_(char *trana, char *tranb, integer *isgn, integer 
00026         *m, integer *n, real *a, integer *lda, real *b, integer *ldb, real *
00027         c__, integer *ldc, real *scale, integer *info)
00028 {
00029     /* System generated locals */
00030     integer a_dim1, a_offset, b_dim1, b_offset, c_dim1, c_offset, i__1, i__2, 
00031             i__3, i__4;
00032     real r__1, r__2;
00033 
00034     /* Local variables */
00035     integer j, k, l;
00036     real x[4]   /* was [2][2] */;
00037     integer k1, k2, l1, l2;
00038     real a11, db, da11, vec[4]  /* was [2][2] */, dum[1], eps, sgn;
00039     integer ierr;
00040     real smin;
00041     extern doublereal sdot_(integer *, real *, integer *, real *, integer *);
00042     real suml, sumr;
00043     extern logical lsame_(char *, char *);
00044     extern /* Subroutine */ int sscal_(integer *, real *, real *, integer *);
00045     integer knext, lnext;
00046     real xnorm;
00047     extern /* Subroutine */ int slaln2_(logical *, integer *, integer *, real 
00048             *, real *, real *, integer *, real *, real *, real *, integer *, 
00049             real *, real *, real *, integer *, real *, real *, integer *), 
00050             slasy2_(logical *, logical *, integer *, integer *, integer *, 
00051             real *, integer *, real *, integer *, real *, integer *, real *, 
00052             real *, integer *, real *, integer *), slabad_(real *, real *);
00053     real scaloc;
00054     extern doublereal slamch_(char *), slange_(char *, integer *, 
00055             integer *, real *, integer *, real *);
00056     extern /* Subroutine */ int xerbla_(char *, integer *);
00057     real bignum;
00058     logical notrna, notrnb;
00059     real smlnum;
00060 
00061 
00062 /*  -- LAPACK routine (version 3.2) -- */
00063 /*  -- LAPACK is a software package provided by Univ. of Tennessee,    -- */
00064 /*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- */
00065 /*     November 2006 */
00066 
00067 /*     .. Scalar Arguments .. */
00068 /*     .. */
00069 /*     .. Array Arguments .. */
00070 /*     .. */
00071 
00072 /*  Purpose */
00073 /*  ======= */
00074 
00075 /*  STRSYL solves the real Sylvester matrix equation: */
00076 
00077 /*     op(A)*X + X*op(B) = scale*C or */
00078 /*     op(A)*X - X*op(B) = scale*C, */
00079 
00080 /*  where op(A) = A or A**T, and  A and B are both upper quasi- */
00081 /*  triangular. A is M-by-M and B is N-by-N; the right hand side C and */
00082 /*  the solution X are M-by-N; and scale is an output scale factor, set */
00083 /*  <= 1 to avoid overflow in X. */
00084 
00085 /*  A and B must be in Schur canonical form (as returned by SHSEQR), that */
00086 /*  is, block upper triangular with 1-by-1 and 2-by-2 diagonal blocks; */
00087 /*  each 2-by-2 diagonal block has its diagonal elements equal and its */
00088 /*  off-diagonal elements of opposite sign. */
00089 
00090 /*  Arguments */
00091 /*  ========= */
00092 
00093 /*  TRANA   (input) CHARACTER*1 */
00094 /*          Specifies the option op(A): */
00095 /*          = 'N': op(A) = A    (No transpose) */
00096 /*          = 'T': op(A) = A**T (Transpose) */
00097 /*          = 'C': op(A) = A**H (Conjugate transpose = Transpose) */
00098 
00099 /*  TRANB   (input) CHARACTER*1 */
00100 /*          Specifies the option op(B): */
00101 /*          = 'N': op(B) = B    (No transpose) */
00102 /*          = 'T': op(B) = B**T (Transpose) */
00103 /*          = 'C': op(B) = B**H (Conjugate transpose = Transpose) */
00104 
00105 /*  ISGN    (input) INTEGER */
00106 /*          Specifies the sign in the equation: */
00107 /*          = +1: solve op(A)*X + X*op(B) = scale*C */
00108 /*          = -1: solve op(A)*X - X*op(B) = scale*C */
00109 
00110 /*  M       (input) INTEGER */
00111 /*          The order of the matrix A, and the number of rows in the */
00112 /*          matrices X and C. M >= 0. */
00113 
00114 /*  N       (input) INTEGER */
00115 /*          The order of the matrix B, and the number of columns in the */
00116 /*          matrices X and C. N >= 0. */
00117 
00118 /*  A       (input) REAL array, dimension (LDA,M) */
00119 /*          The upper quasi-triangular matrix A, in Schur canonical form. */
00120 
00121 /*  LDA     (input) INTEGER */
00122 /*          The leading dimension of the array A. LDA >= max(1,M). */
00123 
00124 /*  B       (input) REAL array, dimension (LDB,N) */
00125 /*          The upper quasi-triangular matrix B, in Schur canonical form. */
00126 
00127 /*  LDB     (input) INTEGER */
00128 /*          The leading dimension of the array B. LDB >= max(1,N). */
00129 
00130 /*  C       (input/output) REAL array, dimension (LDC,N) */
00131 /*          On entry, the M-by-N right hand side matrix C. */
00132 /*          On exit, C is overwritten by the solution matrix X. */
00133 
00134 /*  LDC     (input) INTEGER */
00135 /*          The leading dimension of the array C. LDC >= max(1,M) */
00136 
00137 /*  SCALE   (output) REAL */
00138 /*          The scale factor, scale, set <= 1 to avoid overflow in X. */
00139 
00140 /*  INFO    (output) INTEGER */
00141 /*          = 0: successful exit */
00142 /*          < 0: if INFO = -i, the i-th argument had an illegal value */
00143 /*          = 1: A and B have common or very close eigenvalues; perturbed */
00144 /*               values were used to solve the equation (but the matrices */
00145 /*               A and B are unchanged). */
00146 
00147 /*  ===================================================================== */
00148 
00149 /*     .. Parameters .. */
00150 /*     .. */
00151 /*     .. Local Scalars .. */
00152 /*     .. */
00153 /*     .. Local Arrays .. */
00154 /*     .. */
00155 /*     .. External Functions .. */
00156 /*     .. */
00157 /*     .. External Subroutines .. */
00158 /*     .. */
00159 /*     .. Intrinsic Functions .. */
00160 /*     .. */
00161 /*     .. Executable Statements .. */
00162 
00163 /*     Decode and Test input parameters */
00164 
00165     /* Parameter adjustments */
00166     a_dim1 = *lda;
00167     a_offset = 1 + a_dim1;
00168     a -= a_offset;
00169     b_dim1 = *ldb;
00170     b_offset = 1 + b_dim1;
00171     b -= b_offset;
00172     c_dim1 = *ldc;
00173     c_offset = 1 + c_dim1;
00174     c__ -= c_offset;
00175 
00176     /* Function Body */
00177     notrna = lsame_(trana, "N");
00178     notrnb = lsame_(tranb, "N");
00179 
00180     *info = 0;
00181     if (! notrna && ! lsame_(trana, "T") && ! lsame_(
00182             trana, "C")) {
00183         *info = -1;
00184     } else if (! notrnb && ! lsame_(tranb, "T") && ! 
00185             lsame_(tranb, "C")) {
00186         *info = -2;
00187     } else if (*isgn != 1 && *isgn != -1) {
00188         *info = -3;
00189     } else if (*m < 0) {
00190         *info = -4;
00191     } else if (*n < 0) {
00192         *info = -5;
00193     } else if (*lda < max(1,*m)) {
00194         *info = -7;
00195     } else if (*ldb < max(1,*n)) {
00196         *info = -9;
00197     } else if (*ldc < max(1,*m)) {
00198         *info = -11;
00199     }
00200     if (*info != 0) {
00201         i__1 = -(*info);
00202         xerbla_("STRSYL", &i__1);
00203         return 0;
00204     }
00205 
00206 /*     Quick return if possible */
00207 
00208     *scale = 1.f;
00209     if (*m == 0 || *n == 0) {
00210         return 0;
00211     }
00212 
00213 /*     Set constants to control overflow */
00214 
00215     eps = slamch_("P");
00216     smlnum = slamch_("S");
00217     bignum = 1.f / smlnum;
00218     slabad_(&smlnum, &bignum);
00219     smlnum = smlnum * (real) (*m * *n) / eps;
00220     bignum = 1.f / smlnum;
00221 
00222 /* Computing MAX */
00223     r__1 = smlnum, r__2 = eps * slange_("M", m, m, &a[a_offset], lda, dum), r__1 = max(r__1,r__2), r__2 = eps * slange_("M", n, n, 
00224             &b[b_offset], ldb, dum);
00225     smin = dmax(r__1,r__2);
00226 
00227     sgn = (real) (*isgn);
00228 
00229     if (notrna && notrnb) {
00230 
00231 /*        Solve    A*X + ISGN*X*B = scale*C. */
00232 
00233 /*        The (K,L)th block of X is determined starting from */
00234 /*        bottom-left corner column by column by */
00235 
00236 /*         A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L) */
00237 
00238 /*        Where */
00239 /*                  M                         L-1 */
00240 /*        R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)]. */
00241 /*                I=K+1                       J=1 */
00242 
00243 /*        Start column loop (index = L) */
00244 /*        L1 (L2) : column index of the first (first) row of X(K,L). */
00245 
00246         lnext = 1;
00247         i__1 = *n;
00248         for (l = 1; l <= i__1; ++l) {
00249             if (l < lnext) {
00250                 goto L70;
00251             }
00252             if (l == *n) {
00253                 l1 = l;
00254                 l2 = l;
00255             } else {
00256                 if (b[l + 1 + l * b_dim1] != 0.f) {
00257                     l1 = l;
00258                     l2 = l + 1;
00259                     lnext = l + 2;
00260                 } else {
00261                     l1 = l;
00262                     l2 = l;
00263                     lnext = l + 1;
00264                 }
00265             }
00266 
00267 /*           Start row loop (index = K) */
00268 /*           K1 (K2): row index of the first (last) row of X(K,L). */
00269 
00270             knext = *m;
00271             for (k = *m; k >= 1; --k) {
00272                 if (k > knext) {
00273                     goto L60;
00274                 }
00275                 if (k == 1) {
00276                     k1 = k;
00277                     k2 = k;
00278                 } else {
00279                     if (a[k + (k - 1) * a_dim1] != 0.f) {
00280                         k1 = k - 1;
00281                         k2 = k;
00282                         knext = k - 2;
00283                     } else {
00284                         k1 = k;
00285                         k2 = k;
00286                         knext = k - 1;
00287                     }
00288                 }
00289 
00290                 if (l1 == l2 && k1 == k2) {
00291                     i__2 = *m - k1;
00292 /* Computing MIN */
00293                     i__3 = k1 + 1;
00294 /* Computing MIN */
00295                     i__4 = k1 + 1;
00296                     suml = sdot_(&i__2, &a[k1 + min(i__3, *m)* a_dim1], lda, &
00297                             c__[min(i__4, *m)+ l1 * c_dim1], &c__1);
00298                     i__2 = l1 - 1;
00299                     sumr = sdot_(&i__2, &c__[k1 + c_dim1], ldc, &b[l1 * 
00300                             b_dim1 + 1], &c__1);
00301                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00302                     scaloc = 1.f;
00303 
00304                     a11 = a[k1 + k1 * a_dim1] + sgn * b[l1 + l1 * b_dim1];
00305                     da11 = dabs(a11);
00306                     if (da11 <= smin) {
00307                         a11 = smin;
00308                         da11 = smin;
00309                         *info = 1;
00310                     }
00311                     db = dabs(vec[0]);
00312                     if (da11 < 1.f && db > 1.f) {
00313                         if (db > bignum * da11) {
00314                             scaloc = 1.f / db;
00315                         }
00316                     }
00317                     x[0] = vec[0] * scaloc / a11;
00318 
00319                     if (scaloc != 1.f) {
00320                         i__2 = *n;
00321                         for (j = 1; j <= i__2; ++j) {
00322                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00323 /* L10: */
00324                         }
00325                         *scale *= scaloc;
00326                     }
00327                     c__[k1 + l1 * c_dim1] = x[0];
00328 
00329                 } else if (l1 == l2 && k1 != k2) {
00330 
00331                     i__2 = *m - k2;
00332 /* Computing MIN */
00333                     i__3 = k2 + 1;
00334 /* Computing MIN */
00335                     i__4 = k2 + 1;
00336                     suml = sdot_(&i__2, &a[k1 + min(i__3, *m)* a_dim1], lda, &
00337                             c__[min(i__4, *m)+ l1 * c_dim1], &c__1);
00338                     i__2 = l1 - 1;
00339                     sumr = sdot_(&i__2, &c__[k1 + c_dim1], ldc, &b[l1 * 
00340                             b_dim1 + 1], &c__1);
00341                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00342 
00343                     i__2 = *m - k2;
00344 /* Computing MIN */
00345                     i__3 = k2 + 1;
00346 /* Computing MIN */
00347                     i__4 = k2 + 1;
00348                     suml = sdot_(&i__2, &a[k2 + min(i__3, *m)* a_dim1], lda, &
00349                             c__[min(i__4, *m)+ l1 * c_dim1], &c__1);
00350                     i__2 = l1 - 1;
00351                     sumr = sdot_(&i__2, &c__[k2 + c_dim1], ldc, &b[l1 * 
00352                             b_dim1 + 1], &c__1);
00353                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
00354 
00355                     r__1 = -sgn * b[l1 + l1 * b_dim1];
00356                     slaln2_(&c_false, &c__2, &c__1, &smin, &c_b26, &a[k1 + k1 
00357                             * a_dim1], lda, &c_b26, &c_b26, vec, &c__2, &r__1, 
00358                              &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
00359                     if (ierr != 0) {
00360                         *info = 1;
00361                     }
00362 
00363                     if (scaloc != 1.f) {
00364                         i__2 = *n;
00365                         for (j = 1; j <= i__2; ++j) {
00366                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00367 /* L20: */
00368                         }
00369                         *scale *= scaloc;
00370                     }
00371                     c__[k1 + l1 * c_dim1] = x[0];
00372                     c__[k2 + l1 * c_dim1] = x[1];
00373 
00374                 } else if (l1 != l2 && k1 == k2) {
00375 
00376                     i__2 = *m - k1;
00377 /* Computing MIN */
00378                     i__3 = k1 + 1;
00379 /* Computing MIN */
00380                     i__4 = k1 + 1;
00381                     suml = sdot_(&i__2, &a[k1 + min(i__3, *m)* a_dim1], lda, &
00382                             c__[min(i__4, *m)+ l1 * c_dim1], &c__1);
00383                     i__2 = l1 - 1;
00384                     sumr = sdot_(&i__2, &c__[k1 + c_dim1], ldc, &b[l1 * 
00385                             b_dim1 + 1], &c__1);
00386                     vec[0] = sgn * (c__[k1 + l1 * c_dim1] - (suml + sgn * 
00387                             sumr));
00388 
00389                     i__2 = *m - k1;
00390 /* Computing MIN */
00391                     i__3 = k1 + 1;
00392 /* Computing MIN */
00393                     i__4 = k1 + 1;
00394                     suml = sdot_(&i__2, &a[k1 + min(i__3, *m)* a_dim1], lda, &
00395                             c__[min(i__4, *m)+ l2 * c_dim1], &c__1);
00396                     i__2 = l1 - 1;
00397                     sumr = sdot_(&i__2, &c__[k1 + c_dim1], ldc, &b[l2 * 
00398                             b_dim1 + 1], &c__1);
00399                     vec[1] = sgn * (c__[k1 + l2 * c_dim1] - (suml + sgn * 
00400                             sumr));
00401 
00402                     r__1 = -sgn * a[k1 + k1 * a_dim1];
00403                     slaln2_(&c_true, &c__2, &c__1, &smin, &c_b26, &b[l1 + l1 *
00404                              b_dim1], ldb, &c_b26, &c_b26, vec, &c__2, &r__1, 
00405                             &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
00406                     if (ierr != 0) {
00407                         *info = 1;
00408                     }
00409 
00410                     if (scaloc != 1.f) {
00411                         i__2 = *n;
00412                         for (j = 1; j <= i__2; ++j) {
00413                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00414 /* L40: */
00415                         }
00416                         *scale *= scaloc;
00417                     }
00418                     c__[k1 + l1 * c_dim1] = x[0];
00419                     c__[k1 + l2 * c_dim1] = x[1];
00420 
00421                 } else if (l1 != l2 && k1 != k2) {
00422 
00423                     i__2 = *m - k2;
00424 /* Computing MIN */
00425                     i__3 = k2 + 1;
00426 /* Computing MIN */
00427                     i__4 = k2 + 1;
00428                     suml = sdot_(&i__2, &a[k1 + min(i__3, *m)* a_dim1], lda, &
00429                             c__[min(i__4, *m)+ l1 * c_dim1], &c__1);
00430                     i__2 = l1 - 1;
00431                     sumr = sdot_(&i__2, &c__[k1 + c_dim1], ldc, &b[l1 * 
00432                             b_dim1 + 1], &c__1);
00433                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00434 
00435                     i__2 = *m - k2;
00436 /* Computing MIN */
00437                     i__3 = k2 + 1;
00438 /* Computing MIN */
00439                     i__4 = k2 + 1;
00440                     suml = sdot_(&i__2, &a[k1 + min(i__3, *m)* a_dim1], lda, &
00441                             c__[min(i__4, *m)+ l2 * c_dim1], &c__1);
00442                     i__2 = l1 - 1;
00443                     sumr = sdot_(&i__2, &c__[k1 + c_dim1], ldc, &b[l2 * 
00444                             b_dim1 + 1], &c__1);
00445                     vec[2] = c__[k1 + l2 * c_dim1] - (suml + sgn * sumr);
00446 
00447                     i__2 = *m - k2;
00448 /* Computing MIN */
00449                     i__3 = k2 + 1;
00450 /* Computing MIN */
00451                     i__4 = k2 + 1;
00452                     suml = sdot_(&i__2, &a[k2 + min(i__3, *m)* a_dim1], lda, &
00453                             c__[min(i__4, *m)+ l1 * c_dim1], &c__1);
00454                     i__2 = l1 - 1;
00455                     sumr = sdot_(&i__2, &c__[k2 + c_dim1], ldc, &b[l1 * 
00456                             b_dim1 + 1], &c__1);
00457                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
00458 
00459                     i__2 = *m - k2;
00460 /* Computing MIN */
00461                     i__3 = k2 + 1;
00462 /* Computing MIN */
00463                     i__4 = k2 + 1;
00464                     suml = sdot_(&i__2, &a[k2 + min(i__3, *m)* a_dim1], lda, &
00465                             c__[min(i__4, *m)+ l2 * c_dim1], &c__1);
00466                     i__2 = l1 - 1;
00467                     sumr = sdot_(&i__2, &c__[k2 + c_dim1], ldc, &b[l2 * 
00468                             b_dim1 + 1], &c__1);
00469                     vec[3] = c__[k2 + l2 * c_dim1] - (suml + sgn * sumr);
00470 
00471                     slasy2_(&c_false, &c_false, isgn, &c__2, &c__2, &a[k1 + 
00472                             k1 * a_dim1], lda, &b[l1 + l1 * b_dim1], ldb, vec, 
00473                              &c__2, &scaloc, x, &c__2, &xnorm, &ierr);
00474                     if (ierr != 0) {
00475                         *info = 1;
00476                     }
00477 
00478                     if (scaloc != 1.f) {
00479                         i__2 = *n;
00480                         for (j = 1; j <= i__2; ++j) {
00481                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00482 /* L50: */
00483                         }
00484                         *scale *= scaloc;
00485                     }
00486                     c__[k1 + l1 * c_dim1] = x[0];
00487                     c__[k1 + l2 * c_dim1] = x[2];
00488                     c__[k2 + l1 * c_dim1] = x[1];
00489                     c__[k2 + l2 * c_dim1] = x[3];
00490                 }
00491 
00492 L60:
00493                 ;
00494             }
00495 
00496 L70:
00497             ;
00498         }
00499 
00500     } else if (! notrna && notrnb) {
00501 
00502 /*        Solve    A' *X + ISGN*X*B = scale*C. */
00503 
00504 /*        The (K,L)th block of X is determined starting from */
00505 /*        upper-left corner column by column by */
00506 
00507 /*          A(K,K)'*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L) */
00508 
00509 /*        Where */
00510 /*                   K-1                        L-1 */
00511 /*          R(K,L) = SUM [A(I,K)'*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)] */
00512 /*                   I=1                        J=1 */
00513 
00514 /*        Start column loop (index = L) */
00515 /*        L1 (L2): column index of the first (last) row of X(K,L) */
00516 
00517         lnext = 1;
00518         i__1 = *n;
00519         for (l = 1; l <= i__1; ++l) {
00520             if (l < lnext) {
00521                 goto L130;
00522             }
00523             if (l == *n) {
00524                 l1 = l;
00525                 l2 = l;
00526             } else {
00527                 if (b[l + 1 + l * b_dim1] != 0.f) {
00528                     l1 = l;
00529                     l2 = l + 1;
00530                     lnext = l + 2;
00531                 } else {
00532                     l1 = l;
00533                     l2 = l;
00534                     lnext = l + 1;
00535                 }
00536             }
00537 
00538 /*           Start row loop (index = K) */
00539 /*           K1 (K2): row index of the first (last) row of X(K,L) */
00540 
00541             knext = 1;
00542             i__2 = *m;
00543             for (k = 1; k <= i__2; ++k) {
00544                 if (k < knext) {
00545                     goto L120;
00546                 }
00547                 if (k == *m) {
00548                     k1 = k;
00549                     k2 = k;
00550                 } else {
00551                     if (a[k + 1 + k * a_dim1] != 0.f) {
00552                         k1 = k;
00553                         k2 = k + 1;
00554                         knext = k + 2;
00555                     } else {
00556                         k1 = k;
00557                         k2 = k;
00558                         knext = k + 1;
00559                     }
00560                 }
00561 
00562                 if (l1 == l2 && k1 == k2) {
00563                     i__3 = k1 - 1;
00564                     suml = sdot_(&i__3, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00565                             c_dim1 + 1], &c__1);
00566                     i__3 = l1 - 1;
00567                     sumr = sdot_(&i__3, &c__[k1 + c_dim1], ldc, &b[l1 * 
00568                             b_dim1 + 1], &c__1);
00569                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00570                     scaloc = 1.f;
00571 
00572                     a11 = a[k1 + k1 * a_dim1] + sgn * b[l1 + l1 * b_dim1];
00573                     da11 = dabs(a11);
00574                     if (da11 <= smin) {
00575                         a11 = smin;
00576                         da11 = smin;
00577                         *info = 1;
00578                     }
00579                     db = dabs(vec[0]);
00580                     if (da11 < 1.f && db > 1.f) {
00581                         if (db > bignum * da11) {
00582                             scaloc = 1.f / db;
00583                         }
00584                     }
00585                     x[0] = vec[0] * scaloc / a11;
00586 
00587                     if (scaloc != 1.f) {
00588                         i__3 = *n;
00589                         for (j = 1; j <= i__3; ++j) {
00590                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00591 /* L80: */
00592                         }
00593                         *scale *= scaloc;
00594                     }
00595                     c__[k1 + l1 * c_dim1] = x[0];
00596 
00597                 } else if (l1 == l2 && k1 != k2) {
00598 
00599                     i__3 = k1 - 1;
00600                     suml = sdot_(&i__3, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00601                             c_dim1 + 1], &c__1);
00602                     i__3 = l1 - 1;
00603                     sumr = sdot_(&i__3, &c__[k1 + c_dim1], ldc, &b[l1 * 
00604                             b_dim1 + 1], &c__1);
00605                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00606 
00607                     i__3 = k1 - 1;
00608                     suml = sdot_(&i__3, &a[k2 * a_dim1 + 1], &c__1, &c__[l1 * 
00609                             c_dim1 + 1], &c__1);
00610                     i__3 = l1 - 1;
00611                     sumr = sdot_(&i__3, &c__[k2 + c_dim1], ldc, &b[l1 * 
00612                             b_dim1 + 1], &c__1);
00613                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
00614 
00615                     r__1 = -sgn * b[l1 + l1 * b_dim1];
00616                     slaln2_(&c_true, &c__2, &c__1, &smin, &c_b26, &a[k1 + k1 *
00617                              a_dim1], lda, &c_b26, &c_b26, vec, &c__2, &r__1, 
00618                             &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
00619                     if (ierr != 0) {
00620                         *info = 1;
00621                     }
00622 
00623                     if (scaloc != 1.f) {
00624                         i__3 = *n;
00625                         for (j = 1; j <= i__3; ++j) {
00626                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00627 /* L90: */
00628                         }
00629                         *scale *= scaloc;
00630                     }
00631                     c__[k1 + l1 * c_dim1] = x[0];
00632                     c__[k2 + l1 * c_dim1] = x[1];
00633 
00634                 } else if (l1 != l2 && k1 == k2) {
00635 
00636                     i__3 = k1 - 1;
00637                     suml = sdot_(&i__3, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00638                             c_dim1 + 1], &c__1);
00639                     i__3 = l1 - 1;
00640                     sumr = sdot_(&i__3, &c__[k1 + c_dim1], ldc, &b[l1 * 
00641                             b_dim1 + 1], &c__1);
00642                     vec[0] = sgn * (c__[k1 + l1 * c_dim1] - (suml + sgn * 
00643                             sumr));
00644 
00645                     i__3 = k1 - 1;
00646                     suml = sdot_(&i__3, &a[k1 * a_dim1 + 1], &c__1, &c__[l2 * 
00647                             c_dim1 + 1], &c__1);
00648                     i__3 = l1 - 1;
00649                     sumr = sdot_(&i__3, &c__[k1 + c_dim1], ldc, &b[l2 * 
00650                             b_dim1 + 1], &c__1);
00651                     vec[1] = sgn * (c__[k1 + l2 * c_dim1] - (suml + sgn * 
00652                             sumr));
00653 
00654                     r__1 = -sgn * a[k1 + k1 * a_dim1];
00655                     slaln2_(&c_true, &c__2, &c__1, &smin, &c_b26, &b[l1 + l1 *
00656                              b_dim1], ldb, &c_b26, &c_b26, vec, &c__2, &r__1, 
00657                             &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
00658                     if (ierr != 0) {
00659                         *info = 1;
00660                     }
00661 
00662                     if (scaloc != 1.f) {
00663                         i__3 = *n;
00664                         for (j = 1; j <= i__3; ++j) {
00665                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00666 /* L100: */
00667                         }
00668                         *scale *= scaloc;
00669                     }
00670                     c__[k1 + l1 * c_dim1] = x[0];
00671                     c__[k1 + l2 * c_dim1] = x[1];
00672 
00673                 } else if (l1 != l2 && k1 != k2) {
00674 
00675                     i__3 = k1 - 1;
00676                     suml = sdot_(&i__3, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00677                             c_dim1 + 1], &c__1);
00678                     i__3 = l1 - 1;
00679                     sumr = sdot_(&i__3, &c__[k1 + c_dim1], ldc, &b[l1 * 
00680                             b_dim1 + 1], &c__1);
00681                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00682 
00683                     i__3 = k1 - 1;
00684                     suml = sdot_(&i__3, &a[k1 * a_dim1 + 1], &c__1, &c__[l2 * 
00685                             c_dim1 + 1], &c__1);
00686                     i__3 = l1 - 1;
00687                     sumr = sdot_(&i__3, &c__[k1 + c_dim1], ldc, &b[l2 * 
00688                             b_dim1 + 1], &c__1);
00689                     vec[2] = c__[k1 + l2 * c_dim1] - (suml + sgn * sumr);
00690 
00691                     i__3 = k1 - 1;
00692                     suml = sdot_(&i__3, &a[k2 * a_dim1 + 1], &c__1, &c__[l1 * 
00693                             c_dim1 + 1], &c__1);
00694                     i__3 = l1 - 1;
00695                     sumr = sdot_(&i__3, &c__[k2 + c_dim1], ldc, &b[l1 * 
00696                             b_dim1 + 1], &c__1);
00697                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
00698 
00699                     i__3 = k1 - 1;
00700                     suml = sdot_(&i__3, &a[k2 * a_dim1 + 1], &c__1, &c__[l2 * 
00701                             c_dim1 + 1], &c__1);
00702                     i__3 = l1 - 1;
00703                     sumr = sdot_(&i__3, &c__[k2 + c_dim1], ldc, &b[l2 * 
00704                             b_dim1 + 1], &c__1);
00705                     vec[3] = c__[k2 + l2 * c_dim1] - (suml + sgn * sumr);
00706 
00707                     slasy2_(&c_true, &c_false, isgn, &c__2, &c__2, &a[k1 + k1 
00708                             * a_dim1], lda, &b[l1 + l1 * b_dim1], ldb, vec, &
00709                             c__2, &scaloc, x, &c__2, &xnorm, &ierr);
00710                     if (ierr != 0) {
00711                         *info = 1;
00712                     }
00713 
00714                     if (scaloc != 1.f) {
00715                         i__3 = *n;
00716                         for (j = 1; j <= i__3; ++j) {
00717                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00718 /* L110: */
00719                         }
00720                         *scale *= scaloc;
00721                     }
00722                     c__[k1 + l1 * c_dim1] = x[0];
00723                     c__[k1 + l2 * c_dim1] = x[2];
00724                     c__[k2 + l1 * c_dim1] = x[1];
00725                     c__[k2 + l2 * c_dim1] = x[3];
00726                 }
00727 
00728 L120:
00729                 ;
00730             }
00731 L130:
00732             ;
00733         }
00734 
00735     } else if (! notrna && ! notrnb) {
00736 
00737 /*        Solve    A'*X + ISGN*X*B' = scale*C. */
00738 
00739 /*        The (K,L)th block of X is determined starting from */
00740 /*        top-right corner column by column by */
00741 
00742 /*           A(K,K)'*X(K,L) + ISGN*X(K,L)*B(L,L)' = C(K,L) - R(K,L) */
00743 
00744 /*        Where */
00745 /*                     K-1                          N */
00746 /*            R(K,L) = SUM [A(I,K)'*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)']. */
00747 /*                     I=1                        J=L+1 */
00748 
00749 /*        Start column loop (index = L) */
00750 /*        L1 (L2): column index of the first (last) row of X(K,L) */
00751 
00752         lnext = *n;
00753         for (l = *n; l >= 1; --l) {
00754             if (l > lnext) {
00755                 goto L190;
00756             }
00757             if (l == 1) {
00758                 l1 = l;
00759                 l2 = l;
00760             } else {
00761                 if (b[l + (l - 1) * b_dim1] != 0.f) {
00762                     l1 = l - 1;
00763                     l2 = l;
00764                     lnext = l - 2;
00765                 } else {
00766                     l1 = l;
00767                     l2 = l;
00768                     lnext = l - 1;
00769                 }
00770             }
00771 
00772 /*           Start row loop (index = K) */
00773 /*           K1 (K2): row index of the first (last) row of X(K,L) */
00774 
00775             knext = 1;
00776             i__1 = *m;
00777             for (k = 1; k <= i__1; ++k) {
00778                 if (k < knext) {
00779                     goto L180;
00780                 }
00781                 if (k == *m) {
00782                     k1 = k;
00783                     k2 = k;
00784                 } else {
00785                     if (a[k + 1 + k * a_dim1] != 0.f) {
00786                         k1 = k;
00787                         k2 = k + 1;
00788                         knext = k + 2;
00789                     } else {
00790                         k1 = k;
00791                         k2 = k;
00792                         knext = k + 1;
00793                     }
00794                 }
00795 
00796                 if (l1 == l2 && k1 == k2) {
00797                     i__2 = k1 - 1;
00798                     suml = sdot_(&i__2, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00799                             c_dim1 + 1], &c__1);
00800                     i__2 = *n - l1;
00801 /* Computing MIN */
00802                     i__3 = l1 + 1;
00803 /* Computing MIN */
00804                     i__4 = l1 + 1;
00805                     sumr = sdot_(&i__2, &c__[k1 + min(i__3, *n)* c_dim1], ldc, 
00806                              &b[l1 + min(i__4, *n)* b_dim1], ldb);
00807                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00808                     scaloc = 1.f;
00809 
00810                     a11 = a[k1 + k1 * a_dim1] + sgn * b[l1 + l1 * b_dim1];
00811                     da11 = dabs(a11);
00812                     if (da11 <= smin) {
00813                         a11 = smin;
00814                         da11 = smin;
00815                         *info = 1;
00816                     }
00817                     db = dabs(vec[0]);
00818                     if (da11 < 1.f && db > 1.f) {
00819                         if (db > bignum * da11) {
00820                             scaloc = 1.f / db;
00821                         }
00822                     }
00823                     x[0] = vec[0] * scaloc / a11;
00824 
00825                     if (scaloc != 1.f) {
00826                         i__2 = *n;
00827                         for (j = 1; j <= i__2; ++j) {
00828                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00829 /* L140: */
00830                         }
00831                         *scale *= scaloc;
00832                     }
00833                     c__[k1 + l1 * c_dim1] = x[0];
00834 
00835                 } else if (l1 == l2 && k1 != k2) {
00836 
00837                     i__2 = k1 - 1;
00838                     suml = sdot_(&i__2, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00839                             c_dim1 + 1], &c__1);
00840                     i__2 = *n - l2;
00841 /* Computing MIN */
00842                     i__3 = l2 + 1;
00843 /* Computing MIN */
00844                     i__4 = l2 + 1;
00845                     sumr = sdot_(&i__2, &c__[k1 + min(i__3, *n)* c_dim1], ldc, 
00846                              &b[l1 + min(i__4, *n)* b_dim1], ldb);
00847                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00848 
00849                     i__2 = k1 - 1;
00850                     suml = sdot_(&i__2, &a[k2 * a_dim1 + 1], &c__1, &c__[l1 * 
00851                             c_dim1 + 1], &c__1);
00852                     i__2 = *n - l2;
00853 /* Computing MIN */
00854                     i__3 = l2 + 1;
00855 /* Computing MIN */
00856                     i__4 = l2 + 1;
00857                     sumr = sdot_(&i__2, &c__[k2 + min(i__3, *n)* c_dim1], ldc, 
00858                              &b[l1 + min(i__4, *n)* b_dim1], ldb);
00859                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
00860 
00861                     r__1 = -sgn * b[l1 + l1 * b_dim1];
00862                     slaln2_(&c_true, &c__2, &c__1, &smin, &c_b26, &a[k1 + k1 *
00863                              a_dim1], lda, &c_b26, &c_b26, vec, &c__2, &r__1, 
00864                             &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
00865                     if (ierr != 0) {
00866                         *info = 1;
00867                     }
00868 
00869                     if (scaloc != 1.f) {
00870                         i__2 = *n;
00871                         for (j = 1; j <= i__2; ++j) {
00872                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00873 /* L150: */
00874                         }
00875                         *scale *= scaloc;
00876                     }
00877                     c__[k1 + l1 * c_dim1] = x[0];
00878                     c__[k2 + l1 * c_dim1] = x[1];
00879 
00880                 } else if (l1 != l2 && k1 == k2) {
00881 
00882                     i__2 = k1 - 1;
00883                     suml = sdot_(&i__2, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00884                             c_dim1 + 1], &c__1);
00885                     i__2 = *n - l2;
00886 /* Computing MIN */
00887                     i__3 = l2 + 1;
00888 /* Computing MIN */
00889                     i__4 = l2 + 1;
00890                     sumr = sdot_(&i__2, &c__[k1 + min(i__3, *n)* c_dim1], ldc, 
00891                              &b[l1 + min(i__4, *n)* b_dim1], ldb);
00892                     vec[0] = sgn * (c__[k1 + l1 * c_dim1] - (suml + sgn * 
00893                             sumr));
00894 
00895                     i__2 = k1 - 1;
00896                     suml = sdot_(&i__2, &a[k1 * a_dim1 + 1], &c__1, &c__[l2 * 
00897                             c_dim1 + 1], &c__1);
00898                     i__2 = *n - l2;
00899 /* Computing MIN */
00900                     i__3 = l2 + 1;
00901 /* Computing MIN */
00902                     i__4 = l2 + 1;
00903                     sumr = sdot_(&i__2, &c__[k1 + min(i__3, *n)* c_dim1], ldc, 
00904                              &b[l2 + min(i__4, *n)* b_dim1], ldb);
00905                     vec[1] = sgn * (c__[k1 + l2 * c_dim1] - (suml + sgn * 
00906                             sumr));
00907 
00908                     r__1 = -sgn * a[k1 + k1 * a_dim1];
00909                     slaln2_(&c_false, &c__2, &c__1, &smin, &c_b26, &b[l1 + l1 
00910                             * b_dim1], ldb, &c_b26, &c_b26, vec, &c__2, &r__1, 
00911                              &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
00912                     if (ierr != 0) {
00913                         *info = 1;
00914                     }
00915 
00916                     if (scaloc != 1.f) {
00917                         i__2 = *n;
00918                         for (j = 1; j <= i__2; ++j) {
00919                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00920 /* L160: */
00921                         }
00922                         *scale *= scaloc;
00923                     }
00924                     c__[k1 + l1 * c_dim1] = x[0];
00925                     c__[k1 + l2 * c_dim1] = x[1];
00926 
00927                 } else if (l1 != l2 && k1 != k2) {
00928 
00929                     i__2 = k1 - 1;
00930                     suml = sdot_(&i__2, &a[k1 * a_dim1 + 1], &c__1, &c__[l1 * 
00931                             c_dim1 + 1], &c__1);
00932                     i__2 = *n - l2;
00933 /* Computing MIN */
00934                     i__3 = l2 + 1;
00935 /* Computing MIN */
00936                     i__4 = l2 + 1;
00937                     sumr = sdot_(&i__2, &c__[k1 + min(i__3, *n)* c_dim1], ldc, 
00938                              &b[l1 + min(i__4, *n)* b_dim1], ldb);
00939                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
00940 
00941                     i__2 = k1 - 1;
00942                     suml = sdot_(&i__2, &a[k1 * a_dim1 + 1], &c__1, &c__[l2 * 
00943                             c_dim1 + 1], &c__1);
00944                     i__2 = *n - l2;
00945 /* Computing MIN */
00946                     i__3 = l2 + 1;
00947 /* Computing MIN */
00948                     i__4 = l2 + 1;
00949                     sumr = sdot_(&i__2, &c__[k1 + min(i__3, *n)* c_dim1], ldc, 
00950                              &b[l2 + min(i__4, *n)* b_dim1], ldb);
00951                     vec[2] = c__[k1 + l2 * c_dim1] - (suml + sgn * sumr);
00952 
00953                     i__2 = k1 - 1;
00954                     suml = sdot_(&i__2, &a[k2 * a_dim1 + 1], &c__1, &c__[l1 * 
00955                             c_dim1 + 1], &c__1);
00956                     i__2 = *n - l2;
00957 /* Computing MIN */
00958                     i__3 = l2 + 1;
00959 /* Computing MIN */
00960                     i__4 = l2 + 1;
00961                     sumr = sdot_(&i__2, &c__[k2 + min(i__3, *n)* c_dim1], ldc, 
00962                              &b[l1 + min(i__4, *n)* b_dim1], ldb);
00963                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
00964 
00965                     i__2 = k1 - 1;
00966                     suml = sdot_(&i__2, &a[k2 * a_dim1 + 1], &c__1, &c__[l2 * 
00967                             c_dim1 + 1], &c__1);
00968                     i__2 = *n - l2;
00969 /* Computing MIN */
00970                     i__3 = l2 + 1;
00971 /* Computing MIN */
00972                     i__4 = l2 + 1;
00973                     sumr = sdot_(&i__2, &c__[k2 + min(i__3, *n)* c_dim1], ldc, 
00974                              &b[l2 + min(i__4, *n)* b_dim1], ldb);
00975                     vec[3] = c__[k2 + l2 * c_dim1] - (suml + sgn * sumr);
00976 
00977                     slasy2_(&c_true, &c_true, isgn, &c__2, &c__2, &a[k1 + k1 *
00978                              a_dim1], lda, &b[l1 + l1 * b_dim1], ldb, vec, &
00979                             c__2, &scaloc, x, &c__2, &xnorm, &ierr);
00980                     if (ierr != 0) {
00981                         *info = 1;
00982                     }
00983 
00984                     if (scaloc != 1.f) {
00985                         i__2 = *n;
00986                         for (j = 1; j <= i__2; ++j) {
00987                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
00988 /* L170: */
00989                         }
00990                         *scale *= scaloc;
00991                     }
00992                     c__[k1 + l1 * c_dim1] = x[0];
00993                     c__[k1 + l2 * c_dim1] = x[2];
00994                     c__[k2 + l1 * c_dim1] = x[1];
00995                     c__[k2 + l2 * c_dim1] = x[3];
00996                 }
00997 
00998 L180:
00999                 ;
01000             }
01001 L190:
01002             ;
01003         }
01004 
01005     } else if (notrna && ! notrnb) {
01006 
01007 /*        Solve    A*X + ISGN*X*B' = scale*C. */
01008 
01009 /*        The (K,L)th block of X is determined starting from */
01010 /*        bottom-right corner column by column by */
01011 
01012 /*            A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)' = C(K,L) - R(K,L) */
01013 
01014 /*        Where */
01015 /*                      M                          N */
01016 /*            R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)']. */
01017 /*                    I=K+1                      J=L+1 */
01018 
01019 /*        Start column loop (index = L) */
01020 /*        L1 (L2): column index of the first (last) row of X(K,L) */
01021 
01022         lnext = *n;
01023         for (l = *n; l >= 1; --l) {
01024             if (l > lnext) {
01025                 goto L250;
01026             }
01027             if (l == 1) {
01028                 l1 = l;
01029                 l2 = l;
01030             } else {
01031                 if (b[l + (l - 1) * b_dim1] != 0.f) {
01032                     l1 = l - 1;
01033                     l2 = l;
01034                     lnext = l - 2;
01035                 } else {
01036                     l1 = l;
01037                     l2 = l;
01038                     lnext = l - 1;
01039                 }
01040             }
01041 
01042 /*           Start row loop (index = K) */
01043 /*           K1 (K2): row index of the first (last) row of X(K,L) */
01044 
01045             knext = *m;
01046             for (k = *m; k >= 1; --k) {
01047                 if (k > knext) {
01048                     goto L240;
01049                 }
01050                 if (k == 1) {
01051                     k1 = k;
01052                     k2 = k;
01053                 } else {
01054                     if (a[k + (k - 1) * a_dim1] != 0.f) {
01055                         k1 = k - 1;
01056                         k2 = k;
01057                         knext = k - 2;
01058                     } else {
01059                         k1 = k;
01060                         k2 = k;
01061                         knext = k - 1;
01062                     }
01063                 }
01064 
01065                 if (l1 == l2 && k1 == k2) {
01066                     i__1 = *m - k1;
01067 /* Computing MIN */
01068                     i__2 = k1 + 1;
01069 /* Computing MIN */
01070                     i__3 = k1 + 1;
01071                     suml = sdot_(&i__1, &a[k1 + min(i__2, *m)* a_dim1], lda, &
01072                             c__[min(i__3, *m)+ l1 * c_dim1], &c__1);
01073                     i__1 = *n - l1;
01074 /* Computing MIN */
01075                     i__2 = l1 + 1;
01076 /* Computing MIN */
01077                     i__3 = l1 + 1;
01078                     sumr = sdot_(&i__1, &c__[k1 + min(i__2, *n)* c_dim1], ldc, 
01079                              &b[l1 + min(i__3, *n)* b_dim1], ldb);
01080                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
01081                     scaloc = 1.f;
01082 
01083                     a11 = a[k1 + k1 * a_dim1] + sgn * b[l1 + l1 * b_dim1];
01084                     da11 = dabs(a11);
01085                     if (da11 <= smin) {
01086                         a11 = smin;
01087                         da11 = smin;
01088                         *info = 1;
01089                     }
01090                     db = dabs(vec[0]);
01091                     if (da11 < 1.f && db > 1.f) {
01092                         if (db > bignum * da11) {
01093                             scaloc = 1.f / db;
01094                         }
01095                     }
01096                     x[0] = vec[0] * scaloc / a11;
01097 
01098                     if (scaloc != 1.f) {
01099                         i__1 = *n;
01100                         for (j = 1; j <= i__1; ++j) {
01101                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
01102 /* L200: */
01103                         }
01104                         *scale *= scaloc;
01105                     }
01106                     c__[k1 + l1 * c_dim1] = x[0];
01107 
01108                 } else if (l1 == l2 && k1 != k2) {
01109 
01110                     i__1 = *m - k2;
01111 /* Computing MIN */
01112                     i__2 = k2 + 1;
01113 /* Computing MIN */
01114                     i__3 = k2 + 1;
01115                     suml = sdot_(&i__1, &a[k1 + min(i__2, *m)* a_dim1], lda, &
01116                             c__[min(i__3, *m)+ l1 * c_dim1], &c__1);
01117                     i__1 = *n - l2;
01118 /* Computing MIN */
01119                     i__2 = l2 + 1;
01120 /* Computing MIN */
01121                     i__3 = l2 + 1;
01122                     sumr = sdot_(&i__1, &c__[k1 + min(i__2, *n)* c_dim1], ldc, 
01123                              &b[l1 + min(i__3, *n)* b_dim1], ldb);
01124                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
01125 
01126                     i__1 = *m - k2;
01127 /* Computing MIN */
01128                     i__2 = k2 + 1;
01129 /* Computing MIN */
01130                     i__3 = k2 + 1;
01131                     suml = sdot_(&i__1, &a[k2 + min(i__2, *m)* a_dim1], lda, &
01132                             c__[min(i__3, *m)+ l1 * c_dim1], &c__1);
01133                     i__1 = *n - l2;
01134 /* Computing MIN */
01135                     i__2 = l2 + 1;
01136 /* Computing MIN */
01137                     i__3 = l2 + 1;
01138                     sumr = sdot_(&i__1, &c__[k2 + min(i__2, *n)* c_dim1], ldc, 
01139                              &b[l1 + min(i__3, *n)* b_dim1], ldb);
01140                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
01141 
01142                     r__1 = -sgn * b[l1 + l1 * b_dim1];
01143                     slaln2_(&c_false, &c__2, &c__1, &smin, &c_b26, &a[k1 + k1 
01144                             * a_dim1], lda, &c_b26, &c_b26, vec, &c__2, &r__1, 
01145                              &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
01146                     if (ierr != 0) {
01147                         *info = 1;
01148                     }
01149 
01150                     if (scaloc != 1.f) {
01151                         i__1 = *n;
01152                         for (j = 1; j <= i__1; ++j) {
01153                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
01154 /* L210: */
01155                         }
01156                         *scale *= scaloc;
01157                     }
01158                     c__[k1 + l1 * c_dim1] = x[0];
01159                     c__[k2 + l1 * c_dim1] = x[1];
01160 
01161                 } else if (l1 != l2 && k1 == k2) {
01162 
01163                     i__1 = *m - k1;
01164 /* Computing MIN */
01165                     i__2 = k1 + 1;
01166 /* Computing MIN */
01167                     i__3 = k1 + 1;
01168                     suml = sdot_(&i__1, &a[k1 + min(i__2, *m)* a_dim1], lda, &
01169                             c__[min(i__3, *m)+ l1 * c_dim1], &c__1);
01170                     i__1 = *n - l2;
01171 /* Computing MIN */
01172                     i__2 = l2 + 1;
01173 /* Computing MIN */
01174                     i__3 = l2 + 1;
01175                     sumr = sdot_(&i__1, &c__[k1 + min(i__2, *n)* c_dim1], ldc, 
01176                              &b[l1 + min(i__3, *n)* b_dim1], ldb);
01177                     vec[0] = sgn * (c__[k1 + l1 * c_dim1] - (suml + sgn * 
01178                             sumr));
01179 
01180                     i__1 = *m - k1;
01181 /* Computing MIN */
01182                     i__2 = k1 + 1;
01183 /* Computing MIN */
01184                     i__3 = k1 + 1;
01185                     suml = sdot_(&i__1, &a[k1 + min(i__2, *m)* a_dim1], lda, &
01186                             c__[min(i__3, *m)+ l2 * c_dim1], &c__1);
01187                     i__1 = *n - l2;
01188 /* Computing MIN */
01189                     i__2 = l2 + 1;
01190 /* Computing MIN */
01191                     i__3 = l2 + 1;
01192                     sumr = sdot_(&i__1, &c__[k1 + min(i__2, *n)* c_dim1], ldc, 
01193                              &b[l2 + min(i__3, *n)* b_dim1], ldb);
01194                     vec[1] = sgn * (c__[k1 + l2 * c_dim1] - (suml + sgn * 
01195                             sumr));
01196 
01197                     r__1 = -sgn * a[k1 + k1 * a_dim1];
01198                     slaln2_(&c_false, &c__2, &c__1, &smin, &c_b26, &b[l1 + l1 
01199                             * b_dim1], ldb, &c_b26, &c_b26, vec, &c__2, &r__1, 
01200                              &c_b30, x, &c__2, &scaloc, &xnorm, &ierr);
01201                     if (ierr != 0) {
01202                         *info = 1;
01203                     }
01204 
01205                     if (scaloc != 1.f) {
01206                         i__1 = *n;
01207                         for (j = 1; j <= i__1; ++j) {
01208                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
01209 /* L220: */
01210                         }
01211                         *scale *= scaloc;
01212                     }
01213                     c__[k1 + l1 * c_dim1] = x[0];
01214                     c__[k1 + l2 * c_dim1] = x[1];
01215 
01216                 } else if (l1 != l2 && k1 != k2) {
01217 
01218                     i__1 = *m - k2;
01219 /* Computing MIN */
01220                     i__2 = k2 + 1;
01221 /* Computing MIN */
01222                     i__3 = k2 + 1;
01223                     suml = sdot_(&i__1, &a[k1 + min(i__2, *m)* a_dim1], lda, &
01224                             c__[min(i__3, *m)+ l1 * c_dim1], &c__1);
01225                     i__1 = *n - l2;
01226 /* Computing MIN */
01227                     i__2 = l2 + 1;
01228 /* Computing MIN */
01229                     i__3 = l2 + 1;
01230                     sumr = sdot_(&i__1, &c__[k1 + min(i__2, *n)* c_dim1], ldc, 
01231                              &b[l1 + min(i__3, *n)* b_dim1], ldb);
01232                     vec[0] = c__[k1 + l1 * c_dim1] - (suml + sgn * sumr);
01233 
01234                     i__1 = *m - k2;
01235 /* Computing MIN */
01236                     i__2 = k2 + 1;
01237 /* Computing MIN */
01238                     i__3 = k2 + 1;
01239                     suml = sdot_(&i__1, &a[k1 + min(i__2, *m)* a_dim1], lda, &
01240                             c__[min(i__3, *m)+ l2 * c_dim1], &c__1);
01241                     i__1 = *n - l2;
01242 /* Computing MIN */
01243                     i__2 = l2 + 1;
01244 /* Computing MIN */
01245                     i__3 = l2 + 1;
01246                     sumr = sdot_(&i__1, &c__[k1 + min(i__2, *n)* c_dim1], ldc, 
01247                              &b[l2 + min(i__3, *n)* b_dim1], ldb);
01248                     vec[2] = c__[k1 + l2 * c_dim1] - (suml + sgn * sumr);
01249 
01250                     i__1 = *m - k2;
01251 /* Computing MIN */
01252                     i__2 = k2 + 1;
01253 /* Computing MIN */
01254                     i__3 = k2 + 1;
01255                     suml = sdot_(&i__1, &a[k2 + min(i__2, *m)* a_dim1], lda, &
01256                             c__[min(i__3, *m)+ l1 * c_dim1], &c__1);
01257                     i__1 = *n - l2;
01258 /* Computing MIN */
01259                     i__2 = l2 + 1;
01260 /* Computing MIN */
01261                     i__3 = l2 + 1;
01262                     sumr = sdot_(&i__1, &c__[k2 + min(i__2, *n)* c_dim1], ldc, 
01263                              &b[l1 + min(i__3, *n)* b_dim1], ldb);
01264                     vec[1] = c__[k2 + l1 * c_dim1] - (suml + sgn * sumr);
01265 
01266                     i__1 = *m - k2;
01267 /* Computing MIN */
01268                     i__2 = k2 + 1;
01269 /* Computing MIN */
01270                     i__3 = k2 + 1;
01271                     suml = sdot_(&i__1, &a[k2 + min(i__2, *m)* a_dim1], lda, &
01272                             c__[min(i__3, *m)+ l2 * c_dim1], &c__1);
01273                     i__1 = *n - l2;
01274 /* Computing MIN */
01275                     i__2 = l2 + 1;
01276 /* Computing MIN */
01277                     i__3 = l2 + 1;
01278                     sumr = sdot_(&i__1, &c__[k2 + min(i__2, *n)* c_dim1], ldc, 
01279                              &b[l2 + min(i__3, *n)* b_dim1], ldb);
01280                     vec[3] = c__[k2 + l2 * c_dim1] - (suml + sgn * sumr);
01281 
01282                     slasy2_(&c_false, &c_true, isgn, &c__2, &c__2, &a[k1 + k1 
01283                             * a_dim1], lda, &b[l1 + l1 * b_dim1], ldb, vec, &
01284                             c__2, &scaloc, x, &c__2, &xnorm, &ierr);
01285                     if (ierr != 0) {
01286                         *info = 1;
01287                     }
01288 
01289                     if (scaloc != 1.f) {
01290                         i__1 = *n;
01291                         for (j = 1; j <= i__1; ++j) {
01292                             sscal_(m, &scaloc, &c__[j * c_dim1 + 1], &c__1);
01293 /* L230: */
01294                         }
01295                         *scale *= scaloc;
01296                     }
01297                     c__[k1 + l1 * c_dim1] = x[0];
01298                     c__[k1 + l2 * c_dim1] = x[2];
01299                     c__[k2 + l1 * c_dim1] = x[1];
01300                     c__[k2 + l2 * c_dim1] = x[3];
01301                 }
01302 
01303 L240:
01304                 ;
01305             }
01306 L250:
01307             ;
01308         }
01309 
01310     }
01311 
01312     return 0;
01313 
01314 /*     End of STRSYL */
01315 
01316 } /* strsyl_ */

swiftnav
Author(s):
autogenerated on Sat Jun 8 2019 18:56:15