slantb.c
Go to the documentation of this file.
00001 /* slantb.f -- translated by f2c (version 20061008).
00002    You must link the resulting object file with libf2c:
00003         on Microsoft Windows system, link with libf2c.lib;
00004         on Linux or Unix systems, link with .../path/to/libf2c.a -lm
00005         or, if you install libf2c.a in a standard place, with -lf2c -lm
00006         -- in that order, at the end of the command line, as in
00007                 cc *.o -lf2c -lm
00008         Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
00009 
00010                 http://www.netlib.org/f2c/libf2c.zip
00011 */
00012 
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015 
00016 /* Table of constant values */
00017 
00018 static integer c__1 = 1;
00019 
00020 doublereal slantb_(char *norm, char *uplo, char *diag, integer *n, integer *k, 
00021          real *ab, integer *ldab, real *work)
00022 {
00023     /* System generated locals */
00024     integer ab_dim1, ab_offset, i__1, i__2, i__3, i__4, i__5;
00025     real ret_val, r__1, r__2, r__3;
00026 
00027     /* Builtin functions */
00028     double sqrt(doublereal);
00029 
00030     /* Local variables */
00031     integer i__, j, l;
00032     real sum, scale;
00033     logical udiag;
00034     extern logical lsame_(char *, char *);
00035     real value;
00036     extern /* Subroutine */ int slassq_(integer *, real *, integer *, real *, 
00037             real *);
00038 
00039 
00040 /*  -- LAPACK auxiliary routine (version 3.2) -- */
00041 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00042 /*     November 2006 */
00043 
00044 /*     .. Scalar Arguments .. */
00045 /*     .. */
00046 /*     .. Array Arguments .. */
00047 /*     .. */
00048 
00049 /*  Purpose */
00050 /*  ======= */
00051 
00052 /*  SLANTB  returns the value of the one norm,  or the Frobenius norm, or */
00053 /*  the  infinity norm,  or the element of  largest absolute value  of an */
00054 /*  n by n triangular band matrix A,  with ( k + 1 ) diagonals. */
00055 
00056 /*  Description */
00057 /*  =========== */
00058 
00059 /*  SLANTB returns the value */
00060 
00061 /*     SLANTB = ( max(abs(A(i,j))), NORM = 'M' or 'm' */
00062 /*              ( */
00063 /*              ( norm1(A),         NORM = '1', 'O' or 'o' */
00064 /*              ( */
00065 /*              ( normI(A),         NORM = 'I' or 'i' */
00066 /*              ( */
00067 /*              ( normF(A),         NORM = 'F', 'f', 'E' or 'e' */
00068 
00069 /*  where  norm1  denotes the  one norm of a matrix (maximum column sum), */
00070 /*  normI  denotes the  infinity norm  of a matrix  (maximum row sum) and */
00071 /*  normF  denotes the  Frobenius norm of a matrix (square root of sum of */
00072 /*  squares).  Note that  max(abs(A(i,j)))  is not a consistent matrix norm. */
00073 
00074 /*  Arguments */
00075 /*  ========= */
00076 
00077 /*  NORM    (input) CHARACTER*1 */
00078 /*          Specifies the value to be returned in SLANTB as described */
00079 /*          above. */
00080 
00081 /*  UPLO    (input) CHARACTER*1 */
00082 /*          Specifies whether the matrix A is upper or lower triangular. */
00083 /*          = 'U':  Upper triangular */
00084 /*          = 'L':  Lower triangular */
00085 
00086 /*  DIAG    (input) CHARACTER*1 */
00087 /*          Specifies whether or not the matrix A is unit triangular. */
00088 /*          = 'N':  Non-unit triangular */
00089 /*          = 'U':  Unit triangular */
00090 
00091 /*  N       (input) INTEGER */
00092 /*          The order of the matrix A.  N >= 0.  When N = 0, SLANTB is */
00093 /*          set to zero. */
00094 
00095 /*  K       (input) INTEGER */
00096 /*          The number of super-diagonals of the matrix A if UPLO = 'U', */
00097 /*          or the number of sub-diagonals of the matrix A if UPLO = 'L'. */
00098 /*          K >= 0. */
00099 
00100 /*  AB      (input) REAL array, dimension (LDAB,N) */
00101 /*          The upper or lower triangular band matrix A, stored in the */
00102 /*          first k+1 rows of AB.  The j-th column of A is stored */
00103 /*          in the j-th column of the array AB as follows: */
00104 /*          if UPLO = 'U', AB(k+1+i-j,j) = A(i,j) for max(1,j-k)<=i<=j; */
00105 /*          if UPLO = 'L', AB(1+i-j,j)   = A(i,j) for j<=i<=min(n,j+k). */
00106 /*          Note that when DIAG = 'U', the elements of the array AB */
00107 /*          corresponding to the diagonal elements of the matrix A are */
00108 /*          not referenced, but are assumed to be one. */
00109 
00110 /*  LDAB    (input) INTEGER */
00111 /*          The leading dimension of the array AB.  LDAB >= K+1. */
00112 
00113 /*  WORK    (workspace) REAL array, dimension (MAX(1,LWORK)), */
00114 /*          where LWORK >= N when NORM = 'I'; otherwise, WORK is not */
00115 /*          referenced. */
00116 
00117 /* ===================================================================== */
00118 
00119 /*     .. Parameters .. */
00120 /*     .. */
00121 /*     .. Local Scalars .. */
00122 /*     .. */
00123 /*     .. External Subroutines .. */
00124 /*     .. */
00125 /*     .. External Functions .. */
00126 /*     .. */
00127 /*     .. Intrinsic Functions .. */
00128 /*     .. */
00129 /*     .. Executable Statements .. */
00130 
00131     /* Parameter adjustments */
00132     ab_dim1 = *ldab;
00133     ab_offset = 1 + ab_dim1;
00134     ab -= ab_offset;
00135     --work;
00136 
00137     /* Function Body */
00138     if (*n == 0) {
00139         value = 0.f;
00140     } else if (lsame_(norm, "M")) {
00141 
00142 /*        Find max(abs(A(i,j))). */
00143 
00144         if (lsame_(diag, "U")) {
00145             value = 1.f;
00146             if (lsame_(uplo, "U")) {
00147                 i__1 = *n;
00148                 for (j = 1; j <= i__1; ++j) {
00149 /* Computing MAX */
00150                     i__2 = *k + 2 - j;
00151                     i__3 = *k;
00152                     for (i__ = max(i__2,1); i__ <= i__3; ++i__) {
00153 /* Computing MAX */
00154                         r__2 = value, r__3 = (r__1 = ab[i__ + j * ab_dim1], 
00155                                 dabs(r__1));
00156                         value = dmax(r__2,r__3);
00157 /* L10: */
00158                     }
00159 /* L20: */
00160                 }
00161             } else {
00162                 i__1 = *n;
00163                 for (j = 1; j <= i__1; ++j) {
00164 /* Computing MIN */
00165                     i__2 = *n + 1 - j, i__4 = *k + 1;
00166                     i__3 = min(i__2,i__4);
00167                     for (i__ = 2; i__ <= i__3; ++i__) {
00168 /* Computing MAX */
00169                         r__2 = value, r__3 = (r__1 = ab[i__ + j * ab_dim1], 
00170                                 dabs(r__1));
00171                         value = dmax(r__2,r__3);
00172 /* L30: */
00173                     }
00174 /* L40: */
00175                 }
00176             }
00177         } else {
00178             value = 0.f;
00179             if (lsame_(uplo, "U")) {
00180                 i__1 = *n;
00181                 for (j = 1; j <= i__1; ++j) {
00182 /* Computing MAX */
00183                     i__3 = *k + 2 - j;
00184                     i__2 = *k + 1;
00185                     for (i__ = max(i__3,1); i__ <= i__2; ++i__) {
00186 /* Computing MAX */
00187                         r__2 = value, r__3 = (r__1 = ab[i__ + j * ab_dim1], 
00188                                 dabs(r__1));
00189                         value = dmax(r__2,r__3);
00190 /* L50: */
00191                     }
00192 /* L60: */
00193                 }
00194             } else {
00195                 i__1 = *n;
00196                 for (j = 1; j <= i__1; ++j) {
00197 /* Computing MIN */
00198                     i__3 = *n + 1 - j, i__4 = *k + 1;
00199                     i__2 = min(i__3,i__4);
00200                     for (i__ = 1; i__ <= i__2; ++i__) {
00201 /* Computing MAX */
00202                         r__2 = value, r__3 = (r__1 = ab[i__ + j * ab_dim1], 
00203                                 dabs(r__1));
00204                         value = dmax(r__2,r__3);
00205 /* L70: */
00206                     }
00207 /* L80: */
00208                 }
00209             }
00210         }
00211     } else if (lsame_(norm, "O") || *(unsigned char *)
00212             norm == '1') {
00213 
00214 /*        Find norm1(A). */
00215 
00216         value = 0.f;
00217         udiag = lsame_(diag, "U");
00218         if (lsame_(uplo, "U")) {
00219             i__1 = *n;
00220             for (j = 1; j <= i__1; ++j) {
00221                 if (udiag) {
00222                     sum = 1.f;
00223 /* Computing MAX */
00224                     i__2 = *k + 2 - j;
00225                     i__3 = *k;
00226                     for (i__ = max(i__2,1); i__ <= i__3; ++i__) {
00227                         sum += (r__1 = ab[i__ + j * ab_dim1], dabs(r__1));
00228 /* L90: */
00229                     }
00230                 } else {
00231                     sum = 0.f;
00232 /* Computing MAX */
00233                     i__3 = *k + 2 - j;
00234                     i__2 = *k + 1;
00235                     for (i__ = max(i__3,1); i__ <= i__2; ++i__) {
00236                         sum += (r__1 = ab[i__ + j * ab_dim1], dabs(r__1));
00237 /* L100: */
00238                     }
00239                 }
00240                 value = dmax(value,sum);
00241 /* L110: */
00242             }
00243         } else {
00244             i__1 = *n;
00245             for (j = 1; j <= i__1; ++j) {
00246                 if (udiag) {
00247                     sum = 1.f;
00248 /* Computing MIN */
00249                     i__3 = *n + 1 - j, i__4 = *k + 1;
00250                     i__2 = min(i__3,i__4);
00251                     for (i__ = 2; i__ <= i__2; ++i__) {
00252                         sum += (r__1 = ab[i__ + j * ab_dim1], dabs(r__1));
00253 /* L120: */
00254                     }
00255                 } else {
00256                     sum = 0.f;
00257 /* Computing MIN */
00258                     i__3 = *n + 1 - j, i__4 = *k + 1;
00259                     i__2 = min(i__3,i__4);
00260                     for (i__ = 1; i__ <= i__2; ++i__) {
00261                         sum += (r__1 = ab[i__ + j * ab_dim1], dabs(r__1));
00262 /* L130: */
00263                     }
00264                 }
00265                 value = dmax(value,sum);
00266 /* L140: */
00267             }
00268         }
00269     } else if (lsame_(norm, "I")) {
00270 
00271 /*        Find normI(A). */
00272 
00273         value = 0.f;
00274         if (lsame_(uplo, "U")) {
00275             if (lsame_(diag, "U")) {
00276                 i__1 = *n;
00277                 for (i__ = 1; i__ <= i__1; ++i__) {
00278                     work[i__] = 1.f;
00279 /* L150: */
00280                 }
00281                 i__1 = *n;
00282                 for (j = 1; j <= i__1; ++j) {
00283                     l = *k + 1 - j;
00284 /* Computing MAX */
00285                     i__2 = 1, i__3 = j - *k;
00286                     i__4 = j - 1;
00287                     for (i__ = max(i__2,i__3); i__ <= i__4; ++i__) {
00288                         work[i__] += (r__1 = ab[l + i__ + j * ab_dim1], dabs(
00289                                 r__1));
00290 /* L160: */
00291                     }
00292 /* L170: */
00293                 }
00294             } else {
00295                 i__1 = *n;
00296                 for (i__ = 1; i__ <= i__1; ++i__) {
00297                     work[i__] = 0.f;
00298 /* L180: */
00299                 }
00300                 i__1 = *n;
00301                 for (j = 1; j <= i__1; ++j) {
00302                     l = *k + 1 - j;
00303 /* Computing MAX */
00304                     i__4 = 1, i__2 = j - *k;
00305                     i__3 = j;
00306                     for (i__ = max(i__4,i__2); i__ <= i__3; ++i__) {
00307                         work[i__] += (r__1 = ab[l + i__ + j * ab_dim1], dabs(
00308                                 r__1));
00309 /* L190: */
00310                     }
00311 /* L200: */
00312                 }
00313             }
00314         } else {
00315             if (lsame_(diag, "U")) {
00316                 i__1 = *n;
00317                 for (i__ = 1; i__ <= i__1; ++i__) {
00318                     work[i__] = 1.f;
00319 /* L210: */
00320                 }
00321                 i__1 = *n;
00322                 for (j = 1; j <= i__1; ++j) {
00323                     l = 1 - j;
00324 /* Computing MIN */
00325                     i__4 = *n, i__2 = j + *k;
00326                     i__3 = min(i__4,i__2);
00327                     for (i__ = j + 1; i__ <= i__3; ++i__) {
00328                         work[i__] += (r__1 = ab[l + i__ + j * ab_dim1], dabs(
00329                                 r__1));
00330 /* L220: */
00331                     }
00332 /* L230: */
00333                 }
00334             } else {
00335                 i__1 = *n;
00336                 for (i__ = 1; i__ <= i__1; ++i__) {
00337                     work[i__] = 0.f;
00338 /* L240: */
00339                 }
00340                 i__1 = *n;
00341                 for (j = 1; j <= i__1; ++j) {
00342                     l = 1 - j;
00343 /* Computing MIN */
00344                     i__4 = *n, i__2 = j + *k;
00345                     i__3 = min(i__4,i__2);
00346                     for (i__ = j; i__ <= i__3; ++i__) {
00347                         work[i__] += (r__1 = ab[l + i__ + j * ab_dim1], dabs(
00348                                 r__1));
00349 /* L250: */
00350                     }
00351 /* L260: */
00352                 }
00353             }
00354         }
00355         i__1 = *n;
00356         for (i__ = 1; i__ <= i__1; ++i__) {
00357 /* Computing MAX */
00358             r__1 = value, r__2 = work[i__];
00359             value = dmax(r__1,r__2);
00360 /* L270: */
00361         }
00362     } else if (lsame_(norm, "F") || lsame_(norm, "E")) {
00363 
00364 /*        Find normF(A). */
00365 
00366         if (lsame_(uplo, "U")) {
00367             if (lsame_(diag, "U")) {
00368                 scale = 1.f;
00369                 sum = (real) (*n);
00370                 if (*k > 0) {
00371                     i__1 = *n;
00372                     for (j = 2; j <= i__1; ++j) {
00373 /* Computing MIN */
00374                         i__4 = j - 1;
00375                         i__3 = min(i__4,*k);
00376 /* Computing MAX */
00377                         i__2 = *k + 2 - j;
00378                         slassq_(&i__3, &ab[max(i__2, 1)+ j * ab_dim1], &c__1, 
00379                                 &scale, &sum);
00380 /* L280: */
00381                     }
00382                 }
00383             } else {
00384                 scale = 0.f;
00385                 sum = 1.f;
00386                 i__1 = *n;
00387                 for (j = 1; j <= i__1; ++j) {
00388 /* Computing MIN */
00389                     i__4 = j, i__2 = *k + 1;
00390                     i__3 = min(i__4,i__2);
00391 /* Computing MAX */
00392                     i__5 = *k + 2 - j;
00393                     slassq_(&i__3, &ab[max(i__5, 1)+ j * ab_dim1], &c__1, &
00394                             scale, &sum);
00395 /* L290: */
00396                 }
00397             }
00398         } else {
00399             if (lsame_(diag, "U")) {
00400                 scale = 1.f;
00401                 sum = (real) (*n);
00402                 if (*k > 0) {
00403                     i__1 = *n - 1;
00404                     for (j = 1; j <= i__1; ++j) {
00405 /* Computing MIN */
00406                         i__4 = *n - j;
00407                         i__3 = min(i__4,*k);
00408                         slassq_(&i__3, &ab[j * ab_dim1 + 2], &c__1, &scale, &
00409                                 sum);
00410 /* L300: */
00411                     }
00412                 }
00413             } else {
00414                 scale = 0.f;
00415                 sum = 1.f;
00416                 i__1 = *n;
00417                 for (j = 1; j <= i__1; ++j) {
00418 /* Computing MIN */
00419                     i__4 = *n - j + 1, i__2 = *k + 1;
00420                     i__3 = min(i__4,i__2);
00421                     slassq_(&i__3, &ab[j * ab_dim1 + 1], &c__1, &scale, &sum);
00422 /* L310: */
00423                 }
00424             }
00425         }
00426         value = scale * sqrt(sum);
00427     }
00428 
00429     ret_val = value;
00430     return ret_val;
00431 
00432 /*     End of SLANTB */
00433 
00434 } /* slantb_ */


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