zgelsd.f.html

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

</span><span class="comment">*</span><span class="comment">          where MINMN = MIN( M,N ).
</span><span class="comment">*</span><span class="comment">          On exit, if INFO = 0, IWORK(1) returns the minimum LIWORK.
</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:  the algorithm for computing the SVD failed to converge;
</span><span class="comment">*</span><span class="comment">                if INFO = i, i off-diagonal elements of an intermediate
</span><span class="comment">*</span><span class="comment">                bidiagonal form did not converge to zero.
</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">     Ming Gu and Ren-Cang Li, 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">     Osni Marques, LBNL/NERSC, 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 )
      COMPLEX*16         CZERO
      PARAMETER          ( CZERO = ( 0.0D+0, 0.0D+0 ) )
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Local Scalars ..
</span>      LOGICAL            LQUERY
      INTEGER            IASCL, IBSCL, IE, IL, ITAU, ITAUP, ITAUQ,
     $                   LDWORK, LIWORK, LRWORK, MAXMN, MAXWRK, MINMN,
     $                   MINWRK, MM, MNTHR, NLVL, NRWORK, NWORK, SMLSIZ
      DOUBLE PRECISION   ANRM, BIGNUM, BNRM, EPS, SFMIN, SMLNUM
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Subroutines ..
</span>      EXTERNAL           <a name="DLABAD.164"></a><a href="dlabad.f.html#DLABAD.1">DLABAD</a>, <a name="DLASCL.164"></a><a href="dlascl.f.html#DLASCL.1">DLASCL</a>, <a name="DLASET.164"></a><a href="dlaset.f.html#DLASET.1">DLASET</a>, <a name="XERBLA.164"></a><a href="xerbla.f.html#XERBLA.1">XERBLA</a>, <a name="ZGEBRD.164"></a><a href="zgebrd.f.html#ZGEBRD.1">ZGEBRD</a>, <a name="ZGELQF.164"></a><a href="zgelqf.f.html#ZGELQF.1">ZGELQF</a>,
     $                   <a name="ZGEQRF.165"></a><a href="zgeqrf.f.html#ZGEQRF.1">ZGEQRF</a>, <a name="ZLACPY.165"></a><a href="zlacpy.f.html#ZLACPY.1">ZLACPY</a>, <a name="ZLALSD.165"></a><a href="zlalsd.f.html#ZLALSD.1">ZLALSD</a>, <a name="ZLASCL.165"></a><a href="zlascl.f.html#ZLASCL.1">ZLASCL</a>, <a name="ZLASET.165"></a><a href="zlaset.f.html#ZLASET.1">ZLASET</a>, <a name="ZUNMBR.165"></a><a href="zunmbr.f.html#ZUNMBR.1">ZUNMBR</a>,
     $                   <a name="ZUNMLQ.166"></a><a href="zunmlq.f.html#ZUNMLQ.1">ZUNMLQ</a>, <a name="ZUNMQR.166"></a><a href="zunmqr.f.html#ZUNMQR.1">ZUNMQR</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. External Functions ..
</span>      INTEGER            <a name="ILAENV.169"></a><a href="hfy-index.html#ILAENV">ILAENV</a>
      DOUBLE PRECISION   <a name="DLAMCH.170"></a><a href="dlamch.f.html#DLAMCH.1">DLAMCH</a>, <a name="ZLANGE.170"></a><a href="zlange.f.html#ZLANGE.1">ZLANGE</a>
      EXTERNAL           <a name="ILAENV.171"></a><a href="hfy-index.html#ILAENV">ILAENV</a>, <a name="DLAMCH.171"></a><a href="dlamch.f.html#DLAMCH.1">DLAMCH</a>, <a name="ZLANGE.171"></a><a href="zlange.f.html#ZLANGE.1">ZLANGE</a>
