cgeevx.f.html

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

HTML
557 行
字号
123
         CALL <a name="SLASCL.354"></a><a href="slascl.f.html#SLASCL.1">SLASCL</a>( <span class="string">'G'</span>, 0, 0, CSCALE, ANRM, 1, 1, DUM, 1, IERR )
         ABNRM = DUM( 1 )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Reduce to upper Hessenberg form
</span><span class="comment">*</span><span class="comment">     (CWorkspace: need 2*N, prefer N+N*NB)
</span><span class="comment">*</span><span class="comment">     (RWorkspace: none)
</span><span class="comment">*</span><span class="comment">
</span>      ITAU = 1
      IWRK = ITAU + N
      CALL <a name="CGEHRD.364"></a><a href="cgehrd.f.html#CGEHRD.1">CGEHRD</a>( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ),
     $             LWORK-IWRK+1, IERR )
<span class="comment">*</span><span class="comment">
</span>      IF( WANTVL ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Want left eigenvectors
</span><span class="comment">*</span><span class="comment">        Copy Householder vectors to VL
</span><span class="comment">*</span><span class="comment">
</span>         SIDE = <span class="string">'L'</span>
         CALL <a name="CLACPY.373"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'L'</span>, N, N, A, LDA, VL, LDVL )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Generate unitary matrix in VL
