dtfsm.c
Go to the documentation of this file.
00001 /* dtfsm.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 doublereal c_b23 = -1.;
00019 static doublereal c_b27 = 1.;
00020 
00021 /* Subroutine */ int dtfsm_(char *transr, char *side, char *uplo, char *trans, 
00022          char *diag, integer *m, integer *n, doublereal *alpha, doublereal *a, 
00023          doublereal *b, integer *ldb)
00024 {
00025     /* System generated locals */
00026     integer b_dim1, b_offset, i__1, i__2;
00027 
00028     /* Local variables */
00029     integer i__, j, k, m1, m2, n1, n2, info;
00030     logical normaltransr;
00031     extern /* Subroutine */ int dgemm_(char *, char *, integer *, integer *, 
00032             integer *, doublereal *, doublereal *, integer *, doublereal *, 
00033             integer *, doublereal *, doublereal *, integer *);
00034     logical lside;
00035     extern logical lsame_(char *, char *);
00036     logical lower;
00037     extern /* Subroutine */ int dtrsm_(char *, char *, char *, char *, 
00038             integer *, integer *, doublereal *, doublereal *, integer *, 
00039             doublereal *, integer *), xerbla_(
00040             char *, integer *);
00041     logical misodd, nisodd, notrans;
00042 
00043 
00044 /*  -- LAPACK routine (version 3.2.1)                                    -- */
00045 
00046 /*  -- Contributed by Fred Gustavson of the IBM Watson Research Center -- */
00047 /*  -- April 2009                                                      -- */
00048 
00049 /*  -- LAPACK is a software package provided by Univ. of Tennessee,    -- */
00050 /*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- */
00051 
00052 /*     .. */
00053 /*     .. Scalar Arguments .. */
00054 /*     .. */
00055 /*     .. Array Arguments .. */
00056 /*     .. */
00057 
00058 /*  Purpose */
00059 /*  ======= */
00060 
00061 /*  Level 3 BLAS like routine for A in RFP Format. */
00062 
00063 /*  DTFSM  solves the matrix equation */
00064 
00065 /*     op( A )*X = alpha*B  or  X*op( A ) = alpha*B */
00066 
00067 /*  where alpha is a scalar, X and B are m by n matrices, A is a unit, or */
00068 /*  non-unit,  upper or lower triangular matrix  and  op( A )  is one  of */
00069 
00070 /*     op( A ) = A   or   op( A ) = A'. */
00071 
00072 /*  A is in Rectangular Full Packed (RFP) Format. */
00073 
00074 /*  The matrix X is overwritten on B. */
00075 
00076 /*  Arguments */
00077 /*  ========== */
00078 
00079 /*  TRANSR - (input) CHARACTER */
00080 /*          = 'N':  The Normal Form of RFP A is stored; */
00081 /*          = 'T':  The Transpose Form of RFP A is stored. */
00082 
00083 /*  SIDE   - (input) CHARACTER */
00084 /*           On entry, SIDE specifies whether op( A ) appears on the left */
00085 /*           or right of X as follows: */
00086 
00087 /*              SIDE = 'L' or 'l'   op( A )*X = alpha*B. */
00088 
00089 /*              SIDE = 'R' or 'r'   X*op( A ) = alpha*B. */
00090 
00091 /*           Unchanged on exit. */
00092 
00093 /*  UPLO   - (input) CHARACTER */
00094 /*           On entry, UPLO specifies whether the RFP matrix A came from */
00095 /*           an upper or lower triangular matrix as follows: */
00096 /*           UPLO = 'U' or 'u' RFP A came from an upper triangular matrix */
00097 /*           UPLO = 'L' or 'l' RFP A came from a  lower triangular matrix */
00098 
00099 /*           Unchanged on exit. */
00100 
00101 /*  TRANS  - (input) CHARACTER */
00102 /*           On entry, TRANS  specifies the form of op( A ) to be used */
00103 /*           in the matrix multiplication as follows: */
00104 
00105 /*              TRANS  = 'N' or 'n'   op( A ) = A. */
00106 
00107 /*              TRANS  = 'T' or 't'   op( A ) = A'. */
00108 
00109 /*           Unchanged on exit. */
00110 
00111 /*  DIAG   - (input) CHARACTER */
00112 /*           On entry, DIAG specifies whether or not RFP A is unit */
00113 /*           triangular as follows: */
00114 
00115 /*              DIAG = 'U' or 'u'   A is assumed to be unit triangular. */
00116 
00117 /*              DIAG = 'N' or 'n'   A is not assumed to be unit */
00118 /*                                  triangular. */
00119 
00120 /*           Unchanged on exit. */
00121 
00122 /*  M      - (input) INTEGER. */
00123 /*           On entry, M specifies the number of rows of B. M must be at */
00124 /*           least zero. */
00125 /*           Unchanged on exit. */
00126 
00127 /*  N      - (input) INTEGER. */
00128 /*           On entry, N specifies the number of columns of B.  N must be */
00129 /*           at least zero. */
00130 /*           Unchanged on exit. */
00131 
00132 /*  ALPHA  - (input) DOUBLE PRECISION. */
00133 /*           On entry,  ALPHA specifies the scalar  alpha. When  alpha is */
00134 /*           zero then  A is not referenced and  B need not be set before */
00135 /*           entry. */
00136 /*           Unchanged on exit. */
00137 
00138 /*  A      - (input) DOUBLE PRECISION array, dimension (NT); */
00139 /*           NT = N*(N+1)/2. On entry, the matrix A in RFP Format. */
00140 /*           RFP Format is described by TRANSR, UPLO and N as follows: */
00141 /*           If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even; */
00142 /*           K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If */
00143 /*           TRANSR = 'T' then RFP is the transpose of RFP A as */
00144 /*           defined when TRANSR = 'N'. The contents of RFP A are defined */
00145 /*           by UPLO as follows: If UPLO = 'U' the RFP A contains the NT */
00146 /*           elements of upper packed A either in normal or */
00147 /*           transpose Format. If UPLO = 'L' the RFP A contains */
00148 /*           the NT elements of lower packed A either in normal or */
00149 /*           transpose Format. The LDA of RFP A is (N+1)/2 when */
00150 /*           TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is */
00151 /*           even and is N when is odd. */
00152 /*           See the Note below for more details. Unchanged on exit. */
00153 
00154 /*  B      - (input/ouptut) DOUBLE PRECISION array,  DIMENSION (LDB,N) */
00155 /*           Before entry,  the leading  m by n part of the array  B must */
00156 /*           contain  the  right-hand  side  matrix  B,  and  on exit  is */
00157 /*           overwritten by the solution matrix  X. */
00158 
00159 /*  LDB    - (input) INTEGER. */
00160 /*           On entry, LDB specifies the first dimension of B as declared */
00161 /*           in  the  calling  (sub)  program.   LDB  must  be  at  least */
00162 /*           max( 1, m ). */
00163 /*           Unchanged on exit. */
00164 
00165 /*  Further Details */
00166 /*  =============== */
00167 
00168 /*  We first consider Rectangular Full Packed (RFP) Format when N is */
00169 /*  even. We give an example where N = 6. */
00170 
00171 /*      AP is Upper             AP is Lower */
00172 
00173 /*   00 01 02 03 04 05       00 */
00174 /*      11 12 13 14 15       10 11 */
00175 /*         22 23 24 25       20 21 22 */
00176 /*            33 34 35       30 31 32 33 */
00177 /*               44 45       40 41 42 43 44 */
00178 /*                  55       50 51 52 53 54 55 */
00179 
00180 
00181 /*  Let TRANSR = 'N'. RFP holds AP as follows: */
00182 /*  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last */
00183 /*  three columns of AP upper. The lower triangle A(4:6,0:2) consists of */
00184 /*  the transpose of the first three columns of AP upper. */
00185 /*  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first */
00186 /*  three columns of AP lower. The upper triangle A(0:2,0:2) consists of */
00187 /*  the transpose of the last three columns of AP lower. */
00188 /*  This covers the case N even and TRANSR = 'N'. */
00189 
00190 /*         RFP A                   RFP A */
00191 
00192 /*        03 04 05                33 43 53 */
00193 /*        13 14 15                00 44 54 */
00194 /*        23 24 25                10 11 55 */
00195 /*        33 34 35                20 21 22 */
00196 /*        00 44 45                30 31 32 */
00197 /*        01 11 55                40 41 42 */
00198 /*        02 12 22                50 51 52 */
00199 
00200 /*  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the */
00201 /*  transpose of RFP A above. One therefore gets: */
00202 
00203 
00204 /*           RFP A                   RFP A */
00205 
00206 /*     03 13 23 33 00 01 02    33 00 10 20 30 40 50 */
00207 /*     04 14 24 34 44 11 12    43 44 11 21 31 41 51 */
00208 /*     05 15 25 35 45 55 22    53 54 55 22 32 42 52 */
00209 
00210 
00211 /*  We first consider Rectangular Full Packed (RFP) Format when N is */
00212 /*  odd. We give an example where N = 5. */
00213 
00214 /*     AP is Upper                 AP is Lower */
00215 
00216 /*   00 01 02 03 04              00 */
00217 /*      11 12 13 14              10 11 */
00218 /*         22 23 24              20 21 22 */
00219 /*            33 34              30 31 32 33 */
00220 /*               44              40 41 42 43 44 */
00221 
00222 
00223 /*  Let TRANSR = 'N'. RFP holds AP as follows: */
00224 /*  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last */
00225 /*  three columns of AP upper. The lower triangle A(3:4,0:1) consists of */
00226 /*  the transpose of the first two columns of AP upper. */
00227 /*  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first */
00228 /*  three columns of AP lower. The upper triangle A(0:1,1:2) consists of */
00229 /*  the transpose of the last two columns of AP lower. */
00230 /*  This covers the case N odd and TRANSR = 'N'. */
00231 
00232 /*         RFP A                   RFP A */
00233 
00234 /*        02 03 04                00 33 43 */
00235 /*        12 13 14                10 11 44 */
00236 /*        22 23 24                20 21 22 */
00237 /*        00 33 34                30 31 32 */
00238 /*        01 11 44                40 41 42 */
00239 
00240 /*  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the */
00241 /*  transpose of RFP A above. One therefore gets: */
00242 
00243 /*           RFP A                   RFP A */
00244 
00245 /*     02 12 22 00 01             00 10 20 30 40 50 */
00246 /*     03 13 23 33 11             33 11 21 31 41 51 */
00247 /*     04 14 24 34 44             43 44 22 32 42 52 */
00248 
00249 /*  Reference */
00250 /*  ========= */
00251 
00252 /*  ===================================================================== */
00253 
00254 /*     .. */
00255 /*     .. Parameters .. */
00256 /*     .. */
00257 /*     .. Local Scalars .. */
00258 /*     .. */
00259 /*     .. External Functions .. */
00260 /*     .. */
00261 /*     .. External Subroutines .. */
00262 /*     .. */
00263 /*     .. Intrinsic Functions .. */
00264 /*     .. */
00265 /*     .. Executable Statements .. */
00266 
00267 /*     Test the input parameters. */
00268 
00269     /* Parameter adjustments */
00270     b_dim1 = *ldb - 1 - 0 + 1;
00271     b_offset = 0 + b_dim1 * 0;
00272     b -= b_offset;
00273 
00274     /* Function Body */
00275     info = 0;
00276     normaltransr = lsame_(transr, "N");
00277     lside = lsame_(side, "L");
00278     lower = lsame_(uplo, "L");
00279     notrans = lsame_(trans, "N");
00280     if (! normaltransr && ! lsame_(transr, "T")) {
00281         info = -1;
00282     } else if (! lside && ! lsame_(side, "R")) {
00283         info = -2;
00284     } else if (! lower && ! lsame_(uplo, "U")) {
00285         info = -3;
00286     } else if (! notrans && ! lsame_(trans, "T")) {
00287         info = -4;
00288     } else if (! lsame_(diag, "N") && ! lsame_(diag, 
00289             "U")) {
00290         info = -5;
00291     } else if (*m < 0) {
00292         info = -6;
00293     } else if (*n < 0) {
00294         info = -7;
00295     } else if (*ldb < max(1,*m)) {
00296         info = -11;
00297     }
00298     if (info != 0) {
00299         i__1 = -info;
00300         xerbla_("DTFSM ", &i__1);
00301         return 0;
00302     }
00303 
00304 /*     Quick return when ( (N.EQ.0).OR.(M.EQ.0) ) */
00305 
00306     if (*m == 0 || *n == 0) {
00307         return 0;
00308     }
00309 
00310 /*     Quick return when ALPHA.EQ.(0D+0) */
00311 
00312     if (*alpha == 0.) {
00313         i__1 = *n - 1;
00314         for (j = 0; j <= i__1; ++j) {
00315             i__2 = *m - 1;
00316             for (i__ = 0; i__ <= i__2; ++i__) {
00317                 b[i__ + j * b_dim1] = 0.;
00318 /* L10: */
00319             }
00320 /* L20: */
00321         }
00322         return 0;
00323     }
00324 
00325     if (lside) {
00326 
00327 /*        SIDE = 'L' */
00328 
00329 /*        A is M-by-M. */
00330 /*        If M is odd, set NISODD = .TRUE., and M1 and M2. */
00331 /*        If M is even, NISODD = .FALSE., and M. */
00332 
00333         if (*m % 2 == 0) {
00334             misodd = FALSE_;
00335             k = *m / 2;
00336         } else {
00337             misodd = TRUE_;
00338             if (lower) {
00339                 m2 = *m / 2;
00340                 m1 = *m - m2;
00341             } else {
00342                 m1 = *m / 2;
00343                 m2 = *m - m1;
00344             }
00345         }
00346 
00347 
00348         if (misodd) {
00349 
00350 /*           SIDE = 'L' and N is odd */
00351 
00352             if (normaltransr) {
00353 
00354 /*              SIDE = 'L', N is odd, and TRANSR = 'N' */
00355 
00356                 if (lower) {
00357 
00358 /*                 SIDE  ='L', N is odd, TRANSR = 'N', and UPLO = 'L' */
00359 
00360                     if (notrans) {
00361 
00362 /*                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'L', and */
00363 /*                    TRANS = 'N' */
00364 
00365                         if (*m == 1) {
00366                             dtrsm_("L", "L", "N", diag, &m1, n, alpha, a, m, &
00367                                     b[b_offset], ldb);
00368                         } else {
00369                             dtrsm_("L", "L", "N", diag, &m1, n, alpha, a, m, &
00370                                     b[b_offset], ldb);
00371                             dgemm_("N", "N", &m2, n, &m1, &c_b23, &a[m1], m, &
00372                                     b[b_offset], ldb, alpha, &b[m1], ldb);
00373                             dtrsm_("L", "U", "T", diag, &m2, n, &c_b27, &a[*m]
00374 , m, &b[m1], ldb);
00375                         }
00376 
00377                     } else {
00378 
00379 /*                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'L', and */
00380 /*                    TRANS = 'T' */
00381 
00382                         if (*m == 1) {
00383                             dtrsm_("L", "L", "T", diag, &m1, n, alpha, a, m, &
00384                                     b[b_offset], ldb);
00385                         } else {
00386                             dtrsm_("L", "U", "N", diag, &m2, n, alpha, &a[*m], 
00387                                      m, &b[m1], ldb);
00388                             dgemm_("T", "N", &m1, n, &m2, &c_b23, &a[m1], m, &
00389                                     b[m1], ldb, alpha, &b[b_offset], ldb);
00390                             dtrsm_("L", "L", "T", diag, &m1, n, &c_b27, a, m, 
00391                                     &b[b_offset], ldb);
00392                         }
00393 
00394                     }
00395 
00396                 } else {
00397 
00398 /*                 SIDE  ='L', N is odd, TRANSR = 'N', and UPLO = 'U' */
00399 
00400                     if (! notrans) {
00401 
00402 /*                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'U', and */
00403 /*                    TRANS = 'N' */
00404 
00405                         dtrsm_("L", "L", "N", diag, &m1, n, alpha, &a[m2], m, 
00406                                 &b[b_offset], ldb);
00407                         dgemm_("T", "N", &m2, n, &m1, &c_b23, a, m, &b[
00408                                 b_offset], ldb, alpha, &b[m1], ldb);
00409                         dtrsm_("L", "U", "T", diag, &m2, n, &c_b27, &a[m1], m, 
00410                                  &b[m1], ldb);
00411 
00412                     } else {
00413 
00414 /*                    SIDE  ='L', N is odd, TRANSR = 'N', UPLO = 'U', and */
00415 /*                    TRANS = 'T' */
00416 
00417                         dtrsm_("L", "U", "N", diag, &m2, n, alpha, &a[m1], m, 
00418                                 &b[m1], ldb);
00419                         dgemm_("N", "N", &m1, n, &m2, &c_b23, a, m, &b[m1], 
00420                                 ldb, alpha, &b[b_offset], ldb);
00421                         dtrsm_("L", "L", "T", diag, &m1, n, &c_b27, &a[m2], m, 
00422                                  &b[b_offset], ldb);
00423 
00424                     }
00425 
00426                 }
00427 
00428             } else {
00429 
00430 /*              SIDE = 'L', N is odd, and TRANSR = 'T' */
00431 
00432                 if (lower) {
00433 
00434 /*                 SIDE  ='L', N is odd, TRANSR = 'T', and UPLO = 'L' */
00435 
00436                     if (notrans) {
00437 
00438 /*                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'L', and */
00439 /*                    TRANS = 'N' */
00440 
00441                         if (*m == 1) {
00442                             dtrsm_("L", "U", "T", diag, &m1, n, alpha, a, &m1, 
00443                                      &b[b_offset], ldb);
00444                         } else {
00445                             dtrsm_("L", "U", "T", diag, &m1, n, alpha, a, &m1, 
00446                                      &b[b_offset], ldb);
00447                             dgemm_("T", "N", &m2, n, &m1, &c_b23, &a[m1 * m1], 
00448                                      &m1, &b[b_offset], ldb, alpha, &b[m1], 
00449                                     ldb);
00450                             dtrsm_("L", "L", "N", diag, &m2, n, &c_b27, &a[1], 
00451                                      &m1, &b[m1], ldb);
00452                         }
00453 
00454                     } else {
00455 
00456 /*                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'L', and */
00457 /*                    TRANS = 'T' */
00458 
00459                         if (*m == 1) {
00460                             dtrsm_("L", "U", "N", diag, &m1, n, alpha, a, &m1, 
00461                                      &b[b_offset], ldb);
00462                         } else {
00463                             dtrsm_("L", "L", "T", diag, &m2, n, alpha, &a[1], 
00464                                     &m1, &b[m1], ldb);
00465                             dgemm_("N", "N", &m1, n, &m2, &c_b23, &a[m1 * m1], 
00466                                      &m1, &b[m1], ldb, alpha, &b[b_offset], 
00467                                     ldb);
00468                             dtrsm_("L", "U", "N", diag, &m1, n, &c_b27, a, &
00469                                     m1, &b[b_offset], ldb);
00470                         }
00471 
00472                     }
00473 
00474                 } else {
00475 
00476 /*                 SIDE  ='L', N is odd, TRANSR = 'T', and UPLO = 'U' */
00477 
00478                     if (! notrans) {
00479 
00480 /*                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'U', and */
00481 /*                    TRANS = 'N' */
00482 
00483                         dtrsm_("L", "U", "T", diag, &m1, n, alpha, &a[m2 * m2]
00484 , &m2, &b[b_offset], ldb);
00485                         dgemm_("N", "N", &m2, n, &m1, &c_b23, a, &m2, &b[
00486                                 b_offset], ldb, alpha, &b[m1], ldb);
00487                         dtrsm_("L", "L", "N", diag, &m2, n, &c_b27, &a[m1 * 
00488                                 m2], &m2, &b[m1], ldb);
00489 
00490                     } else {
00491 
00492 /*                    SIDE  ='L', N is odd, TRANSR = 'T', UPLO = 'U', and */
00493 /*                    TRANS = 'T' */
00494 
00495                         dtrsm_("L", "L", "T", diag, &m2, n, alpha, &a[m1 * m2]
00496 , &m2, &b[m1], ldb);
00497                         dgemm_("T", "N", &m1, n, &m2, &c_b23, a, &m2, &b[m1], 
00498                                 ldb, alpha, &b[b_offset], ldb);
00499                         dtrsm_("L", "U", "N", diag, &m1, n, &c_b27, &a[m2 * 
00500                                 m2], &m2, &b[b_offset], ldb);
00501 
00502                     }
00503 
00504                 }
00505 
00506             }
00507 
00508         } else {
00509 
00510 /*           SIDE = 'L' and N is even */
00511 
00512             if (normaltransr) {
00513 
00514 /*              SIDE = 'L', N is even, and TRANSR = 'N' */
00515 
00516                 if (lower) {
00517 
00518 /*                 SIDE  ='L', N is even, TRANSR = 'N', and UPLO = 'L' */
00519 
00520                     if (notrans) {
00521 
00522 /*                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'L', */
00523 /*                    and TRANS = 'N' */
00524 
00525                         i__1 = *m + 1;
00526                         dtrsm_("L", "L", "N", diag, &k, n, alpha, &a[1], &
00527                                 i__1, &b[b_offset], ldb);
00528                         i__1 = *m + 1;
00529                         dgemm_("N", "N", &k, n, &k, &c_b23, &a[k + 1], &i__1, 
00530                                 &b[b_offset], ldb, alpha, &b[k], ldb);
00531                         i__1 = *m + 1;
00532                         dtrsm_("L", "U", "T", diag, &k, n, &c_b27, a, &i__1, &
00533                                 b[k], ldb);
00534 
00535                     } else {
00536 
00537 /*                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'L', */
00538 /*                    and TRANS = 'T' */
00539 
00540                         i__1 = *m + 1;
00541                         dtrsm_("L", "U", "N", diag, &k, n, alpha, a, &i__1, &
00542                                 b[k], ldb);
00543                         i__1 = *m + 1;
00544                         dgemm_("T", "N", &k, n, &k, &c_b23, &a[k + 1], &i__1, 
00545                                 &b[k], ldb, alpha, &b[b_offset], ldb);
00546                         i__1 = *m + 1;
00547                         dtrsm_("L", "L", "T", diag, &k, n, &c_b27, &a[1], &
00548                                 i__1, &b[b_offset], ldb);
00549 
00550                     }
00551 
00552                 } else {
00553 
00554 /*                 SIDE  ='L', N is even, TRANSR = 'N', and UPLO = 'U' */
00555 
00556                     if (! notrans) {
00557 
00558 /*                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'U', */
00559 /*                    and TRANS = 'N' */
00560 
00561                         i__1 = *m + 1;
00562                         dtrsm_("L", "L", "N", diag, &k, n, alpha, &a[k + 1], &
00563                                 i__1, &b[b_offset], ldb);
00564                         i__1 = *m + 1;
00565                         dgemm_("T", "N", &k, n, &k, &c_b23, a, &i__1, &b[
00566                                 b_offset], ldb, alpha, &b[k], ldb);
00567                         i__1 = *m + 1;
00568                         dtrsm_("L", "U", "T", diag, &k, n, &c_b27, &a[k], &
00569                                 i__1, &b[k], ldb);
00570 
00571                     } else {
00572 
00573 /*                    SIDE  ='L', N is even, TRANSR = 'N', UPLO = 'U', */
00574 /*                    and TRANS = 'T' */
00575                         i__1 = *m + 1;
00576                         dtrsm_("L", "U", "N", diag, &k, n, alpha, &a[k], &
00577                                 i__1, &b[k], ldb);
00578                         i__1 = *m + 1;
00579                         dgemm_("N", "N", &k, n, &k, &c_b23, a, &i__1, &b[k], 
00580                                 ldb, alpha, &b[b_offset], ldb);
00581                         i__1 = *m + 1;
00582                         dtrsm_("L", "L", "T", diag, &k, n, &c_b27, &a[k + 1], 
00583                                 &i__1, &b[b_offset], ldb);
00584 
00585                     }
00586 
00587                 }
00588 
00589             } else {
00590 
00591 /*              SIDE = 'L', N is even, and TRANSR = 'T' */
00592 
00593                 if (lower) {
00594 
00595 /*                 SIDE  ='L', N is even, TRANSR = 'T', and UPLO = 'L' */
00596 
00597                     if (notrans) {
00598 
00599 /*                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'L', */
00600 /*                    and TRANS = 'N' */
00601 
00602                         dtrsm_("L", "U", "T", diag, &k, n, alpha, &a[k], &k, &
00603                                 b[b_offset], ldb);
00604                         dgemm_("T", "N", &k, n, &k, &c_b23, &a[k * (k + 1)], &
00605                                 k, &b[b_offset], ldb, alpha, &b[k], ldb);
00606                         dtrsm_("L", "L", "N", diag, &k, n, &c_b27, a, &k, &b[
00607                                 k], ldb);
00608 
00609                     } else {
00610 
00611 /*                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'L', */
00612 /*                    and TRANS = 'T' */
00613 
00614                         dtrsm_("L", "L", "T", diag, &k, n, alpha, a, &k, &b[k]
00615 , ldb);
00616                         dgemm_("N", "N", &k, n, &k, &c_b23, &a[k * (k + 1)], &
00617                                 k, &b[k], ldb, alpha, &b[b_offset], ldb);
00618                         dtrsm_("L", "U", "N", diag, &k, n, &c_b27, &a[k], &k, 
00619                                 &b[b_offset], ldb);
00620 
00621                     }
00622 
00623                 } else {
00624 
00625 /*                 SIDE  ='L', N is even, TRANSR = 'T', and UPLO = 'U' */
00626 
00627                     if (! notrans) {
00628 
00629 /*                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'U', */
00630 /*                    and TRANS = 'N' */
00631 
00632                         dtrsm_("L", "U", "T", diag, &k, n, alpha, &a[k * (k + 
00633                                 1)], &k, &b[b_offset], ldb);
00634                         dgemm_("N", "N", &k, n, &k, &c_b23, a, &k, &b[
00635                                 b_offset], ldb, alpha, &b[k], ldb);
00636                         dtrsm_("L", "L", "N", diag, &k, n, &c_b27, &a[k * k], 
00637                                 &k, &b[k], ldb);
00638 
00639                     } else {
00640 
00641 /*                    SIDE  ='L', N is even, TRANSR = 'T', UPLO = 'U', */
00642 /*                    and TRANS = 'T' */
00643 
00644                         dtrsm_("L", "L", "T", diag, &k, n, alpha, &a[k * k], &
00645                                 k, &b[k], ldb);
00646                         dgemm_("T", "N", &k, n, &k, &c_b23, a, &k, &b[k], ldb, 
00647                                  alpha, &b[b_offset], ldb);
00648                         dtrsm_("L", "U", "N", diag, &k, n, &c_b27, &a[k * (k 
00649                                 + 1)], &k, &b[b_offset], ldb);
00650 
00651                     }
00652 
00653                 }
00654 
00655             }
00656 
00657         }
00658 
00659     } else {
00660 
00661 /*        SIDE = 'R' */
00662 
00663 /*        A is N-by-N. */
00664 /*        If N is odd, set NISODD = .TRUE., and N1 and N2. */
00665 /*        If N is even, NISODD = .FALSE., and K. */
00666 
00667         if (*n % 2 == 0) {
00668             nisodd = FALSE_;
00669             k = *n / 2;
00670         } else {
00671             nisodd = TRUE_;
00672             if (lower) {
00673                 n2 = *n / 2;
00674                 n1 = *n - n2;
00675             } else {
00676                 n1 = *n / 2;
00677                 n2 = *n - n1;
00678             }
00679         }
00680 
00681         if (nisodd) {
00682 
00683 /*           SIDE = 'R' and N is odd */
00684 
00685             if (normaltransr) {
00686 
00687 /*              SIDE = 'R', N is odd, and TRANSR = 'N' */
00688 
00689                 if (lower) {
00690 
00691 /*                 SIDE  ='R', N is odd, TRANSR = 'N', and UPLO = 'L' */
00692 
00693                     if (notrans) {
00694 
00695 /*                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'L', and */
00696 /*                    TRANS = 'N' */
00697 
00698                         dtrsm_("R", "U", "T", diag, m, &n2, alpha, &a[*n], n, 
00699                                 &b[n1 * b_dim1], ldb);
00700                         dgemm_("N", "N", m, &n1, &n2, &c_b23, &b[n1 * b_dim1], 
00701                                  ldb, &a[n1], n, alpha, b, ldb);
00702                         dtrsm_("R", "L", "N", diag, m, &n1, &c_b27, a, n, b, 
00703                                 ldb);
00704 
00705                     } else {
00706 
00707 /*                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'L', and */
00708 /*                    TRANS = 'T' */
00709 
00710                         dtrsm_("R", "L", "T", diag, m, &n1, alpha, a, n, b, 
00711                                 ldb);
00712                         dgemm_("N", "T", m, &n2, &n1, &c_b23, b, ldb, &a[n1], 
00713                                 n, alpha, &b[n1 * b_dim1], ldb);
00714                         dtrsm_("R", "U", "N", diag, m, &n2, &c_b27, &a[*n], n, 
00715                                  &b[n1 * b_dim1], ldb);
00716 
00717                     }
00718 
00719                 } else {
00720 
00721 /*                 SIDE  ='R', N is odd, TRANSR = 'N', and UPLO = 'U' */
00722 
00723                     if (notrans) {
00724 
00725 /*                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'U', and */
00726 /*                    TRANS = 'N' */
00727 
00728                         dtrsm_("R", "L", "T", diag, m, &n1, alpha, &a[n2], n, 
00729                                 b, ldb);
00730                         dgemm_("N", "N", m, &n2, &n1, &c_b23, b, ldb, a, n, 
00731                                 alpha, &b[n1 * b_dim1], ldb);
00732                         dtrsm_("R", "U", "N", diag, m, &n2, &c_b27, &a[n1], n, 
00733                                  &b[n1 * b_dim1], ldb);
00734 
00735                     } else {
00736 
00737 /*                    SIDE  ='R', N is odd, TRANSR = 'N', UPLO = 'U', and */
00738 /*                    TRANS = 'T' */
00739 
00740                         dtrsm_("R", "U", "T", diag, m, &n2, alpha, &a[n1], n, 
00741                                 &b[n1 * b_dim1], ldb);
00742                         dgemm_("N", "T", m, &n1, &n2, &c_b23, &b[n1 * b_dim1], 
00743                                  ldb, a, n, alpha, b, ldb);
00744                         dtrsm_("R", "L", "N", diag, m, &n1, &c_b27, &a[n2], n, 
00745                                  b, ldb);
00746 
00747                     }
00748 
00749                 }
00750 
00751             } else {
00752 
00753 /*              SIDE = 'R', N is odd, and TRANSR = 'T' */
00754 
00755                 if (lower) {
00756 
00757 /*                 SIDE  ='R', N is odd, TRANSR = 'T', and UPLO = 'L' */
00758 
00759                     if (notrans) {
00760 
00761 /*                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'L', and */
00762 /*                    TRANS = 'N' */
00763 
00764                         dtrsm_("R", "L", "N", diag, m, &n2, alpha, &a[1], &n1, 
00765                                  &b[n1 * b_dim1], ldb);
00766                         dgemm_("N", "T", m, &n1, &n2, &c_b23, &b[n1 * b_dim1], 
00767                                  ldb, &a[n1 * n1], &n1, alpha, b, ldb);
00768                         dtrsm_("R", "U", "T", diag, m, &n1, &c_b27, a, &n1, b, 
00769                                  ldb);
00770 
00771                     } else {
00772 
00773 /*                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'L', and */
00774 /*                    TRANS = 'T' */
00775 
00776                         dtrsm_("R", "U", "N", diag, m, &n1, alpha, a, &n1, b, 
00777                                 ldb);
00778                         dgemm_("N", "N", m, &n2, &n1, &c_b23, b, ldb, &a[n1 * 
00779                                 n1], &n1, alpha, &b[n1 * b_dim1], ldb);
00780                         dtrsm_("R", "L", "T", diag, m, &n2, &c_b27, &a[1], &
00781                                 n1, &b[n1 * b_dim1], ldb);
00782 
00783                     }
00784 
00785                 } else {
00786 
00787 /*                 SIDE  ='R', N is odd, TRANSR = 'T', and UPLO = 'U' */
00788 
00789                     if (notrans) {
00790 
00791 /*                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'U', and */
00792 /*                    TRANS = 'N' */
00793 
00794                         dtrsm_("R", "U", "N", diag, m, &n1, alpha, &a[n2 * n2]
00795 , &n2, b, ldb);
00796                         dgemm_("N", "T", m, &n2, &n1, &c_b23, b, ldb, a, &n2, 
00797                                 alpha, &b[n1 * b_dim1], ldb);
00798                         dtrsm_("R", "L", "T", diag, m, &n2, &c_b27, &a[n1 * 
00799                                 n2], &n2, &b[n1 * b_dim1], ldb);
00800 
00801                     } else {
00802 
00803 /*                    SIDE  ='R', N is odd, TRANSR = 'T', UPLO = 'U', and */
00804 /*                    TRANS = 'T' */
00805 
00806                         dtrsm_("R", "L", "N", diag, m, &n2, alpha, &a[n1 * n2]
00807 , &n2, &b[n1 * b_dim1], ldb);
00808                         dgemm_("N", "N", m, &n1, &n2, &c_b23, &b[n1 * b_dim1], 
00809                                  ldb, a, &n2, alpha, b, ldb);
00810                         dtrsm_("R", "U", "T", diag, m, &n1, &c_b27, &a[n2 * 
00811                                 n2], &n2, b, ldb);
00812 
00813                     }
00814 
00815                 }
00816 
00817             }
00818 
00819         } else {
00820 
00821 /*           SIDE = 'R' and N is even */
00822 
00823             if (normaltransr) {
00824 
00825 /*              SIDE = 'R', N is even, and TRANSR = 'N' */
00826 
00827                 if (lower) {
00828 
00829 /*                 SIDE  ='R', N is even, TRANSR = 'N', and UPLO = 'L' */
00830 
00831                     if (notrans) {
00832 
00833 /*                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'L', */
00834 /*                    and TRANS = 'N' */
00835 
00836                         i__1 = *n + 1;
00837                         dtrsm_("R", "U", "T", diag, m, &k, alpha, a, &i__1, &
00838                                 b[k * b_dim1], ldb);
00839                         i__1 = *n + 1;
00840                         dgemm_("N", "N", m, &k, &k, &c_b23, &b[k * b_dim1], 
00841                                 ldb, &a[k + 1], &i__1, alpha, b, ldb);
00842                         i__1 = *n + 1;
00843                         dtrsm_("R", "L", "N", diag, m, &k, &c_b27, &a[1], &
00844                                 i__1, b, ldb);
00845 
00846                     } else {
00847 
00848 /*                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'L', */
00849 /*                    and TRANS = 'T' */
00850 
00851                         i__1 = *n + 1;
00852                         dtrsm_("R", "L", "T", diag, m, &k, alpha, &a[1], &
00853                                 i__1, b, ldb);
00854                         i__1 = *n + 1;
00855                         dgemm_("N", "T", m, &k, &k, &c_b23, b, ldb, &a[k + 1], 
00856                                  &i__1, alpha, &b[k * b_dim1], ldb);
00857                         i__1 = *n + 1;
00858                         dtrsm_("R", "U", "N", diag, m, &k, &c_b27, a, &i__1, &
00859                                 b[k * b_dim1], ldb);
00860 
00861                     }
00862 
00863                 } else {
00864 
00865 /*                 SIDE  ='R', N is even, TRANSR = 'N', and UPLO = 'U' */
00866 
00867                     if (notrans) {
00868 
00869 /*                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'U', */
00870 /*                    and TRANS = 'N' */
00871 
00872                         i__1 = *n + 1;
00873                         dtrsm_("R", "L", "T", diag, m, &k, alpha, &a[k + 1], &
00874                                 i__1, b, ldb);
00875                         i__1 = *n + 1;
00876                         dgemm_("N", "N", m, &k, &k, &c_b23, b, ldb, a, &i__1, 
00877                                 alpha, &b[k * b_dim1], ldb);
00878                         i__1 = *n + 1;
00879                         dtrsm_("R", "U", "N", diag, m, &k, &c_b27, &a[k], &
00880                                 i__1, &b[k * b_dim1], ldb);
00881 
00882                     } else {
00883 
00884 /*                    SIDE  ='R', N is even, TRANSR = 'N', UPLO = 'U', */
00885 /*                    and TRANS = 'T' */
00886 
00887                         i__1 = *n + 1;
00888                         dtrsm_("R", "U", "T", diag, m, &k, alpha, &a[k], &
00889                                 i__1, &b[k * b_dim1], ldb);
00890                         i__1 = *n + 1;
00891                         dgemm_("N", "T", m, &k, &k, &c_b23, &b[k * b_dim1], 
00892                                 ldb, a, &i__1, alpha, b, ldb);
00893                         i__1 = *n + 1;
00894                         dtrsm_("R", "L", "N", diag, m, &k, &c_b27, &a[k + 1], 
00895                                 &i__1, b, ldb);
00896 
00897                     }
00898 
00899                 }
00900 
00901             } else {
00902 
00903 /*              SIDE = 'R', N is even, and TRANSR = 'T' */
00904 
00905                 if (lower) {
00906 
00907 /*                 SIDE  ='R', N is even, TRANSR = 'T', and UPLO = 'L' */
00908 
00909                     if (notrans) {
00910 
00911 /*                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'L', */
00912 /*                    and TRANS = 'N' */
00913 
00914                         dtrsm_("R", "L", "N", diag, m, &k, alpha, a, &k, &b[k 
00915                                 * b_dim1], ldb);
00916                         dgemm_("N", "T", m, &k, &k, &c_b23, &b[k * b_dim1], 
00917                                 ldb, &a[(k + 1) * k], &k, alpha, b, ldb);
00918                         dtrsm_("R", "U", "T", diag, m, &k, &c_b27, &a[k], &k, 
00919                                 b, ldb);
00920 
00921                     } else {
00922 
00923 /*                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'L', */
00924 /*                    and TRANS = 'T' */
00925 
00926                         dtrsm_("R", "U", "N", diag, m, &k, alpha, &a[k], &k, 
00927                                 b, ldb);
00928                         dgemm_("N", "N", m, &k, &k, &c_b23, b, ldb, &a[(k + 1)
00929                                  * k], &k, alpha, &b[k * b_dim1], ldb);
00930                         dtrsm_("R", "L", "T", diag, m, &k, &c_b27, a, &k, &b[
00931                                 k * b_dim1], ldb);
00932 
00933                     }
00934 
00935                 } else {
00936 
00937 /*                 SIDE  ='R', N is even, TRANSR = 'T', and UPLO = 'U' */
00938 
00939                     if (notrans) {
00940 
00941 /*                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'U', */
00942 /*                    and TRANS = 'N' */
00943 
00944                         dtrsm_("R", "U", "N", diag, m, &k, alpha, &a[(k + 1) *
00945                                  k], &k, b, ldb);
00946                         dgemm_("N", "T", m, &k, &k, &c_b23, b, ldb, a, &k, 
00947                                 alpha, &b[k * b_dim1], ldb);
00948                         dtrsm_("R", "L", "T", diag, m, &k, &c_b27, &a[k * k], 
00949                                 &k, &b[k * b_dim1], ldb);
00950 
00951                     } else {
00952 
00953 /*                    SIDE  ='R', N is even, TRANSR = 'T', UPLO = 'U', */
00954 /*                    and TRANS = 'T' */
00955 
00956                         dtrsm_("R", "L", "N", diag, m, &k, alpha, &a[k * k], &
00957                                 k, &b[k * b_dim1], ldb);
00958                         dgemm_("N", "N", m, &k, &k, &c_b23, &b[k * b_dim1], 
00959                                 ldb, a, &k, alpha, b, ldb);
00960                         dtrsm_("R", "U", "T", diag, m, &k, &c_b27, &a[(k + 1) 
00961                                 * k], &k, b, ldb);
00962 
00963                     }
00964 
00965                 }
00966 
00967             }
00968 
00969         }
00970     }
00971 
00972     return 0;
00973 
00974 /*     End of DTFSM */
00975 
00976 } /* dtfsm_ */

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