<span class="comment">*</span><span class="comment">     ..
</span><span class="comment">*</span><span class="comment">     .. Intrinsic Functions ..
</span>      INTRINSIC          INT, LOG, MAX, MIN, DBLE
<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 arguments.
</span><span class="comment">*</span><span class="comment">
</span>      INFO = 0
      MINMN = MIN( M, N )
      MAXMN = MAX( M, N )
      LQUERY = ( LWORK.EQ.-1 )
      IF( M.LT.0 ) THEN
         INFO = -1
      ELSE IF( N.LT.0 ) THEN
         INFO = -2
      ELSE IF( NRHS.LT.0 ) THEN
         INFO = -3
      ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
         INFO = -5
      ELSE IF( LDB.LT.MAX( 1, MAXMN ) ) THEN
         INFO = -7
      END IF
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">     Compute workspace.
</span><span class="comment">*</span><span class="comment">     (Note: Comments in the code beginning &quot;Workspace:&quot; describe the
</span><span class="comment">*</span><span class="comment">     minimal amount of workspace needed at that point in the code,
</span><span class="comment">*</span><span class="comment">     as well as the preferred amount for good performance.
</span><span class="comment">*</span><span class="comment">     NB refers to the optimal block size for the immediately
</span><span class="comment">*</span><span class="comment">     following subroutine, as returned by <a name="ILAENV.201"></a><a href="hfy-index.html#ILAENV">ILAENV</a>.)
</span><span class="comment">*</span><span class="comment">
</span>      IF( INFO.EQ.0 ) THEN
         MINWRK = 1
         MAXWRK = 1
         LIWORK = 1
         LRWORK = 1
         IF( MINMN.GT.0 ) THEN
            SMLSIZ = <a name="ILAENV.209"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 9, <span class="string">'<a name="ZGELSD.209"></a><a href="zgelsd.f.html#ZGELSD.1">ZGELSD</a>'</span>, <span class="string">' '</span>, 0, 0, 0, 0 )
            MNTHR = <a name="ILAENV.210"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 6, <span class="string">'<a name="ZGELSD.210"></a><a href="zgelsd.f.html#ZGELSD.1">ZGELSD</a>'</span>, <span class="string">' '</span>, M, N, NRHS, -1 )
            NLVL = MAX( INT( LOG( DBLE( MINMN ) / DBLE( SMLSIZ + 1 ) ) /
     $                  LOG( TWO ) ) + 1, 0 )
            LIWORK = 3*MINMN*NLVL + 11*MINMN
            MM = M
            IF( M.GE.N .AND. M.GE.MNTHR ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Path 1a - overdetermined, with many more rows than
</span><span class="comment">*</span><span class="comment">                        columns.
</span><span class="comment">*</span><span class="comment">
</span>               MM = N
               MAXWRK = MAX( MAXWRK, N*<a name="ILAENV.221"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZGEQRF.221"></a><a href="zgeqrf.f.html#ZGEQRF.1">ZGEQRF</a>'</span>, <span class="string">' '</span>, M, N,
     $                       -1, -1 ) )
               MAXWRK = MAX( MAXWRK, NRHS*<a name="ILAENV.223"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZUNMQR.223"></a><a href="zunmqr.f.html#ZUNMQR.1">ZUNMQR</a>'</span>, <span class="string">'LC'</span>, M,
     $                       NRHS, N, -1 ) )
            END IF
            IF( M.GE.N ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">              Path 1 - overdetermined or exactly determined.
</span><span class="comment">*</span><span class="comment">
</span>               LRWORK = 10*N + 2*N*SMLSIZ + 8*N*NLVL + 3*SMLSIZ*NRHS +
     $                  ( SMLSIZ + 1 )**2
               MAXWRK = MAX( MAXWRK, 2*N + ( MM + N )*<a name="ILAENV.232"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1,
     $                       <span class="string">'<a name="ZGEBRD.233"></a><a href="zgebrd.f.html#ZGEBRD.1">ZGEBRD</a>'</span>, <span class="string">' '</span>, MM, N, -1, -1 ) )
               MAXWRK = MAX( MAXWRK, 2*N + NRHS*<a name="ILAENV.234"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZUNMBR.234"></a><a href="zunmbr.f.html#ZUNMBR.1">ZUNMBR</a>'</span>,
     $                       <span class="string">'QLC'</span>, MM, NRHS, N, -1 ) )
               MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*<a name="ILAENV.236"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1,
     $                       <span class="string">'<a name="ZUNMBR.237"></a><a href="zunmbr.f.html#ZUNMBR.1">ZUNMBR</a>'</span>, <span class="string">'PLN'</span>, N, NRHS, N, -1 ) )
               MAXWRK = MAX( MAXWRK, 2*N + N*NRHS )
               MINWRK = MAX( 2*N + MM, 2*N + N*NRHS )
            END IF
            IF( N.GT.M ) THEN
               LRWORK = 10*M + 2*M*SMLSIZ + 8*M*NLVL + 3*SMLSIZ*NRHS +
     $                  ( SMLSIZ + 1 )**2
               IF( N.GE.MNTHR ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Path 2a - underdetermined, with many more columns
</span><span class="comment">*</span><span class="comment">                           than rows.
</span><span class="comment">*</span><span class="comment">
</span>                  MAXWRK = M + M*<a name="ILAENV.249"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZGELQF.249"></a><a href="zgelqf.f.html#ZGELQF.1">ZGELQF</a>'</span>, <span class="string">' '</span>, M, N, -1,
     $                     -1 )
                  MAXWRK = MAX( MAXWRK, M*M + 4*M + 2*M*<a name="ILAENV.251"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1,
     $                          <span class="string">'<a name="ZGEBRD.252"></a><a href="zgebrd.f.html#ZGEBRD.1">ZGEBRD</a>'</span>, <span class="string">' '</span>, M, M, -1, -1 ) )
                  MAXWRK = MAX( MAXWRK, M*M + 4*M + NRHS*<a name="ILAENV.253"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1,
     $                          <span class="string">'<a name="ZUNMBR.254"></a><a href="zunmbr.f.html#ZUNMBR.1">ZUNMBR</a>'</span>, <span class="string">'QLC'</span>, M, NRHS, M, -1 ) )
                  MAXWRK = MAX( MAXWRK, M*M + 4*M + ( M - 1 )*<a name="ILAENV.255"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1,
     $                          <span class="string">'<a name="ZUNMLQ.256"></a><a href="zunmlq.f.html#ZUNMLQ.1">ZUNMLQ</a>'</span>, <span class="string">'LC'</span>, N, NRHS, M, -1 ) )
                  IF( NRHS.GT.1 ) THEN
                     MAXWRK = MAX( MAXWRK, M*M + M + M*NRHS )
                  ELSE
                     MAXWRK = MAX( MAXWRK, M*M + 2*M )
                  END IF
                  MAXWRK = MAX( MAXWRK, M*M + 4*M + M*NRHS )
               ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment">                 Path 2 - underdetermined.