</span><span class="comment">*</span><span class="comment">        (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
</span><span class="comment">*</span><span class="comment">        (RWorkspace: none)
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="CUNGHR.379"></a><a href="cunghr.f.html#CUNGHR.1">CUNGHR</a>( N, ILO, IHI, VL, LDVL, WORK( ITAU ), WORK( IWRK ),
     $                LWORK-IWRK+1, IERR )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Perform QR iteration, accumulating Schur vectors in VL
</span><span class="comment">*</span><span class="comment">        (CWorkspace: need 1, prefer HSWORK (see comments) )
</span><span class="comment">*</span><span class="comment">        (RWorkspace: none)
</span><span class="comment">*</span><span class="comment">
</span>         IWRK = ITAU
         CALL <a name="CHSEQR.387"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a>( <span class="string">'S'</span>, <span class="string">'V'</span>, N, ILO, IHI, A, LDA, W, VL, LDVL,
     $                WORK( IWRK ), LWORK-IWRK+1, INFO )
<span class="comment">*</span><span class="comment">
</span>         IF( WANTVR ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">           Want left and right eigenvectors
</span><span class="comment">*</span><span class="comment">           Copy Schur vectors to VR
</span><span class="comment">*</span><span class="comment">
</span>            SIDE = <span class="string">'B'</span>
            CALL <a name="CLACPY.396"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'F'</span>, N, N, VL, LDVL, VR, LDVR )
         END IF
<span class="comment">*</span><span class="comment">
</span>      ELSE IF( WANTVR ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Want right eigenvectors
</span><span class="comment">*</span><span class="comment">        Copy Householder vectors to VR
</span><span class="comment">*</span><span class="comment">
</span>         SIDE = <span class="string">'R'</span>
         CALL <a name="CLACPY.405"></a><a href="clacpy.f.html#CLACPY.1">CLACPY</a>( <span class="string">'L'</span>, N, N, A, LDA, VR, LDVR )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Generate unitary matrix in VR
</span><span class="comment">*</span><span class="comment">        (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
</span><span class="comment">*</span><span class="comment">        (RWorkspace: none)
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="CUNGHR.411"></a><a href="cunghr.f.html#CUNGHR.1">CUNGHR</a>( N, ILO, IHI, VR, LDVR, WORK( ITAU ), WORK( IWRK ),
     $                LWORK-IWRK+1, IERR )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Perform QR iteration, accumulating Schur vectors in VR
</span><span class="comment">*</span><span class="comment">        (CWorkspace: need 1, prefer HSWORK (see comments) )
</span><span class="comment">*</span><span class="comment">        (RWorkspace: none)
</span><span class="comment">*</span><span class="comment">
</span>         IWRK = ITAU
         CALL <a name="CHSEQR.419"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a>( <span class="string">'S'</span>, <span class="string">'V'</span>, N, ILO, IHI, A, LDA, W, VR, LDVR,
     $                WORK( IWRK ), LWORK-IWRK+1, INFO )
<span class="comment">*</span><span class="comment">
</span>      ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Compute eigenvalues only
</span><span class="comment">*</span><span class="comment">        If condition numbers desired, compute Schur form
</span><span class="comment">*</span><span class="comment">
</span>         IF( WNTSNN ) THEN
            JOB = <span class="string">'E'</span>
         ELSE
            JOB = <span class="string">'S'</span>
         END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        (CWorkspace: need 1, prefer HSWORK (see comments) )
</span><span class="comment">*</span><span class="comment">        (RWorkspace: none)
</span><span class="comment">*</span><span class="comment">
</span>         IWRK = ITAU
         CALL <a name="CHSEQR.437"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a>( JOB, <span class="string">'N'</span>, N, ILO, IHI, A, LDA, W, VR, LDVR,
     $                WORK( IWRK ), LWORK-IWRK+1, INFO )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     If INFO &gt; 0 from <a name="CHSEQR.441"></a><a href="chseqr.f.html#CHSEQR.1">CHSEQR</a>, then quit
</span><span class="comment">*</span><span class="comment">
</span>      IF( INFO.GT.0 )
     $   GO TO 50
<span class="comment">*</span><span class="comment">
</span>      IF( WANTVL .OR. WANTVR ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Compute left and/or right eigenvectors
</span><span class="comment">*</span><span class="comment">        (CWorkspace: need 2*N)
</span><span class="comment">*</span><span class="comment">        (RWorkspace: need N)
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="CTREVC.452"></a><a href="ctrevc.f.html#CTREVC.1">CTREVC</a>( SIDE, <span class="string">'B'</span>, SELECT, N, A, LDA, VL, LDVL, VR, LDVR,
     $                N, NOUT, WORK( IWRK ), RWORK, IERR )
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Compute condition numbers if desired
</span><span class="comment">*</span><span class="comment">     (CWorkspace: need N*N+2*N unless SENSE = 'E')
</span><span class="comment">*</span><span class="comment">     (RWorkspace: need 2*N unless SENSE = 'E')
</span><span class="comment">*</span><span class="comment">
</span>      IF( .NOT.WNTSNN ) THEN
         CALL <a name="CTRSNA.461"></a><a href="ctrsna.f.html#CTRSNA.1">CTRSNA</a>( SENSE, <span class="string">'A'</span>, SELECT, N, A, LDA, VL, LDVL, VR, LDVR,
     $                RCONDE, RCONDV, N, NOUT, WORK( IWRK ), N, RWORK,
     $                ICOND )
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( WANTVL ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Undo balancing of left eigenvectors
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="CGEBAK.470"></a><a href="cgebak.f.html#CGEBAK.1">CGEBAK</a>( BALANC, <span class="string">'L'</span>, N, ILO, IHI, SCALE, N, VL, LDVL,
     $                IERR )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Normalize left eigenvectors and make largest component real
</span><span class="comment">*</span><span class="comment">
</span>         DO 20 I = 1, N
            SCL = ONE / SCNRM2( N, VL( 1, I ), 1 )
            CALL CSSCAL( N, SCL, VL( 1, I ), 1 )
            DO 10 K = 1, N
               RWORK( K ) = REAL( VL( K, I ) )**2 +
     $                      AIMAG( VL( K, I ) )**2
   10       CONTINUE
            K = ISAMAX( N, RWORK, 1 )
            TMP = CONJG( VL( K, I ) ) / SQRT( RWORK( K ) )
            CALL CSCAL( N, TMP, VL( 1, I ), 1 )
            VL( K, I ) = CMPLX( REAL( VL( K, I ) ), ZERO )
   20    CONTINUE
      END IF
<span class="comment">*</span><span class="comment">
</span>      IF( WANTVR ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Undo balancing of right eigenvectors
</span><span class="comment">*</span><span class="comment">
</span>         CALL <a name="CGEBAK.493"></a><a href="cgebak.f.html#CGEBAK.1">CGEBAK</a>( BALANC, <span class="string">'R'</span>, N, ILO, IHI, SCALE, N, VR, LDVR,
     $                IERR )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">        Normalize right eigenvectors and make largest component real
</span><span class="comment">*</span><span class="comment">
</span>         DO 40 I = 1, N
            SCL = ONE / SCNRM2( N, VR( 1, I ), 1 )
            CALL CSSCAL( N, SCL, VR( 1, I ), 1 )
            DO 30 K = 1, N
               RWORK( K ) = REAL( VR( K, I ) )**2 +
     $                      AIMAG( VR( K, I ) )**2
   30       CONTINUE
            K = ISAMAX( N, RWORK, 1 )
            TMP = CONJG( VR( K, I ) ) / SQRT( RWORK( K ) )
            CALL CSCAL( N, TMP, VR( 1, I ), 1 )
            VR( K, I ) = CMPLX( REAL( VR( K, I ) ), ZERO )
   40    CONTINUE
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Undo scaling if necessary
</span><span class="comment">*</span><span class="comment">
</span>   50 CONTINUE
      IF( SCALEA ) THEN
         CALL <a name="CLASCL.516"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>( <span class="string">'G'</span>, 0, 0, CSCALE, ANRM, N-INFO, 1, W( INFO+1 ),
     $                MAX( N-INFO, 1 ), IERR )
         IF( INFO.EQ.0 ) THEN
            IF( ( WNTSNV .OR. WNTSNB ) .AND. ICOND.EQ.0 )
     $         CALL <a name="SLASCL.520"></a><a href="slascl.f.html#SLASCL.1">SLASCL</a>( <span class="string">'G'</span>, 0, 0, CSCALE, ANRM, N, 1, RCONDV, N,
     $                      IERR )
         ELSE
            CALL <a name="CLASCL.523"></a><a href="clascl.f.html#CLASCL.1">CLASCL</a>( <span class="string">'G'</span>, 0, 0, CSCALE, ANRM, ILO-1, 1, W, N, IERR )
         END IF
      END IF
<span class="comment">*</span><span class="comment">
</span>      WORK( 1 ) = MAXWRK
      RETURN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     End of <a name="CGEEVX.530"></a><a href="cgeevx.f.html#CGEEVX.1">CGEEVX</a>
</span><span class="comment">*</span><span class="comment">
</span>      END

</pre>

 </body>
</html>
123

cgeevx.f.html - 源码说明

本页面展示了「famous linear algebra library (LAPACK) ports to windows」中的 cgeevx.f.html 源码文件，采用 HTML 编程语言编写，共 557 行代码。您可以在线阅读完整代码内容，也可以返回资源详情页下载完整源码包进行本地学习和开发。

虫虫下载站收录了大量与algebra相关的技术资源，包括源代码、技术文档、电路图等，是电子工程师和嵌入式开发者的专业学习平台。

⌨️ 快捷键说明

复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?