slangb.c
Go to the documentation of this file.
00001 /* slangb.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 slangb_(char *norm, integer *n, integer *kl, integer *ku, real *ab, 
00021          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, i__6;
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, k, l;
00032     real sum, scale;
00033     extern logical lsame_(char *, char *);
00034     real value;
00035     extern /* Subroutine */ int slassq_(integer *, real *, integer *, real *, 
00036             real *);
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 /*  SLANGB  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 band matrix  A,  with kl sub-diagonals and ku super-diagonals. */
00054 
00055 /*  Description */
00056 /*  =========== */
00057 
00058 /*  SLANGB returns the value */
00059 
00060 /*     SLANGB = ( 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 SLANGB as described */
00078 /*          above. */
00079 
00080 /*  N       (input) INTEGER */
00081 /*          The order of the matrix A.  N >= 0.  When N = 0, SLANGB is */
00082 /*          set to zero. */
00083 
00084 /*  KL      (input) INTEGER */
00085 /*          The number of sub-diagonals of the matrix A.  KL >= 0. */
00086 
00087 /*  KU      (input) INTEGER */
00088 /*          The number of super-diagonals of the matrix A.  KU >= 0. */
00089 
00090 /*  AB      (input) REAL array, dimension (LDAB,N) */
00091 /*          The band matrix A, stored in rows 1 to KL+KU+1.  The j-th */
00092 /*          column of A is stored in the j-th column of the array AB as */
00093 /*          follows: */
00094 /*          AB(ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(n,j+kl). */
00095 
00096 /*  LDAB    (input) INTEGER */
00097 /*          The leading dimension of the array AB.  LDAB >= KL+KU+1. */
00098 
00099 /*  WORK    (workspace) REAL array, dimension (MAX(1,LWORK)), */
00100 /*          where LWORK >= N when NORM = 'I'; otherwise, WORK is not */
00101 /*          referenced. */
00102 
00103 /* ===================================================================== */
00104 
00105 
00106 /*     .. Parameters .. */
00107 /*     .. */
00108 /*     .. Local Scalars .. */
00109 /*     .. */
00110 /*     .. External Subroutines .. */
00111 /*     .. */
00112 /*     .. External Functions .. */
00113 /*     .. */
00114 /*     .. Intrinsic Functions .. */
00115 /*     .. */
00116 /*     .. Executable Statements .. */
00117 
00118     /* Parameter adjustments */
00119     ab_dim1 = *ldab;
00120     ab_offset = 1 + ab_dim1;
00121     ab -= ab_offset;
00122     --work;
00123 
00124     /* Function Body */
00125     if (*n == 0) {
00126         value = 0.f;
00127     } else if (lsame_(norm, "M")) {
00128 
00129 /*        Find max(abs(A(i,j))). */
00130 
00131         value = 0.f;
00132         i__1 = *n;
00133         for (j = 1; j <= i__1; ++j) {
00134 /* Computing MAX */
00135             i__2 = *ku + 2 - j;
00136 /* Computing MIN */
00137             i__4 = *n + *ku + 1 - j, i__5 = *kl + *ku + 1;
00138             i__3 = min(i__4,i__5);
00139             for (i__ = max(i__2,1); i__ <= i__3; ++i__) {
00140 /* Computing MAX */
00141                 r__2 = value, r__3 = (r__1 = ab[i__ + j * ab_dim1], dabs(r__1)
00142                         );
00143                 value = dmax(r__2,r__3);
00144 /* L10: */
00145             }
00146 /* L20: */
00147         }
00148     } else if (lsame_(norm, "O") || *(unsigned char *)
00149             norm == '1') {
00150 
00151 /*        Find norm1(A). */
00152 
00153         value = 0.f;
00154         i__1 = *n;
00155         for (j = 1; j <= i__1; ++j) {
00156             sum = 0.f;
00157 /* Computing MAX */
00158             i__3 = *ku + 2 - j;
00159 /* Computing MIN */
00160             i__4 = *n + *ku + 1 - j, i__5 = *kl + *ku + 1;
00161             i__2 = min(i__4,i__5);
00162             for (i__ = max(i__3,1); i__ <= i__2; ++i__) {
00163                 sum += (r__1 = ab[i__ + j * ab_dim1], dabs(r__1));
00164 /* L30: */
00165             }
00166             value = dmax(value,sum);
00167 /* L40: */
00168         }
00169     } else if (lsame_(norm, "I")) {
00170 
00171 /*        Find normI(A). */
00172 
00173         i__1 = *n;
00174         for (i__ = 1; i__ <= i__1; ++i__) {
00175             work[i__] = 0.f;
00176 /* L50: */
00177         }
00178         i__1 = *n;
00179         for (j = 1; j <= i__1; ++j) {
00180             k = *ku + 1 - j;
00181 /* Computing MAX */
00182             i__2 = 1, i__3 = j - *ku;
00183 /* Computing MIN */
00184             i__5 = *n, i__6 = j + *kl;
00185             i__4 = min(i__5,i__6);
00186             for (i__ = max(i__2,i__3); i__ <= i__4; ++i__) {
00187                 work[i__] += (r__1 = ab[k + i__ + j * ab_dim1], dabs(r__1));
00188 /* L60: */
00189             }
00190 /* L70: */
00191         }
00192         value = 0.f;
00193         i__1 = *n;
00194         for (i__ = 1; i__ <= i__1; ++i__) {
00195 /* Computing MAX */
00196             r__1 = value, r__2 = work[i__];
00197             value = dmax(r__1,r__2);
00198 /* L80: */
00199         }
00200     } else if (lsame_(norm, "F") || lsame_(norm, "E")) {
00201 
00202 /*        Find normF(A). */
00203 
00204         scale = 0.f;
00205         sum = 1.f;
00206         i__1 = *n;
00207         for (j = 1; j <= i__1; ++j) {
00208 /* Computing MAX */
00209             i__4 = 1, i__2 = j - *ku;
00210             l = max(i__4,i__2);
00211             k = *ku + 1 - j + l;
00212 /* Computing MIN */
00213             i__2 = *n, i__3 = j + *kl;
00214             i__4 = min(i__2,i__3) - l + 1;
00215             slassq_(&i__4, &ab[k + j * ab_dim1], &c__1, &scale, &sum);
00216 /* L90: */
00217         }
00218         value = scale * sqrt(sum);
00219     }
00220 
00221     ret_val = value;
00222     return ret_val;
00223 
00224 /*     End of SLANGB */
00225 
00226 } /* slangb_ */


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