Diff for /rpl/lapack/lapack/dlasq5.f between versions 1.9 and 1.18

version 1.9, 2011/11/21 22:19:35 version 1.18, 2018/05/29 06:55:20
Line 1 Line 1
 *> \brief \b DLASQ5  *> \brief \b DLASQ5 computes one dqds transform in ping-pong form. Used by sbdsqr and sstegr.
 *  *
 *  =========== DOCUMENTATION ===========  *  =========== DOCUMENTATION ===========
 *  *
 * Online html documentation available at   * Online html documentation available at
 *            http://www.netlib.org/lapack/explore-html/   *            http://www.netlib.org/lapack/explore-html/
 *  *
 *> \htmlonly  *> \htmlonly
 *> Download DLASQ5 + dependencies   *> Download DLASQ5 + dependencies
 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq5.f">   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq5.f">
 *> [TGZ]</a>   *> [TGZ]</a>
 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq5.f">   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq5.f">
 *> [ZIP]</a>   *> [ZIP]</a>
 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq5.f">   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq5.f">
 *> [TXT]</a>  *> [TXT]</a>
 *> \endhtmlonly   *> \endhtmlonly
 *  *
 *  Definition:  *  Definition:
 *  ===========  *  ===========
 *  *
 *       SUBROUTINE DLASQ5( I0, N0, Z, PP, TAU, DMIN, DMIN1, DMIN2, DN,  *       SUBROUTINE DLASQ5( I0, N0, Z, PP, TAU, SIGMA, DMIN, DMIN1, DMIN2, DN,
 *                          DNM1, DNM2, IEEE )  *                          DNM1, DNM2, IEEE, EPS )
 *   *
 *       .. Scalar Arguments ..  *       .. Scalar Arguments ..
 *       LOGICAL            IEEE  *       LOGICAL            IEEE
 *       INTEGER            I0, N0, PP  *       INTEGER            I0, N0, PP
 *       DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2, TAU  *       DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2, TAU, SIGMA, EPS
 *       ..  *       ..
 *       .. Array Arguments ..  *       .. Array Arguments ..
 *       DOUBLE PRECISION   Z( * )  *       DOUBLE PRECISION   Z( * )
 *       ..  *       ..
 *    *
 *  *
 *> \par Purpose:  *> \par Purpose:
 *  =============  *  =============
Line 74 Line 74
 *>        This is the shift.  *>        This is the shift.
 *> \endverbatim  *> \endverbatim
 *>  *>
   *> \param[in] SIGMA
   *> \verbatim
   *>          SIGMA is DOUBLE PRECISION
   *>        This is the accumulated shift up to this step.
   *> \endverbatim
   *>
 *> \param[out] DMIN  *> \param[out] DMIN
 *> \verbatim  *> \verbatim
 *>          DMIN is DOUBLE PRECISION  *>          DMIN is DOUBLE PRECISION
Line 115 Line 121
 *>          IEEE is LOGICAL  *>          IEEE is LOGICAL
 *>        Flag for IEEE or non IEEE arithmetic.  *>        Flag for IEEE or non IEEE arithmetic.
 *> \endverbatim  *> \endverbatim
 *  *>
   *> \param[in] EPS
   *> \verbatim
   *>          EPS is DOUBLE PRECISION
   *>        This is the value of epsilon used.
   *> \endverbatim
   *>
 *  Authors:  *  Authors:
 *  ========  *  ========
 *  *
 *> \author Univ. of Tennessee   *> \author Univ. of Tennessee
 *> \author Univ. of California Berkeley   *> \author Univ. of California Berkeley
 *> \author Univ. of Colorado Denver   *> \author Univ. of Colorado Denver
 *> \author NAG Ltd.   *> \author NAG Ltd.
 *  *
 *> \date November 2011  *> \date June 2017
 *  *
 *> \ingroup auxOTHERcomputational  *> \ingroup auxOTHERcomputational
 *  *
 *  =====================================================================  *  =====================================================================
       SUBROUTINE DLASQ5( I0, N0, Z, PP, TAU, DMIN, DMIN1, DMIN2, DN,        SUBROUTINE DLASQ5( I0, N0, Z, PP, TAU, SIGMA, DMIN, DMIN1, DMIN2,
      $                   DNM1, DNM2, IEEE )       $                   DN, DNM1, DNM2, IEEE, EPS )
 *  *
 *  -- LAPACK computational routine (version 3.4.0) --  *  -- LAPACK computational routine (version 3.7.1) --
 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --  *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--  *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 *     November 2011  *     June 2017
 *  *
 *     .. Scalar Arguments ..  *     .. Scalar Arguments ..
       LOGICAL            IEEE        LOGICAL            IEEE
       INTEGER            I0, N0, PP        INTEGER            I0, N0, PP
       DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2, TAU        DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2, TAU,
        $                   SIGMA, EPS
 *     ..  *     ..
 *     .. Array Arguments ..  *     .. Array Arguments ..
       DOUBLE PRECISION   Z( * )        DOUBLE PRECISION   Z( * )
