zppequ.c
Go to the documentation of this file.
00001 /* zppequ.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 zppequ_(char *uplo, integer *n, doublecomplex *ap, 
00017         doublereal *s, doublereal *scond, doublereal *amax, integer *info)
00018 {
00019     /* System generated locals */
00020     integer i__1, i__2;
00021     doublereal d__1, d__2;
00022 
00023     /* Builtin functions */
00024     double sqrt(doublereal);
00025 
00026     /* Local variables */
00027     integer i__, jj;
00028     doublereal smin;
00029     extern logical lsame_(char *, char *);
00030     logical upper;
00031     extern /* Subroutine */ int xerbla_(char *, integer *);
00032 
00033 
00034 /*  -- LAPACK routine (version 3.2) -- */
00035 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00036 /*     November 2006 */
00037 
00038 /*     .. Scalar Arguments .. */
00039 /*     .. */
00040 /*     .. Array Arguments .. */
00041 /*     .. */
00042 
00043 /*  Purpose */
00044 /*  ======= */
00045 
00046 /*  ZPPEQU computes row and column scalings intended to equilibrate a */
00047 /*  Hermitian positive definite matrix A in packed storage and reduce */
00048 /*  its condition number (with respect to the two-norm).  S contains the */
00049 /*  scale factors, S(i)=1/sqrt(A(i,i)), chosen so that the scaled matrix */
00050 /*  B with elements B(i,j)=S(i)*A(i,j)*S(j) has ones on the diagonal. */
00051 /*  This choice of S puts the condition number of B within a factor N of */
00052 /*  the smallest possible condition number over all possible diagonal */
00053 /*  scalings. */
00054 
00055 /*  Arguments */
00056 /*  ========= */
00057 
00058 /*  UPLO    (input) CHARACTER*1 */
00059 /*          = 'U':  Upper triangle of A is stored; */
00060 /*          = 'L':  Lower triangle of A is stored. */
00061 
00062 /*  N       (input) INTEGER */
00063 /*          The order of the matrix A.  N >= 0. */
00064 
00065 /*  AP      (input) COMPLEX*16 array, dimension (N*(N+1)/2) */
00066 /*          The upper or lower triangle of the Hermitian matrix A, packed */
00067 /*          columnwise in a linear array.  The j-th column of A is stored */
00068 /*          in the array AP as follows: */
00069 /*          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; */
00070 /*          if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. */
00071 
00072 /*  S       (output) DOUBLE PRECISION array, dimension (N) */
00073 /*          If INFO = 0, S contains the scale factors for A. */
00074 
00075 /*  SCOND   (output) DOUBLE PRECISION */
00076 /*          If INFO = 0, S contains the ratio of the smallest S(i) to */
00077 /*          the largest S(i).  If SCOND >= 0.1 and AMAX is neither too */
00078 /*          large nor too small, it is not worth scaling by S. */
00079 
00080 /*  AMAX    (output) DOUBLE PRECISION */
00081 /*          Absolute value of largest matrix element.  If AMAX is very */
00082 /*          close to overflow or very close to underflow, the matrix */
00083 /*          should be scaled. */
00084 
00085 /*  INFO    (output) INTEGER */
00086 /*          = 0:  successful exit */
00087 /*          < 0:  if INFO = -i, the i-th argument had an illegal value */
00088 /*          > 0:  if INFO = i, the i-th diagonal element is nonpositive. */
00089 
00090 /*  ===================================================================== */
00091 
00092 /*     .. Parameters .. */
00093 /*     .. */
00094 /*     .. Local Scalars .. */
00095 /*     .. */
00096 /*     .. External Functions .. */
00097 /*     .. */
00098 /*     .. External Subroutines .. */
00099 /*     .. */
00100 /*     .. Intrinsic Functions .. */
00101 /*     .. */
00102 /*     .. Executable Statements .. */
00103 
00104 /*     Test the input parameters. */
00105 
00106     /* Parameter adjustments */
00107     --s;
00108     --ap;
00109 
00110     /* Function Body */
00111     *info = 0;
00112     upper = lsame_(uplo, "U");
00113     if (! upper && ! lsame_(uplo, "L")) {
00114         *info = -1;
00115     } else if (*n < 0) {
00116         *info = -2;
00117     }
00118     if (*info != 0) {
00119         i__1 = -(*info);
00120         xerbla_("ZPPEQU", &i__1);
00121         return 0;
00122     }
00123 
00124 /*     Quick return if possible */
00125 
00126     if (*n == 0) {
00127         *scond = 1.;
00128         *amax = 0.;
00129         return 0;
00130     }
00131 
00132 /*     Initialize SMIN and AMAX. */
00133 
00134     s[1] = ap[1].r;
00135     smin = s[1];
00136     *amax = s[1];
00137 
00138     if (upper) {
00139 
00140 /*        UPLO = 'U':  Upper triangle of A is stored. */
00141 /*        Find the minimum and maximum diagonal elements. */
00142 
00143         jj = 1;
00144         i__1 = *n;
00145         for (i__ = 2; i__ <= i__1; ++i__) {
00146             jj += i__;
00147             i__2 = jj;
00148             s[i__] = ap[i__2].r;
00149 /* Computing MIN */
00150             d__1 = smin, d__2 = s[i__];
00151             smin = min(d__1,d__2);
00152 /* Computing MAX */
00153             d__1 = *amax, d__2 = s[i__];
00154             *amax = max(d__1,d__2);
00155 /* L10: */
00156         }
00157 
00158     } else {
00159 
00160 /*        UPLO = 'L':  Lower triangle of A is stored. */
00161 /*        Find the minimum and maximum diagonal elements. */
00162 
00163         jj = 1;
00164         i__1 = *n;
00165         for (i__ = 2; i__ <= i__1; ++i__) {
00166             jj = jj + *n - i__ + 2;
00167             i__2 = jj;
00168             s[i__] = ap[i__2].r;
00169 /* Computing MIN */
00170             d__1 = smin, d__2 = s[i__];
00171             smin = min(d__1,d__2);
00172 /* Computing MAX */
00173             d__1 = *amax, d__2 = s[i__];
00174             *amax = max(d__1,d__2);
00175 /* L20: */
00176         }
00177     }
00178 
00179     if (smin <= 0.) {
00180 
00181 /*        Find the first non-positive diagonal element and return. */
00182 
00183         i__1 = *n;
00184         for (i__ = 1; i__ <= i__1; ++i__) {
00185             if (s[i__] <= 0.) {
00186                 *info = i__;
00187                 return 0;
00188             }
00189 /* L30: */
00190         }
00191     } else {
00192 
00193 /*        Set the scale factors to the reciprocals */
00194 /*        of the diagonal elements. */
00195 
00196         i__1 = *n;
00197         for (i__ = 1; i__ <= i__1; ++i__) {
00198             s[i__] = 1. / sqrt(s[i__]);
00199 /* L40: */
00200         }
00201 
00202 /*        Compute SCOND = min(S(I)) / max(S(I)) */
00203 
00204         *scond = sqrt(smin) / sqrt(*amax);
00205     }
00206     return 0;
00207 
00208 /*     End of ZPPEQU */
00209 
00210 } /* zppequ_ */


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