</span><span class="comment">*</span><span class="comment">
</span>                  MAXWRK = 2*M + ( N + M )*<a name="ILAENV.267"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZGEBRD.267"></a><a href="zgebrd.f.html#ZGEBRD.1">ZGEBRD</a>'</span>, <span class="string">' '</span>, M,
     $                     N, -1, -1 )
                  MAXWRK = MAX( MAXWRK, 2*M + NRHS*<a name="ILAENV.269"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZUNMBR.269"></a><a href="zunmbr.f.html#ZUNMBR.1">ZUNMBR</a>'</span>,
     $                          <span class="string">'QLC'</span>, M, NRHS, M, -1 ) )
                  MAXWRK = MAX( MAXWRK, 2*M + M*<a name="ILAENV.271"></a><a href="hfy-index.html#ILAENV">ILAENV</a>( 1, <span class="string">'<a name="ZUNMBR.271"></a><a href="zunmbr.f.html#ZUNMBR.1">ZUNMBR</a>'</span>,
     $                          <span class="string">'PLN'</span>, N, NRHS, M, -1 ) )
                  MAXWRK = MAX( MAXWRK, 2*M + M*NRHS )
               END IF
               MINWRK = MAX( 2*M + N, 2*M + M*NRHS )
            END IF
         END IF