cgtts2.c
Go to the documentation of this file.
00001 /* cgtts2.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 /* Subroutine */ int cgtts2_(integer *itrans, integer *n, integer *nrhs, 
00017         complex *dl, complex *d__, complex *du, complex *du2, integer *ipiv, 
00018         complex *b, integer *ldb)
00019 {
00020     /* System generated locals */
00021     integer b_dim1, b_offset, i__1, i__2, i__3, i__4, i__5, i__6, i__7, i__8;
00022     complex q__1, q__2, q__3, q__4, q__5, q__6, q__7, q__8;
00023 
00024     /* Builtin functions */
00025     void c_div(complex *, complex *, complex *), r_cnjg(complex *, complex *);
00026 
00027     /* Local variables */
00028     integer i__, j;
00029     complex temp;
00030 
00031 
00032 /*  -- LAPACK auxiliary routine (version 3.2) -- */
00033 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00034 /*     November 2006 */
00035 
00036 /*     .. Scalar Arguments .. */
00037 /*     .. */
00038 /*     .. Array Arguments .. */
00039 /*     .. */
00040 
00041 /*  Purpose */
00042 /*  ======= */
00043 
00044 /*  CGTTS2 solves one of the systems of equations */
00045 /*     A * X = B,  A**T * X = B,  or  A**H * X = B, */
00046 /*  with a tridiagonal matrix A using the LU factorization computed */
00047 /*  by CGTTRF. */
00048 
00049 /*  Arguments */
00050 /*  ========= */
00051 
00052 /*  ITRANS  (input) INTEGER */
00053 /*          Specifies the form of the system of equations. */
00054 /*          = 0:  A * X = B     (No transpose) */
00055 /*          = 1:  A**T * X = B  (Transpose) */
00056 /*          = 2:  A**H * X = B  (Conjugate transpose) */
00057 
00058 /*  N       (input) INTEGER */
00059 /*          The order of the matrix A. */
00060 
00061 /*  NRHS    (input) INTEGER */
00062 /*          The number of right hand sides, i.e., the number of columns */
00063 /*          of the matrix B.  NRHS >= 0. */
00064 
00065 /*  DL      (input) COMPLEX array, dimension (N-1) */
00066 /*          The (n-1) multipliers that define the matrix L from the */
00067 /*          LU factorization of A. */
00068 
00069 /*  D       (input) COMPLEX array, dimension (N) */
00070 /*          The n diagonal elements of the upper triangular matrix U from */
00071 /*          the LU factorization of A. */
00072 
00073 /*  DU      (input) COMPLEX array, dimension (N-1) */
00074 /*          The (n-1) elements of the first super-diagonal of U. */
00075 
00076 /*  DU2     (input) COMPLEX array, dimension (N-2) */
00077 /*          The (n-2) elements of the second super-diagonal of U. */
00078 
00079 /*  IPIV    (input) INTEGER array, dimension (N) */
00080 /*          The pivot indices; for 1 <= i <= n, row i of the matrix was */
00081 /*          interchanged with row IPIV(i).  IPIV(i) will always be either */
00082 /*          i or i+1; IPIV(i) = i indicates a row interchange was not */
00083 /*          required. */
00084 
00085 /*  B       (input/output) COMPLEX array, dimension (LDB,NRHS) */
00086 /*          On entry, the matrix of right hand side vectors B. */
00087 /*          On exit, B is overwritten by the solution vectors X. */
00088 
00089 /*  LDB     (input) INTEGER */
00090 /*          The leading dimension of the array B.  LDB >= max(1,N). */
00091 
00092 /*  ===================================================================== */
00093 
00094 /*     .. Local Scalars .. */
00095 /*     .. */
00096 /*     .. Intrinsic Functions .. */
00097 /*     .. */
00098 /*     .. Executable Statements .. */
00099 
00100 /*     Quick return if possible */
00101 
00102     /* Parameter adjustments */
00103     --dl;
00104     --d__;
00105     --du;
00106     --du2;
00107     --ipiv;
00108     b_dim1 = *ldb;
00109     b_offset = 1 + b_dim1;
00110     b -= b_offset;
00111 
00112     /* Function Body */
00113     if (*n == 0 || *nrhs == 0) {
00114         return 0;
00115     }
00116 
00117     if (*itrans == 0) {
00118 
00119 /*        Solve A*X = B using the LU factorization of A, */
00120 /*        overwriting each right hand side vector with its solution. */
00121 
00122         if (*nrhs <= 1) {
00123             j = 1;
00124 L10:
00125 
00126 /*           Solve L*x = b. */
00127 
00128             i__1 = *n - 1;
00129             for (i__ = 1; i__ <= i__1; ++i__) {
00130                 if (ipiv[i__] == i__) {
00131                     i__2 = i__ + 1 + j * b_dim1;
00132                     i__3 = i__ + 1 + j * b_dim1;
00133                     i__4 = i__;
00134                     i__5 = i__ + j * b_dim1;
00135                     q__2.r = dl[i__4].r * b[i__5].r - dl[i__4].i * b[i__5].i, 
00136                             q__2.i = dl[i__4].r * b[i__5].i + dl[i__4].i * b[
00137                             i__5].r;
00138                     q__1.r = b[i__3].r - q__2.r, q__1.i = b[i__3].i - q__2.i;
00139                     b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00140                 } else {
00141                     i__2 = i__ + j * b_dim1;
00142                     temp.r = b[i__2].r, temp.i = b[i__2].i;
00143                     i__2 = i__ + j * b_dim1;
00144                     i__3 = i__ + 1 + j * b_dim1;
00145                     b[i__2].r = b[i__3].r, b[i__2].i = b[i__3].i;
00146                     i__2 = i__ + 1 + j * b_dim1;
00147                     i__3 = i__;
00148                     i__4 = i__ + j * b_dim1;
00149                     q__2.r = dl[i__3].r * b[i__4].r - dl[i__3].i * b[i__4].i, 
00150                             q__2.i = dl[i__3].r * b[i__4].i + dl[i__3].i * b[
00151                             i__4].r;
00152                     q__1.r = temp.r - q__2.r, q__1.i = temp.i - q__2.i;
00153                     b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00154                 }
00155 /* L20: */
00156             }
00157 
00158 /*           Solve U*x = b. */
00159 
00160             i__1 = *n + j * b_dim1;
00161             c_div(&q__1, &b[*n + j * b_dim1], &d__[*n]);
00162             b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00163             if (*n > 1) {
00164                 i__1 = *n - 1 + j * b_dim1;
00165                 i__2 = *n - 1 + j * b_dim1;
00166                 i__3 = *n - 1;
00167                 i__4 = *n + j * b_dim1;
00168                 q__3.r = du[i__3].r * b[i__4].r - du[i__3].i * b[i__4].i, 
00169                         q__3.i = du[i__3].r * b[i__4].i + du[i__3].i * b[i__4]
00170                         .r;
00171                 q__2.r = b[i__2].r - q__3.r, q__2.i = b[i__2].i - q__3.i;
00172                 c_div(&q__1, &q__2, &d__[*n - 1]);
00173                 b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00174             }
00175             for (i__ = *n - 2; i__ >= 1; --i__) {
00176                 i__1 = i__ + j * b_dim1;
00177                 i__2 = i__ + j * b_dim1;
00178                 i__3 = i__;
00179                 i__4 = i__ + 1 + j * b_dim1;
00180                 q__4.r = du[i__3].r * b[i__4].r - du[i__3].i * b[i__4].i, 
00181                         q__4.i = du[i__3].r * b[i__4].i + du[i__3].i * b[i__4]
00182                         .r;
00183                 q__3.r = b[i__2].r - q__4.r, q__3.i = b[i__2].i - q__4.i;
00184                 i__5 = i__;
00185                 i__6 = i__ + 2 + j * b_dim1;
00186                 q__5.r = du2[i__5].r * b[i__6].r - du2[i__5].i * b[i__6].i, 
00187                         q__5.i = du2[i__5].r * b[i__6].i + du2[i__5].i * b[
00188                         i__6].r;
00189                 q__2.r = q__3.r - q__5.r, q__2.i = q__3.i - q__5.i;
00190                 c_div(&q__1, &q__2, &d__[i__]);
00191                 b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00192 /* L30: */
00193             }
00194             if (j < *nrhs) {
00195                 ++j;
00196                 goto L10;
00197             }
00198         } else {
00199             i__1 = *nrhs;
00200             for (j = 1; j <= i__1; ++j) {
00201 
00202 /*           Solve L*x = b. */
00203 
00204                 i__2 = *n - 1;
00205                 for (i__ = 1; i__ <= i__2; ++i__) {
00206                     if (ipiv[i__] == i__) {
00207                         i__3 = i__ + 1 + j * b_dim1;
00208                         i__4 = i__ + 1 + j * b_dim1;
00209                         i__5 = i__;
00210                         i__6 = i__ + j * b_dim1;
00211                         q__2.r = dl[i__5].r * b[i__6].r - dl[i__5].i * b[i__6]
00212                                 .i, q__2.i = dl[i__5].r * b[i__6].i + dl[i__5]
00213                                 .i * b[i__6].r;
00214                         q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4].i - 
00215                                 q__2.i;
00216                         b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00217                     } else {
00218                         i__3 = i__ + j * b_dim1;
00219                         temp.r = b[i__3].r, temp.i = b[i__3].i;
00220                         i__3 = i__ + j * b_dim1;
00221                         i__4 = i__ + 1 + j * b_dim1;
00222                         b[i__3].r = b[i__4].r, b[i__3].i = b[i__4].i;
00223                         i__3 = i__ + 1 + j * b_dim1;
00224                         i__4 = i__;
00225                         i__5 = i__ + j * b_dim1;
00226                         q__2.r = dl[i__4].r * b[i__5].r - dl[i__4].i * b[i__5]
00227                                 .i, q__2.i = dl[i__4].r * b[i__5].i + dl[i__4]
00228                                 .i * b[i__5].r;
00229                         q__1.r = temp.r - q__2.r, q__1.i = temp.i - q__2.i;
00230                         b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00231                     }
00232 /* L40: */
00233                 }
00234 
00235 /*           Solve U*x = b. */
00236 
00237                 i__2 = *n + j * b_dim1;
00238                 c_div(&q__1, &b[*n + j * b_dim1], &d__[*n]);
00239                 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00240                 if (*n > 1) {
00241                     i__2 = *n - 1 + j * b_dim1;
00242                     i__3 = *n - 1 + j * b_dim1;
00243                     i__4 = *n - 1;
00244                     i__5 = *n + j * b_dim1;
00245                     q__3.r = du[i__4].r * b[i__5].r - du[i__4].i * b[i__5].i, 
00246                             q__3.i = du[i__4].r * b[i__5].i + du[i__4].i * b[
00247                             i__5].r;
00248                     q__2.r = b[i__3].r - q__3.r, q__2.i = b[i__3].i - q__3.i;
00249                     c_div(&q__1, &q__2, &d__[*n - 1]);
00250                     b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00251                 }
00252                 for (i__ = *n - 2; i__ >= 1; --i__) {
00253                     i__2 = i__ + j * b_dim1;
00254                     i__3 = i__ + j * b_dim1;
00255                     i__4 = i__;
00256                     i__5 = i__ + 1 + j * b_dim1;
00257                     q__4.r = du[i__4].r * b[i__5].r - du[i__4].i * b[i__5].i, 
00258                             q__4.i = du[i__4].r * b[i__5].i + du[i__4].i * b[
00259                             i__5].r;
00260                     q__3.r = b[i__3].r - q__4.r, q__3.i = b[i__3].i - q__4.i;
00261                     i__6 = i__;
00262                     i__7 = i__ + 2 + j * b_dim1;
00263                     q__5.r = du2[i__6].r * b[i__7].r - du2[i__6].i * b[i__7]
00264                             .i, q__5.i = du2[i__6].r * b[i__7].i + du2[i__6]
00265                             .i * b[i__7].r;
00266                     q__2.r = q__3.r - q__5.r, q__2.i = q__3.i - q__5.i;
00267                     c_div(&q__1, &q__2, &d__[i__]);
00268                     b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00269 /* L50: */
00270                 }
00271 /* L60: */
00272             }
00273         }
00274     } else if (*itrans == 1) {
00275 
00276 /*        Solve A**T * X = B. */
00277 
00278         if (*nrhs <= 1) {
00279             j = 1;
00280 L70:
00281 
00282 /*           Solve U**T * x = b. */
00283 
00284             i__1 = j * b_dim1 + 1;
00285             c_div(&q__1, &b[j * b_dim1 + 1], &d__[1]);
00286             b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00287             if (*n > 1) {
00288                 i__1 = j * b_dim1 + 2;
00289                 i__2 = j * b_dim1 + 2;
00290                 i__3 = j * b_dim1 + 1;
00291                 q__3.r = du[1].r * b[i__3].r - du[1].i * b[i__3].i, q__3.i = 
00292                         du[1].r * b[i__3].i + du[1].i * b[i__3].r;
00293                 q__2.r = b[i__2].r - q__3.r, q__2.i = b[i__2].i - q__3.i;
00294                 c_div(&q__1, &q__2, &d__[2]);
00295                 b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00296             }
00297             i__1 = *n;
00298             for (i__ = 3; i__ <= i__1; ++i__) {
00299                 i__2 = i__ + j * b_dim1;
00300                 i__3 = i__ + j * b_dim1;
00301                 i__4 = i__ - 1;
00302                 i__5 = i__ - 1 + j * b_dim1;
00303                 q__4.r = du[i__4].r * b[i__5].r - du[i__4].i * b[i__5].i, 
00304                         q__4.i = du[i__4].r * b[i__5].i + du[i__4].i * b[i__5]
00305                         .r;
00306                 q__3.r = b[i__3].r - q__4.r, q__3.i = b[i__3].i - q__4.i;
00307                 i__6 = i__ - 2;
00308                 i__7 = i__ - 2 + j * b_dim1;
00309                 q__5.r = du2[i__6].r * b[i__7].r - du2[i__6].i * b[i__7].i, 
00310                         q__5.i = du2[i__6].r * b[i__7].i + du2[i__6].i * b[
00311                         i__7].r;
00312                 q__2.r = q__3.r - q__5.r, q__2.i = q__3.i - q__5.i;
00313                 c_div(&q__1, &q__2, &d__[i__]);
00314                 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00315 /* L80: */
00316             }
00317 
00318 /*           Solve L**T * x = b. */
00319 
00320             for (i__ = *n - 1; i__ >= 1; --i__) {
00321                 if (ipiv[i__] == i__) {
00322                     i__1 = i__ + j * b_dim1;
00323                     i__2 = i__ + j * b_dim1;
00324                     i__3 = i__;
00325                     i__4 = i__ + 1 + j * b_dim1;
00326                     q__2.r = dl[i__3].r * b[i__4].r - dl[i__3].i * b[i__4].i, 
00327                             q__2.i = dl[i__3].r * b[i__4].i + dl[i__3].i * b[
00328                             i__4].r;
00329                     q__1.r = b[i__2].r - q__2.r, q__1.i = b[i__2].i - q__2.i;
00330                     b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00331                 } else {
00332                     i__1 = i__ + 1 + j * b_dim1;
00333                     temp.r = b[i__1].r, temp.i = b[i__1].i;
00334                     i__1 = i__ + 1 + j * b_dim1;
00335                     i__2 = i__ + j * b_dim1;
00336                     i__3 = i__;
00337                     q__2.r = dl[i__3].r * temp.r - dl[i__3].i * temp.i, 
00338                             q__2.i = dl[i__3].r * temp.i + dl[i__3].i * 
00339                             temp.r;
00340                     q__1.r = b[i__2].r - q__2.r, q__1.i = b[i__2].i - q__2.i;
00341                     b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00342                     i__1 = i__ + j * b_dim1;
00343                     b[i__1].r = temp.r, b[i__1].i = temp.i;
00344                 }
00345 /* L90: */
00346             }
00347             if (j < *nrhs) {
00348                 ++j;
00349                 goto L70;
00350             }
00351         } else {
00352             i__1 = *nrhs;
00353             for (j = 1; j <= i__1; ++j) {
00354 
00355 /*           Solve U**T * x = b. */
00356 
00357                 i__2 = j * b_dim1 + 1;
00358                 c_div(&q__1, &b[j * b_dim1 + 1], &d__[1]);
00359                 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00360                 if (*n > 1) {
00361                     i__2 = j * b_dim1 + 2;
00362                     i__3 = j * b_dim1 + 2;
00363                     i__4 = j * b_dim1 + 1;
00364                     q__3.r = du[1].r * b[i__4].r - du[1].i * b[i__4].i, 
00365                             q__3.i = du[1].r * b[i__4].i + du[1].i * b[i__4]
00366                             .r;
00367                     q__2.r = b[i__3].r - q__3.r, q__2.i = b[i__3].i - q__3.i;
00368                     c_div(&q__1, &q__2, &d__[2]);
00369                     b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00370                 }
00371                 i__2 = *n;
00372                 for (i__ = 3; i__ <= i__2; ++i__) {
00373                     i__3 = i__ + j * b_dim1;
00374                     i__4 = i__ + j * b_dim1;
00375                     i__5 = i__ - 1;
00376                     i__6 = i__ - 1 + j * b_dim1;
00377                     q__4.r = du[i__5].r * b[i__6].r - du[i__5].i * b[i__6].i, 
00378                             q__4.i = du[i__5].r * b[i__6].i + du[i__5].i * b[
00379                             i__6].r;
00380                     q__3.r = b[i__4].r - q__4.r, q__3.i = b[i__4].i - q__4.i;
00381                     i__7 = i__ - 2;
00382                     i__8 = i__ - 2 + j * b_dim1;
00383                     q__5.r = du2[i__7].r * b[i__8].r - du2[i__7].i * b[i__8]
00384                             .i, q__5.i = du2[i__7].r * b[i__8].i + du2[i__7]
00385                             .i * b[i__8].r;
00386                     q__2.r = q__3.r - q__5.r, q__2.i = q__3.i - q__5.i;
00387                     c_div(&q__1, &q__2, &d__[i__]);
00388                     b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00389 /* L100: */
00390                 }
00391 
00392 /*           Solve L**T * x = b. */
00393 
00394                 for (i__ = *n - 1; i__ >= 1; --i__) {
00395                     if (ipiv[i__] == i__) {
00396                         i__2 = i__ + j * b_dim1;
00397                         i__3 = i__ + j * b_dim1;
00398                         i__4 = i__;
00399                         i__5 = i__ + 1 + j * b_dim1;
00400                         q__2.r = dl[i__4].r * b[i__5].r - dl[i__4].i * b[i__5]
00401                                 .i, q__2.i = dl[i__4].r * b[i__5].i + dl[i__4]
00402                                 .i * b[i__5].r;
00403                         q__1.r = b[i__3].r - q__2.r, q__1.i = b[i__3].i - 
00404                                 q__2.i;
00405                         b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00406                     } else {
00407                         i__2 = i__ + 1 + j * b_dim1;
00408                         temp.r = b[i__2].r, temp.i = b[i__2].i;
00409                         i__2 = i__ + 1 + j * b_dim1;
00410                         i__3 = i__ + j * b_dim1;
00411                         i__4 = i__;
00412                         q__2.r = dl[i__4].r * temp.r - dl[i__4].i * temp.i, 
00413                                 q__2.i = dl[i__4].r * temp.i + dl[i__4].i * 
00414                                 temp.r;
00415                         q__1.r = b[i__3].r - q__2.r, q__1.i = b[i__3].i - 
00416                                 q__2.i;
00417                         b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00418                         i__2 = i__ + j * b_dim1;
00419                         b[i__2].r = temp.r, b[i__2].i = temp.i;
00420                     }
00421 /* L110: */
00422                 }
00423 /* L120: */
00424             }
00425         }
00426     } else {
00427 
00428 /*        Solve A**H * X = B. */
00429 
00430         if (*nrhs <= 1) {
00431             j = 1;
00432 L130:
00433 
00434 /*           Solve U**H * x = b. */
00435 
00436             i__1 = j * b_dim1 + 1;
00437             r_cnjg(&q__2, &d__[1]);
00438             c_div(&q__1, &b[j * b_dim1 + 1], &q__2);
00439             b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00440             if (*n > 1) {
00441                 i__1 = j * b_dim1 + 2;
00442                 i__2 = j * b_dim1 + 2;
00443                 r_cnjg(&q__4, &du[1]);
00444                 i__3 = j * b_dim1 + 1;
00445                 q__3.r = q__4.r * b[i__3].r - q__4.i * b[i__3].i, q__3.i = 
00446                         q__4.r * b[i__3].i + q__4.i * b[i__3].r;
00447                 q__2.r = b[i__2].r - q__3.r, q__2.i = b[i__2].i - q__3.i;
00448                 r_cnjg(&q__5, &d__[2]);
00449                 c_div(&q__1, &q__2, &q__5);
00450                 b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00451             }
00452             i__1 = *n;
00453             for (i__ = 3; i__ <= i__1; ++i__) {
00454                 i__2 = i__ + j * b_dim1;
00455                 i__3 = i__ + j * b_dim1;
00456                 r_cnjg(&q__5, &du[i__ - 1]);
00457                 i__4 = i__ - 1 + j * b_dim1;
00458                 q__4.r = q__5.r * b[i__4].r - q__5.i * b[i__4].i, q__4.i = 
00459                         q__5.r * b[i__4].i + q__5.i * b[i__4].r;
00460                 q__3.r = b[i__3].r - q__4.r, q__3.i = b[i__3].i - q__4.i;
00461                 r_cnjg(&q__7, &du2[i__ - 2]);
00462                 i__5 = i__ - 2 + j * b_dim1;
00463                 q__6.r = q__7.r * b[i__5].r - q__7.i * b[i__5].i, q__6.i = 
00464                         q__7.r * b[i__5].i + q__7.i * b[i__5].r;
00465                 q__2.r = q__3.r - q__6.r, q__2.i = q__3.i - q__6.i;
00466                 r_cnjg(&q__8, &d__[i__]);
00467                 c_div(&q__1, &q__2, &q__8);
00468                 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00469 /* L140: */
00470             }
00471 
00472 /*           Solve L**H * x = b. */
00473 
00474             for (i__ = *n - 1; i__ >= 1; --i__) {
00475                 if (ipiv[i__] == i__) {
00476                     i__1 = i__ + j * b_dim1;
00477                     i__2 = i__ + j * b_dim1;
00478                     r_cnjg(&q__3, &dl[i__]);
00479                     i__3 = i__ + 1 + j * b_dim1;
00480                     q__2.r = q__3.r * b[i__3].r - q__3.i * b[i__3].i, q__2.i =
00481                              q__3.r * b[i__3].i + q__3.i * b[i__3].r;
00482                     q__1.r = b[i__2].r - q__2.r, q__1.i = b[i__2].i - q__2.i;
00483                     b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00484                 } else {
00485                     i__1 = i__ + 1 + j * b_dim1;
00486                     temp.r = b[i__1].r, temp.i = b[i__1].i;
00487                     i__1 = i__ + 1 + j * b_dim1;
00488                     i__2 = i__ + j * b_dim1;
00489                     r_cnjg(&q__3, &dl[i__]);
00490                     q__2.r = q__3.r * temp.r - q__3.i * temp.i, q__2.i = 
00491                             q__3.r * temp.i + q__3.i * temp.r;
00492                     q__1.r = b[i__2].r - q__2.r, q__1.i = b[i__2].i - q__2.i;
00493                     b[i__1].r = q__1.r, b[i__1].i = q__1.i;
00494                     i__1 = i__ + j * b_dim1;
00495                     b[i__1].r = temp.r, b[i__1].i = temp.i;
00496                 }
00497 /* L150: */
00498             }
00499             if (j < *nrhs) {
00500                 ++j;
00501                 goto L130;
00502             }
00503         } else {
00504             i__1 = *nrhs;
00505             for (j = 1; j <= i__1; ++j) {
00506 
00507 /*           Solve U**H * x = b. */
00508 
00509                 i__2 = j * b_dim1 + 1;
00510                 r_cnjg(&q__2, &d__[1]);
00511                 c_div(&q__1, &b[j * b_dim1 + 1], &q__2);
00512                 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00513                 if (*n > 1) {
00514                     i__2 = j * b_dim1 + 2;
00515                     i__3 = j * b_dim1 + 2;
00516                     r_cnjg(&q__4, &du[1]);
00517                     i__4 = j * b_dim1 + 1;
00518                     q__3.r = q__4.r * b[i__4].r - q__4.i * b[i__4].i, q__3.i =
00519                              q__4.r * b[i__4].i + q__4.i * b[i__4].r;
00520                     q__2.r = b[i__3].r - q__3.r, q__2.i = b[i__3].i - q__3.i;
00521                     r_cnjg(&q__5, &d__[2]);
00522                     c_div(&q__1, &q__2, &q__5);
00523                     b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00524                 }
00525                 i__2 = *n;
00526                 for (i__ = 3; i__ <= i__2; ++i__) {
00527                     i__3 = i__ + j * b_dim1;
00528                     i__4 = i__ + j * b_dim1;
00529                     r_cnjg(&q__5, &du[i__ - 1]);
00530                     i__5 = i__ - 1 + j * b_dim1;
00531                     q__4.r = q__5.r * b[i__5].r - q__5.i * b[i__5].i, q__4.i =
00532                              q__5.r * b[i__5].i + q__5.i * b[i__5].r;
00533                     q__3.r = b[i__4].r - q__4.r, q__3.i = b[i__4].i - q__4.i;
00534                     r_cnjg(&q__7, &du2[i__ - 2]);
00535                     i__6 = i__ - 2 + j * b_dim1;
00536                     q__6.r = q__7.r * b[i__6].r - q__7.i * b[i__6].i, q__6.i =
00537                              q__7.r * b[i__6].i + q__7.i * b[i__6].r;
00538                     q__2.r = q__3.r - q__6.r, q__2.i = q__3.i - q__6.i;
00539                     r_cnjg(&q__8, &d__[i__]);
00540                     c_div(&q__1, &q__2, &q__8);
00541                     b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00542 /* L160: */
00543                 }
00544 
00545 /*           Solve L**H * x = b. */
00546 
00547                 for (i__ = *n - 1; i__ >= 1; --i__) {
00548                     if (ipiv[i__] == i__) {
00549                         i__2 = i__ + j * b_dim1;
00550                         i__3 = i__ + j * b_dim1;
00551                         r_cnjg(&q__3, &dl[i__]);
00552                         i__4 = i__ + 1 + j * b_dim1;
00553                         q__2.r = q__3.r * b[i__4].r - q__3.i * b[i__4].i, 
00554                                 q__2.i = q__3.r * b[i__4].i + q__3.i * b[i__4]
00555                                 .r;
00556                         q__1.r = b[i__3].r - q__2.r, q__1.i = b[i__3].i - 
00557                                 q__2.i;
00558                         b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00559                     } else {
00560                         i__2 = i__ + 1 + j * b_dim1;
00561                         temp.r = b[i__2].r, temp.i = b[i__2].i;
00562                         i__2 = i__ + 1 + j * b_dim1;
00563                         i__3 = i__ + j * b_dim1;
00564                         r_cnjg(&q__3, &dl[i__]);
00565                         q__2.r = q__3.r * temp.r - q__3.i * temp.i, q__2.i = 
00566                                 q__3.r * temp.i + q__3.i * temp.r;
00567                         q__1.r = b[i__3].r - q__2.r, q__1.i = b[i__3].i - 
00568                                 q__2.i;
00569                         b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00570                         i__2 = i__ + j * b_dim1;
00571                         b[i__2].r = temp.r, b[i__2].i = temp.i;
00572                     }
00573 /* L170: */
00574                 }
00575 /* L180: */
00576             }
00577         }
00578     }
00579 
00580 /*     End of CGTTS2 */
00581 
00582     return 0;
00583 } /* cgtts2_ */

swiftnav
Author(s):
autogenerated on Sat Jun 8 2019 18:55:28