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

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