zheevr.f.html

来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 613 行 · 第 1/4 页

</span><span class="comment">*</span><span class="comment">                  ABSTOL + EPS *   max( |a|,|b| ) ,
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">          where EPS is the machine precision.  If ABSTOL is less than
</span><span class="comment">*</span><span class="comment">          or equal to zero, then  EPS*|T|  will be used in its place,
</span><span class="comment">*</span><span class="comment">          where |T| is the 1-norm of the tridiagonal matrix obtained
</span><span class="comment">*</span><span class="comment">          by reducing A to tridiagonal form.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">          See &quot;Computing Small Singular Values of Bidiagonal Matrices
</span><span class="comment">*</span><span class="comment">          with Guaranteed High Relative Accuracy,&quot; by Demmel and
</span><span class="comment">*</span><span class="comment">          Kahan, LAPACK Working Note #3.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">          If high relative accuracy is important, set ABSTOL to
</span><span class="comment">*</span><span class="comment">          <a name="DLAMCH.148"></a><a href="dlamch.f.html#DLAMCH.1">DLAMCH</a>( 'Safe minimum' ).  Doing so will guarantee that
</span><span class="comment">*</span><span class="comment">          eigenvalues are computed to high relative accuracy when
</span><span class="comment">*</span><span class="comment">          possible in future releases.  The current code does not
</span><span class="comment">*</span><span class="comment">          make any guarantees about high relative accuracy, but
</span><span class="comment">*</span><span class="comment">          furutre releases will. See J. Barlow and J. Demmel,
</span><span class="comment">*</span><span class="comment">          &quot;Computing Accurate Eigensystems of Scaled Diagonally
</span><span class="comment">*</span><span class="comment">          Dominant Matrices&quot;, LAPACK Working Note #7, for a discussion
</span><span class="comment">*</span><span class="comment">          of which matrices define their eigenvalues to high relative
</span><span class="comment">*</span><span class="comment">          accuracy.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  M       (output) INTEGER
</span><span class="comment">*</span><span class="comment">          The total number of eigenvalues found.  0 &lt;= M &lt;= N.
</span><span class="comment">*</span><span class="comment">          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  W       (output) DOUBLE PRECISION array, dimension (N)
</span><span class="comment">*</span><span class="comment">          The first M elements contain the selected eigenvalues in
</span><span class="comment">*</span><span class="comment">          ascending order.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Z       (output) COMPLEX*16 array, dimension (LDZ, max(1,M))
</span><span class="comment">*</span><span class="comment">          If JOBZ = 'V', then if INFO = 0, the first M columns of Z
</span><span class="comment">*</span><span class="comment">          contain the orthonormal eigenvectors of the matrix A
</span><span class="comment">*</span><span class="comment">          corresponding to the selected eigenvalues, with the i-th
</span><span class="comment">*</span><span class="comment">          column of Z holding the eigenvector associated with W(i).
</span><span class="comment">*</span><span class="comment">          If JOBZ = 'N', then Z is not referenced.
</span><span class="comment">*</span><span class="comment">          Note: the user must ensure that at least max(1,M) columns are
</span><span class="comment">*</span><span class="comment">          supplied in the array Z; if RANGE = 'V', the exact value of M
</span><span class="comment">*</span><span class="comment">          is not known in advance and an upper bound must be used.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  LDZ     (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The leading dimension of the array Z.  LDZ &gt;= 1, and if
</span><span class="comment">*</span><span class="comment">          JOBZ = 'V', LDZ &gt;= max(1,N).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  ISUPPZ  (output) INTEGER array, dimension ( 2*max(1,M) )
</span><span class="comment">*</span><span class="comment">          The support of the eigenvectors in Z, i.e., the indices
</span><span class="comment">*</span><span class="comment">          indicating the nonzero elements in Z. The i-th eigenvector
</span><span class="comment">*</span><span class="comment">          is nonzero only in elements ISUPPZ( 2*i-1 ) through
</span><span class="comment">*</span><span class="comment">          ISUPPZ( 2*i ).
</span><span class="comment">*</span><span class="comment">********* Implemented only for RANGE = 'A' or 'I' and IU - IL = N - 1
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  WORK    (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
</span><span class="comment">*</span><span class="comment">          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  LWORK   (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The length of the array WORK.  LWORK &gt;= max(1,2*N).
</span><span class="comment">*</span><span class="comment">          For optimal efficiency, LWORK &gt;= (NB+1)*N,
</span><span class="comment">*</span><span class="comment">          where NB is the max of the blocksize for <a name="ZHETRD.193"></a><a href="zhetrd.f.html#ZHETRD.1">ZHETRD</a> and for
</span><span class="comment">*</span><span class="comment">          <a name="ZUNMTR.194"></a><a href="zunmtr.f.html#ZUNMTR.1">ZUNMTR</a> as returned by <a name="ILAENV.194"></a><a href="hfy-index.html#ILAENV">ILAENV</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">          If LWORK = -1, then a workspace query is assumed; the routine
</span><span class="comment">*</span><span class="comment">          only calculates the optimal sizes of the WORK, RWORK and
</span><span class="comment">*</span><span class="comment">          IWORK arrays, returns these values as the first entries of
</span><span class="comment">*</span><span class="comment">          the WORK, RWORK and IWORK arrays, and no error message
</span><span class="comment">*</span><span class="comment">          related to LWORK or LRWORK or LIWORK is issued by <a name="XERBLA.200"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  RWORK   (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LRWORK))
</span><span class="comment">*</span><span class="comment">          On exit, if INFO = 0, RWORK(1) returns the optimal
</span><span class="comment">*</span><span class="comment">          (and minimal) LRWORK.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> LRWORK   (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The length of the array RWORK.  LRWORK &gt;= max(1,24*N).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">          If LRWORK = -1, then a workspace query is assumed; the
</span><span class="comment">*</span><span class="comment">          routine only calculates the optimal sizes of the WORK, RWORK
</span><span class="comment">*</span><span class="comment">          and IWORK arrays, returns these values as the first entries
</span><span class="comment">*</span><span class="comment">          of the WORK, RWORK and IWORK arrays, and no error message
</span><span class="comment">*</span><span class="comment">          related to LWORK or LRWORK or LIWORK is issued by <a name="XERBLA.213"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  IWORK   (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))
</span><span class="comment">*</span><span class="comment">          On exit, if INFO = 0, IWORK(1) returns the optimal
</span><span class="comment">*</span><span class="comment">          (and minimal) LIWORK.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> LIWORK   (input) INTEGER
</span><span class="comment">*</span><span class="comment">          The dimension of the array IWORK.  LIWORK &gt;= max(1,10*N).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">          If LIWORK = -1, then a workspace query is assumed; the
</span><span class="comment">*</span><span class="comment">          routine only calculates the optimal sizes of the WORK, RWORK
</span><span class="comment">*</span><span class="comment">          and IWORK arrays, returns these values as the first entries
</span><span class="comment">*</span><span class="comment">          of the WORK, RWORK and IWORK arrays, and no error message
</span><span class="comment">*</span><span class="comment">          related to LWORK or LRWORK or LIWORK is issued by <a name="XERBLA.226"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  INFO    (output) INTEGER
</span><span class="comment">*</span><span class="comment">          = 0:  successful exit
</span><span class="comment">*</span><span class="comment">          &lt; 0:  if INFO = -i, the i-th argument had an illegal value
</span><span class="comment">*</span><span class="comment">          &gt; 0:  Internal error
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Further Details
</span><span class="comment">*</span><span class="comment">  ===============
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  Based on contributions by
</span><span class="comment">*</span><span class="comment">     Inderjit Dhillon, IBM Almaden, USA
</span><span class="comment">*</span><span class="comment">     Osni Marques, LBNL/NERSC, USA
</span><span class="comment">*</span><span class="comment">     Ken Stanley, Computer Science Division, University of
</span><span class="comment">*</span><span class="comment">       California at Berkeley, USA
</span><span class="comment">*</span><span class="comment">     Jason Riedy, Computer Science Division, University of
</span><span class="comment">*</span><span class="comment">       California at Berkeley, USA
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> =====================================================================
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     .. Parameters ..
</span>      DOUBLE PRECISION   ZERO, ONE, TWO
      PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0 )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      LOGICAL            ALLEIG, INDEIG, LOWER, LQUERY, TEST, VALEIG,
     $                   WANTZ, TRYRAC
      CHARACTER          ORDER
      INTEGER            I, IEEEOK, IINFO, IMAX, INDIBL, INDIFL, INDISP,
     $                   INDIWO, INDRD, INDRDD, INDRE, INDREE, INDRWK,
     $                   INDTAU, INDWK, INDWKN, ISCALE, ITMP1, J, JJ,
     $                   LIWMIN, LLWORK, LLRWORK, LLWRKN, LRWMIN,
     $                   LWKOPT, LWMIN, NB, NSPLIT
      DOUBLE PRECISION   ABSTLL, ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN,
     $                   SIGMA, SMLNUM, TMP1, VLL, VUU
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      LOGICAL            <a name="LSAME.263"></a><a href="lsame.f.html#LSAME.1">LSAME</a>
      INTEGER            <a name="ILAENV.264"></a><a href="hfy-index.html#ILAENV">ILAENV</a>
      DOUBLE PRECISION   <a name="DLAMCH.265"></a><a href="dlamch.f.html#DLAMCH.1">DLAMCH</a>, <a name="ZLANSY.265"></a><a href="zlansy.f.html#ZLANSY.1">ZLANSY</a>
      EXTERNAL           <a name="LSAME.266"></a><a href="lsame.f.html#LSAME.1">LSAME</a>, <a name="ILAENV.266"></a><a href="hfy-index.html#ILAENV">ILAENV</a>, <a name="DLAMCH.266"></a><a href="dlamch.f.html#DLAMCH.1">DLAMCH</a>, <a name="ZLANSY.266"></a><a href="zlansy.f.html#ZLANSY.1">ZLANSY</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           DCOPY, DSCAL, <a name="DSTEBZ.269"></a><a href="dstebz.f.html#DSTEBZ.1">DSTEBZ</a>, <a name="DSTERF.269"></a><a href="dsterf.f.html#DSTERF.1">DSTERF</a>, <a name="XERBLA.269"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>, ZDSCAL,
     $                   <a name="ZHETRD.270"></a><a href="zhetrd.f.html#ZHETRD.1">ZHETRD</a>, <a name="ZSTEMR.270"></a><a href="zstemr.f.html#ZSTEMR.1">ZSTEMR</a>, <a name="ZSTEIN.270"></a><a href="zstein.f.html#ZSTEIN.1">ZSTEIN</a>, ZSWAP, <a name="ZUNMTR.270"></a><a href="zunmtr.f.html#ZUNMTR.1">ZUNMTR</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          DBLE, MAX, MIN, SQRT
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Executable Statements ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Test the input parameters.
</span><span class="comment">*</span><span class="comment">
</span>      IEEEOK = <a name="ILAENV.279"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 10, <span class="string">'<a name="ZHEEVR.279"></a><a href="zheevr.f.html#ZHEEVR.1">ZHEEVR</a>'</span>, <span class="string">'N'</span>, 1, 2, 3, 4 )
<span class="comment">*</span><span class="comment">
</span>      LOWER = <a name="LSAME.281"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'L'</span> )
      WANTZ = <a name="LSAME.282"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( JOBZ, <span class="string">'V'</span> )
      ALLEIG = <a name="LSAME.283"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( RANGE, <span class="string">'A'</span> )
      VALEIG = <a name="LSAME.284"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( RANGE, <span class="string">'V'</span> )
      INDEIG = <a name="LSAME.285"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( RANGE, <span class="string">'I'</span> )
<span class="comment">*</span><span class="comment">
</span>      LQUERY = ( ( LWORK.EQ.-1 ) .OR. ( LRWORK.EQ.-1 ) .OR.
     $         ( LIWORK.EQ.-1 ) )
<span class="comment">*</span><span class="comment">