ssyevx.f.html

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

</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  IWORK   (workspace) INTEGER array, dimension (5*N)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">  IFAIL   (output) INTEGER array, dimension (N)
</span><span class="comment">*</span><span class="comment">          If JOBZ = 'V', then if INFO = 0, the first M elements of
</span><span class="comment">*</span><span class="comment">          IFAIL are zero.  If INFO &gt; 0, then IFAIL contains the
</span><span class="comment">*</span><span class="comment">          indices of the eigenvectors that failed to converge.
</span><span class="comment">*</span><span class="comment">          If JOBZ = 'N', then IFAIL is not referenced.
</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:  if INFO = i, then i eigenvectors failed to converge.
</span><span class="comment">*</span><span class="comment">                Their indices are stored in array IFAIL.
</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>      REAL               ZERO, ONE
      PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+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
      CHARACTER          ORDER
      INTEGER            I, IINFO, IMAX, INDD, INDE, INDEE, INDIBL,
     $                   INDISP, INDIWO, INDTAU, INDWKN, INDWRK, ISCALE,
     $                   ITMP1, J, JJ, LLWORK, LLWRKN, LWKMIN,
     $                   LWKOPT, NB, NSPLIT
      REAL               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.168"></a><a href="lsame.f.html#LSAME.1">LSAME</a>
      INTEGER            <a name="ILAENV.169"></a><a href="hfy-index.html#ILAENV">ILAENV</a>
      REAL               <a name="SLAMCH.170"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>, <a name="SLANSY.170"></a><a href="slansy.f.html#SLANSY.1">SLANSY</a>
      EXTERNAL           <a name="LSAME.171"></a><a href="lsame.f.html#LSAME.1">LSAME</a>, <a name="ILAENV.171"></a><a href="hfy-index.html#ILAENV">ILAENV</a>, <a name="SLAMCH.171"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>, <a name="SLANSY.171"></a><a href="slansy.f.html#SLANSY.1">SLANSY</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           SCOPY, <a name="SLACPY.174"></a><a href="slacpy.f.html#SLACPY.1">SLACPY</a>, <a name="SORGTR.174"></a><a href="sorgtr.f.html#SORGTR.1">SORGTR</a>, <a name="SORMTR.174"></a><a href="sormtr.f.html#SORMTR.1">SORMTR</a>, SSCAL, <a name="SSTEBZ.174"></a><a href="sstebz.f.html#SSTEBZ.1">SSTEBZ</a>,
     $                   <a name="SSTEIN.175"></a><a href="sstein.f.html#SSTEIN.1">SSTEIN</a>, <a name="SSTEQR.175"></a><a href="ssteqr.f.html#SSTEQR.1">SSTEQR</a>, <a name="SSTERF.175"></a><a href="ssterf.f.html#SSTERF.1">SSTERF</a>, SSWAP, <a name="SSYTRD.175"></a><a href="ssytrd.f.html#SSYTRD.1">SSYTRD</a>, <a name="XERBLA.175"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          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>      LOWER = <a name="LSAME.184"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'L'</span> )
      WANTZ = <a name="LSAME.185"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( JOBZ, <span class="string">'V'</span> )
      ALLEIG = <a name="LSAME.186"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( RANGE, <span class="string">'A'</span> )
      VALEIG = <a name="LSAME.187"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( RANGE, <span class="string">'V'</span> )
      INDEIG = <a name="LSAME.188"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( RANGE, <span class="string">'I'</span> )
      LQUERY = ( LWORK.EQ.-1 )