Line 149 Line 162
 *  =====================================================================  *  =====================================================================
 *  *
 *     .. Parameter ..  *     .. Parameter ..
       DOUBLE PRECISION   ZERO        DOUBLE PRECISION   ZERO, HALF
       PARAMETER          ( ZERO = 0.0D0 )        PARAMETER          ( ZERO = 0.0D0, HALF = 0.5 )
 *     ..  *     ..
 *     .. Local Scalars ..  *     .. Local Scalars ..
       INTEGER            J4, J4P2        INTEGER            J4, J4P2
       DOUBLE PRECISION   D, EMIN, TEMP        DOUBLE PRECISION   D, EMIN, TEMP, DTHRESH
 *     ..  *     ..
 *     .. Intrinsic Functions ..  *     .. Intrinsic Functions ..
       INTRINSIC          MIN        INTRINSIC          MIN
Line 164 Line 177
       IF( ( N0-I0-1 ).LE.0 )        IF( ( N0-I0-1 ).LE.0 )
      $   RETURN       $   RETURN
 *  *
         DTHRESH = EPS*(SIGMA+TAU)
         IF( TAU.LT.DTHRESH*HALF ) TAU = ZERO
         IF( TAU.NE.ZERO ) THEN
       J4 = 4*I0 + PP - 3        J4 = 4*I0 + PP - 3
       EMIN = Z( J4+4 )         EMIN = Z( J4+4 )
       D = Z( J4 ) - TAU        D = Z( J4 ) - TAU
       DMIN = D        DMIN = D
       DMIN1 = -Z( J4 )        DMIN1 = -Z( J4 )
Line 176 Line 192
 *  *
          IF( PP.EQ.0 ) THEN           IF( PP.EQ.0 ) THEN
             DO 10 J4 = 4*I0, 4*( N0-3 ), 4              DO 10 J4 = 4*I0, 4*( N0-3 ), 4
                Z( J4-2 ) = D + Z( J4-1 )                  Z( J4-2 ) = D + Z( J4-1 )
                TEMP = Z( J4+1 ) / Z( J4-2 )                 TEMP = Z( J4+1 ) / Z( J4-2 )
                D = D*TEMP - TAU                 D = D*TEMP - TAU
                DMIN = MIN( DMIN, D )                 DMIN = MIN( DMIN, D )
Line 185 Line 201
    10       CONTINUE     10       CONTINUE
          ELSE           ELSE
             DO 20 J4 = 4*I0, 4*( N0-3 ), 4              DO 20 J4 = 4*I0, 4*( N0-3 ), 4
                Z( J4-3 ) = D + Z( J4 )                  Z( J4-3 ) = D + Z( J4 )
                TEMP = Z( J4+2 ) / Z( J4-3 )                 TEMP = Z( J4+2 ) / Z( J4-3 )
                D = D*TEMP - TAU                 D = D*TEMP - TAU
                DMIN = MIN( DMIN, D )                 DMIN = MIN( DMIN, D )
Line 194 Line 210
    20       CONTINUE     20       CONTINUE
          END IF           END IF
 *  *
 *        Unroll last two steps.   *        Unroll last two steps.
 *  *
          DNM2 = D           DNM2 = D
          DMIN2 = DMIN           DMIN2 = DMIN
Line 219 Line 235
 *  *
          IF( PP.EQ.0 ) THEN           IF( PP.EQ.0 ) THEN
             DO 30 J4 = 4*I0, 4*( N0-3 ), 4              DO 30 J4 = 4*I0, 4*( N0-3 ), 4
                Z( J4-2 ) = D + Z( J4-1 )                  Z( J4-2 ) = D + Z( J4-1 )
                IF( D.LT.ZERO ) THEN                 IF( D.LT.ZERO ) THEN
                   RETURN                    RETURN
                ELSE                  ELSE
                   Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) )                    Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) )
                   D = Z( J4+1 )*( D / Z( J4-2 ) ) - TAU                    D = Z( J4+1 )*( D / Z( J4-2 ) ) - TAU
                END IF                 END IF
Line 231 Line 247
    30       CONTINUE     30       CONTINUE
          ELSE           ELSE
             DO 40 J4 = 4*I0, 4*( N0-3 ), 4              DO 40 J4 = 4*I0, 4*( N0-3 ), 4
                Z( J4-3 ) = D + Z( J4 )                  Z( J4-3 ) = D + Z( J4 )
                IF( D.LT.ZERO ) THEN                 IF( D.LT.ZERO ) THEN
                   RETURN                    RETURN
                ELSE                  ELSE
                   Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) )                    Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) )
                   D = Z( J4+2 )*( D / Z( J4-3 ) ) - TAU                    D = Z( J4+2 )*( D / Z( J4-3 ) ) - TAU
                END IF                 END IF
