zlansb.c
Go to the documentation of this file.
00001 /* zlansb.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 zlansb_(char *norm, char *uplo, integer *n, integer *k, 
00021         doublecomplex *ab, integer *ldab, doublereal *work)
00022 {
00023     /* System generated locals */
00024     integer ab_dim1, ab_offset, i__1, i__2, i__3, i__4;
00025     doublereal ret_val, d__1, d__2;
00026 
00027     /* Builtin functions */
00028     double z_abs(doublecomplex *), sqrt(doublereal);
00029 
00030     /* Local variables */
00031     integer i__, j, l;
00032     doublereal sum, absa, scale;
00033     extern logical lsame_(char *, char *);
00034     doublereal value;
00035     extern /* Subroutine */ int zlassq_(integer *, doublecomplex *, integer *, 
00036              doublereal *, doublereal *);
00037 
00038 
00039 /*  -- LAPACK auxiliary routine (version 3.2) -- */
00040 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00041 /*     November 2006 */
00042 
00043 /*     .. Scalar Arguments .. */
00044 /*     .. */
00045 /*     .. Array Arguments .. */
00046 /*     .. */
00047 
00048 /*  Purpose */
00049 /*  ======= */
00050 
00051 /*  ZLANSB  returns the value of the one norm,  or the Frobenius norm, or */
00052 /*  the  infinity norm,  or the element of  largest absolute value  of an */
00053 /*  n by n symmetric band matrix A,  with k super-diagonals. */
00054 
00055 /*  Description */
00056 /*  =========== */
00057 
00058 /*  ZLANSB returns the value */
00059 
00060 /*     ZLANSB = ( max(abs(A(i,j))), NORM = 'M' or 'm' */
00061 /*              ( */
00062 /*              ( norm1(A),         NORM = '1', 'O' or 'o' */
00063 /*              ( */
00064 /*              ( normI(A),         NORM = 'I' or 'i' */
00065 /*              ( */
00066 /*              ( normF(A),         NORM = 'F', 'f', 'E' or 'e' */
00067 
00068 /*  where  norm1  denotes the  one norm of a matrix (maximum column sum), */
00069 /*  normI  denotes the  infinity norm  of a matrix  (maximum row sum) and */
00070 /*  normF  denotes the  Frobenius norm of a matrix (square root of sum of */
00071 /*  squares).  Note that  max(abs(A(i,j)))  is not a consistent matrix norm. */
00072 
00073 /*  Arguments */
00074 /*  ========= */
00075 
00076 /*  NORM    (input) CHARACTER*1 */
00077 /*          Specifies the value to be returned in ZLANSB as described */
00078 /*          above. */
00079 
00080 /*  UPLO    (input) CHARACTER*1 */
00081 /*          Specifies whether the upper or lower triangular part of the */
00082 /*          band matrix A is supplied. */
00083 /*          = 'U':  Upper triangular part is supplied */
00084 /*          = 'L':  Lower triangular part is supplied */
00085 
00086 /*  N       (input) INTEGER */
00087 /*          The order of the matrix A.  N >= 0.  When N = 0, ZLANSB is */
00088 /*          set to zero. */
00089 
00090 /*  K       (input) INTEGER */
00091 /*          The number of super-diagonals or sub-diagonals of the */
00092 /*          band matrix A.  K >= 0. */
00093 
00094 /*  AB      (input) COMPLEX*16 array, dimension (LDAB,N) */
00095 /*          The upper or lower triangle of the symmetric band matrix A, */
00096 /*          stored in the first K+1 rows of AB.  The j-th column of A is */
00097 /*          stored in the j-th column of the array AB as follows: */
00098 /*          if UPLO = 'U', AB(k+1+i-j,j) = A(i,j) for max(1,j-k)<=i<=j; */
00099 /*          if UPLO = 'L', AB(1+i-j,j)   = A(i,j) for j<=i<=min(n,j+k). */
00100 
00101 /*  LDAB    (input) INTEGER */
00102 /*          The leading dimension of the array AB.  LDAB >= K+1. */
00103 
00104 /*  WORK    (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK)), */
00105 /*          where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise, */
00106 /*          WORK is not referenced. */
00107 
00108 /* ===================================================================== */
00109 
00110 /*     .. Parameters .. */
00111 /*     .. */
00112 /*     .. Local Scalars .. */
00113 /*     .. */
00114 /*     .. External Functions .. */
00115 /*     .. */
00116 /*     .. External Subroutines .. */
00117 /*     .. */
00118 /*     .. Intrinsic Functions .. */
00119 /*     .. */
00120 /*     .. Executable Statements .. */
00121 
00122     /* Parameter adjustments */
00123     ab_dim1 = *ldab;
00124     ab_offset = 1 + ab_dim1;
00125     ab -= ab_offset;
00126     --work;
00127 
00128     /* Function Body */
00129     if (*n == 0) {
00130         value = 0.;
00131     } else if (lsame_(norm, "M")) {
00132 
00133 /*        Find max(abs(A(i,j))). */
00134 
00135         value = 0.;
00136         if (lsame_(uplo, "U")) {
00137             i__1 = *n;
00138             for (j = 1; j <= i__1; ++j) {
00139 /* Computing MAX */
00140                 i__2 = *k + 2 - j;
00141                 i__3 = *k + 1;
00142                 for (i__ = max(i__2,1); i__ <= i__3; ++i__) {
00143 /* Computing MAX */
00144                     d__1 = value, d__2 = z_abs(&ab[i__ + j * ab_dim1]);
00145                     value = max(d__1,d__2);
00146 /* L10: */
00147                 }
00148 /* L20: */
00149             }
00150         } else {
00151             i__1 = *n;
00152             for (j = 1; j <= i__1; ++j) {
00153 /* Computing MIN */
00154                 i__2 = *n + 1 - j, i__4 = *k + 1;
00155                 i__3 = min(i__2,i__4);
00156                 for (i__ = 1; i__ <= i__3; ++i__) {
00157 /* Computing MAX */
00158                     d__1 = value, d__2 = z_abs(&ab[i__ + j * ab_dim1]);
00159                     value = max(d__1,d__2);
00160 /* L30: */
00161                 }
00162 /* L40: */
00163             }
00164         }
00165     } else if (lsame_(norm, "I") || lsame_(norm, "O") || *(unsigned char *)norm == '1') {
00166 
00167 /*        Find normI(A) ( = norm1(A), since A is symmetric). */
00168 
00169         value = 0.;
00170         if (lsame_(uplo, "U")) {
00171             i__1 = *n;
00172             for (j = 1; j <= i__1; ++j) {
00173                 sum = 0.;
00174                 l = *k + 1 - j;
00175 /* Computing MAX */
00176                 i__3 = 1, i__2 = j - *k;
00177                 i__4 = j - 1;
00178                 for (i__ = max(i__3,i__2); i__ <= i__4; ++i__) {
00179                     absa = z_abs(&ab[l + i__ + j * ab_dim1]);
00180                     sum += absa;
00181                     work[i__] += absa;
00182 /* L50: */
00183                 }
00184                 work[j] = sum + z_abs(&ab[*k + 1 + j * ab_dim1]);
00185 /* L60: */
00186             }
00187             i__1 = *n;
00188             for (i__ = 1; i__ <= i__1; ++i__) {
00189 /* Computing MAX */
00190                 d__1 = value, d__2 = work[i__];
00191                 value = max(d__1,d__2);
00192 /* L70: */
00193             }
00194         } else {
00195             i__1 = *n;
00196             for (i__ = 1; i__ <= i__1; ++i__) {
00197                 work[i__] = 0.;
00198 /* L80: */
00199             }
00200             i__1 = *n;
00201             for (j = 1; j <= i__1; ++j) {
00202                 sum = work[j] + z_abs(&ab[j * ab_dim1 + 1]);
00203                 l = 1 - j;
00204 /* Computing MIN */
00205                 i__3 = *n, i__2 = j + *k;
00206                 i__4 = min(i__3,i__2);
00207                 for (i__ = j + 1; i__ <= i__4; ++i__) {
00208                     absa = z_abs(&ab[l + i__ + j * ab_dim1]);
00209                     sum += absa;
00210                     work[i__] += absa;
00211 /* L90: */
00212                 }
00213                 value = max(value,sum);
00214 /* L100: */
00215             }
00216         }
00217     } else if (lsame_(norm, "F") || lsame_(norm, "E")) {
00218 
00219 /*        Find normF(A). */
00220 
00221         scale = 0.;
00222         sum = 1.;
00223         if (*k > 0) {
00224             if (lsame_(uplo, "U")) {
00225                 i__1 = *n;
00226                 for (j = 2; j <= i__1; ++j) {
00227 /* Computing MIN */
00228                     i__3 = j - 1;
00229                     i__4 = min(i__3,*k);
00230 /* Computing MAX */
00231                     i__2 = *k + 2 - j;
00232                     zlassq_(&i__4, &ab[max(i__2, 1)+ j * ab_dim1], &c__1, &
00233                             scale, &sum);
00234 /* L110: */
00235                 }
00236                 l = *k + 1;
00237             } else {
00238                 i__1 = *n - 1;
00239                 for (j = 1; j <= i__1; ++j) {
00240 /* Computing MIN */
00241                     i__3 = *n - j;
00242                     i__4 = min(i__3,*k);
00243                     zlassq_(&i__4, &ab[j * ab_dim1 + 2], &c__1, &scale, &sum);
00244 /* L120: */
00245                 }
00246                 l = 1;
00247             }
00248             sum *= 2;
00249         } else {
00250             l = 1;
00251         }
00252         zlassq_(n, &ab[l + ab_dim1], ldab, &scale, &sum);
00253         value = scale * sqrt(sum);
00254     }
00255 
00256     ret_val = value;
00257     return ret_val;
00258 
00259 /*     End of ZLANSB */
00260 
00261 } /* zlansb_ */


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