<span class="comment">*</span><span class="comment">
</span>      INFO = 0
      IF( .NOT.( WANTZ .OR. <a name="LSAME.192"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( JOBZ, <span class="string">'N'</span> ) ) ) THEN
         INFO = -1
      ELSE IF( .NOT.( ALLEIG .OR. VALEIG .OR. INDEIG ) ) THEN
         INFO = -2
      ELSE IF( .NOT.( LOWER .OR. <a name="LSAME.196"></a><a href="lsame.f.html#LSAME.1">LSAME</a>( UPLO, <span class="string">'U'</span> ) ) ) THEN
         INFO = -3
      ELSE IF( N.LT.0 ) THEN
         INFO = -4
      ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
         INFO = -6
      ELSE
         IF( VALEIG ) THEN
            IF( N.GT.0 .AND. VU.LE.VL )
     $         INFO = -8
         ELSE IF( INDEIG ) THEN
            IF( IL.LT.1 .OR. IL.GT.MAX( 1, N ) ) THEN
               INFO = -9
            ELSE IF( IU.LT.MIN( N, IL ) .OR. IU.GT.N ) THEN
               INFO = -10
            END IF
         END IF
      END IF
      IF( INFO.EQ.0 ) THEN
         IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
            INFO = -15
         END IF
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( INFO.EQ.0 ) THEN
         IF( N.LE.1 ) THEN
            LWKMIN = 1
            WORK( 1 ) = LWKMIN
         ELSE
            LWKMIN = 8*N
            NB = <a name="ILAENV.226"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="SSYTRD.226"></a><a href="ssytrd.f.html#SSYTRD.1">SSYTRD</a>'</span>, UPLO, N, -1, -1, -1 )
            NB = MAX( NB, <a name="ILAENV.227"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="SORMTR.227"></a><a href="sormtr.f.html#SORMTR.1">SORMTR</a>'</span>, UPLO, N, -1, -1, -1 ) )
            LWKOPT = MAX( LWKMIN, ( NB + 3 )*N )
            WORK( 1 ) = LWKOPT
         END IF
<span class="comment">*</span><span class="comment">
</span>         IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY )
     $      INFO = -17
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( INFO.NE.0 ) THEN
         CALL <a name="XERBLA.237"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>( <span class="string">'<a name="SSYEVX.237"></a><a href="ssyevx.f.html#SSYEVX.1">SSYEVX</a>'</span>, -INFO )
         RETURN
      ELSE IF( LQUERY ) THEN
         RETURN
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Quick return if possible
</span><span class="comment">*</span><span class="comment">
</span>      M = 0
      IF( N.EQ.0 ) THEN
         RETURN
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( N.EQ.1 ) THEN
         IF( ALLEIG .OR. INDEIG ) THEN
            M = 1
            W( 1 ) = A( 1, 1 )
         ELSE
            IF( VL.LT.A( 1, 1 ) .AND. VU.GE.A( 1, 1 ) ) THEN
               M = 1
               W( 1 ) = A( 1, 1 )
            END IF
         END IF
         IF( WANTZ )
     $      Z( 1, 1 ) = ONE
         RETURN
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Get machine constants.
</span><span class="comment">*</span><span class="comment">
</span>      SAFMIN = <a name="SLAMCH.267"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'Safe minimum'</span> )
      EPS = <a name="SLAMCH.268"></a><a href="slamch.f.html#SLAMCH.1">SLAMCH</a>( <span class="string">'Precision'</span> )
      SMLNUM = SAFMIN / EPS
      BIGNUM = ONE / SMLNUM
      RMIN = SQRT( SMLNUM )
      RMAX = MIN( SQRT( BIGNUM ), ONE / SQRT( SQRT( SAFMIN ) ) )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Scale matrix to allowable range, if necessary.
</span><span class="comment">*</span><span class="comment">
</span>      ISCALE = 0
      ABSTLL = ABSTOL
      IF( VALEIG ) THEN
         VLL = VL
         VUU = VU
      END IF
      ANRM = <a name="SLANSY.282"></a><a href="slansy.f.html#SLANSY.1">SLANSY</a>( <span class="string">'M'</span>, UPLO, N, A, LDA, WORK )
      IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
         ISCALE = 1
         SIGMA = RMIN / ANRM
      ELSE IF( ANRM.GT.RMAX ) THEN
         ISCALE = 1