Line 243 Line 259
    40       CONTINUE     40       CONTINUE
          END IF           END IF
 *  *
 *        Unroll last two steps.   *        Unroll last two steps.
 *  *
          DNM2 = D           DNM2 = D
          DMIN2 = DMIN           DMIN2 = DMIN
Line 271 Line 287
          DMIN = MIN( DMIN, DN )           DMIN = MIN( DMIN, DN )
 *  *
       END IF        END IF
         ELSE
   *     This is the version that sets d's to zero if they are small enough
            J4 = 4*I0 + PP - 3
            EMIN = Z( J4+4 )
            D = Z( J4 ) - TAU
            DMIN = D
            DMIN1 = -Z( J4 )
            IF( IEEE ) THEN
   *
   *     Code for IEEE arithmetic.
   *
               IF( PP.EQ.0 ) THEN
                  DO 50 J4 = 4*I0, 4*( N0-3 ), 4
                     Z( J4-2 ) = D + Z( J4-1 )
                     TEMP = Z( J4+1 ) / Z( J4-2 )
                     D = D*TEMP - TAU
                     IF( D.LT.DTHRESH ) D = ZERO
                     DMIN = MIN( DMIN, D )
                     Z( J4 ) = Z( J4-1 )*TEMP
                     EMIN = MIN( Z( J4 ), EMIN )
    50            CONTINUE
               ELSE
                  DO 60 J4 = 4*I0, 4*( N0-3 ), 4
                     Z( J4-3 ) = D + Z( J4 )
                     TEMP = Z( J4+2 ) / Z( J4-3 )
                     D = D*TEMP - TAU
                     IF( D.LT.DTHRESH ) D = ZERO
                     DMIN = MIN( DMIN, D )
                     Z( J4-1 ) = Z( J4 )*TEMP
                     EMIN = MIN( Z( J4-1 ), EMIN )
    60            CONTINUE
               END IF
   *
   *     Unroll last two steps.
   *
               DNM2 = D
               DMIN2 = DMIN
               J4 = 4*( N0-2 ) - PP
               J4P2 = J4 + 2*PP - 1
               Z( J4-2 ) = DNM2 + Z( J4P2 )
               Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
               DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) - TAU
               DMIN = MIN( DMIN, DNM1 )
   *
               DMIN1 = DMIN
               J4 = J4 + 4
               J4P2 = J4 + 2*PP - 1
               Z( J4-2 ) = DNM1 + Z( J4P2 )
               Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
               DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) - TAU
               DMIN = MIN( DMIN, DN )
   *
            ELSE
   *
   *     Code for non IEEE arithmetic.
   *
               IF( PP.EQ.0 ) THEN
                  DO 70 J4 = 4*I0, 4*( N0-3 ), 4
                     Z( J4-2 ) = D + Z( J4-1 )
                     IF( D.LT.ZERO ) THEN
                        RETURN
                     ELSE
                        Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) )
                        D = Z( J4+1 )*( D / Z( J4-2 ) ) - TAU
                     END IF
                     IF( D.LT.DTHRESH) D = ZERO
                     DMIN = MIN( DMIN, D )
                     EMIN = MIN( EMIN, Z( J4 ) )
    70            CONTINUE
               ELSE
                  DO 80 J4 = 4*I0, 4*( N0-3 ), 4
                     Z( J4-3 ) = D + Z( J4 )
                     IF( D.LT.ZERO ) THEN
                        RETURN
                     ELSE
                        Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) )
                        D = Z( J4+2 )*( D / Z( J4-3 ) ) - TAU
                     END IF
                     IF( D.LT.DTHRESH) D = ZERO
                     DMIN = MIN( DMIN, D )
                     EMIN = MIN( EMIN, Z( J4-1 ) )
    80            CONTINUE
               END IF
   *
   *     Unroll last two steps.
   *
               DNM2 = D
               DMIN2 = DMIN
               J4 = 4*( N0-2 ) - PP
               J4P2 = J4 + 2*PP - 1
               Z( J4-2 ) = DNM2 + Z( J4P2 )
               IF( DNM2.LT.ZERO ) THEN
                  RETURN
               ELSE
                  Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
                  DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) ) - TAU
               END IF
               DMIN = MIN( DMIN, DNM1 )
   *
               DMIN1 = DMIN
               J4 = J4 + 4
               J4P2 = J4 + 2*PP - 1
               Z( J4-2 ) = DNM1 + Z( J4P2 )
               IF( DNM1.LT.ZERO ) THEN
                  RETURN
               ELSE
                  Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
                  DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) ) - TAU
               END IF
               DMIN = MIN( DMIN, DN )
   *
            END IF
         END IF
 *  *
       Z( J4+2 ) = DN        Z( J4+2 ) = DN
       Z( 4*N0-PP ) = EMIN        Z( 4*N0-PP ) = EMIN

Removed from v.1.9  
changed lines
  Added in v.1.18


CVSweb interface <joel.bertrand@systella.fr>