File:  [local] / rpl / lapack / lapack / dlasq6.f
Revision 1.18: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:38:59 2023 UTC (9 months, 1 week ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_35, rpl-4_1_34, HEAD
Première mise à jour de lapack et blas.

    1: *> \brief \b DLASQ6 computes one dqd transform in ping-pong form. Used by sbdsqr and sstegr.
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *> \htmlonly
    9: *> Download DLASQ6 + dependencies
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq6.f">
   11: *> [TGZ]</a>
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq6.f">
   13: *> [ZIP]</a>
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq6.f">
   15: *> [TXT]</a>
   16: *> \endhtmlonly
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN,
   22: *                          DNM1, DNM2 )
   23: *
   24: *       .. Scalar Arguments ..
   25: *       INTEGER            I0, N0, PP
   26: *       DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2
   27: *       ..
   28: *       .. Array Arguments ..
   29: *       DOUBLE PRECISION   Z( * )
   30: *       ..
   31: *
   32: *
   33: *> \par Purpose:
   34: *  =============
   35: *>
   36: *> \verbatim
   37: *>
   38: *> DLASQ6 computes one dqd (shift equal to zero) transform in
   39: *> ping-pong form, with protection against underflow and overflow.
   40: *> \endverbatim
   41: *
   42: *  Arguments:
   43: *  ==========
   44: *
   45: *> \param[in] I0
   46: *> \verbatim
   47: *>          I0 is INTEGER
   48: *>        First index.
   49: *> \endverbatim
   50: *>
   51: *> \param[in] N0
   52: *> \verbatim
   53: *>          N0 is INTEGER
   54: *>        Last index.
   55: *> \endverbatim
   56: *>
   57: *> \param[in] Z
   58: *> \verbatim
   59: *>          Z is DOUBLE PRECISION array, dimension ( 4*N )
   60: *>        Z holds the qd array. EMIN is stored in Z(4*N0) to avoid
   61: *>        an extra argument.
   62: *> \endverbatim
   63: *>
   64: *> \param[in] PP
   65: *> \verbatim
   66: *>          PP is INTEGER
   67: *>        PP=0 for ping, PP=1 for pong.
   68: *> \endverbatim
   69: *>
   70: *> \param[out] DMIN
   71: *> \verbatim
   72: *>          DMIN is DOUBLE PRECISION
   73: *>        Minimum value of d.
   74: *> \endverbatim
   75: *>
   76: *> \param[out] DMIN1
   77: *> \verbatim
   78: *>          DMIN1 is DOUBLE PRECISION
   79: *>        Minimum value of d, excluding D( N0 ).
   80: *> \endverbatim
   81: *>
   82: *> \param[out] DMIN2
   83: *> \verbatim
   84: *>          DMIN2 is DOUBLE PRECISION
   85: *>        Minimum value of d, excluding D( N0 ) and D( N0-1 ).
   86: *> \endverbatim
   87: *>
   88: *> \param[out] DN
   89: *> \verbatim
   90: *>          DN is DOUBLE PRECISION
   91: *>        d(N0), the last value of d.
   92: *> \endverbatim
   93: *>
   94: *> \param[out] DNM1
   95: *> \verbatim
   96: *>          DNM1 is DOUBLE PRECISION
   97: *>        d(N0-1).
   98: *> \endverbatim
   99: *>
  100: *> \param[out] DNM2
  101: *> \verbatim
  102: *>          DNM2 is DOUBLE PRECISION
  103: *>        d(N0-2).
  104: *> \endverbatim
  105: *
  106: *  Authors:
  107: *  ========
  108: *
  109: *> \author Univ. of Tennessee
  110: *> \author Univ. of California Berkeley
  111: *> \author Univ. of Colorado Denver
  112: *> \author NAG Ltd.
  113: *
  114: *> \ingroup auxOTHERcomputational
  115: *
  116: *  =====================================================================
  117:       SUBROUTINE DLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN,
  118:      $                   DNM1, DNM2 )
  119: *
  120: *  -- LAPACK computational routine --
  121: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  122: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  123: *
  124: *     .. Scalar Arguments ..
  125:       INTEGER            I0, N0, PP
  126:       DOUBLE PRECISION   DMIN, DMIN1, DMIN2, DN, DNM1, DNM2
  127: *     ..
  128: *     .. Array Arguments ..
  129:       DOUBLE PRECISION   Z( * )
  130: *     ..
  131: *
  132: *  =====================================================================
  133: *
  134: *     .. Parameter ..
  135:       DOUBLE PRECISION   ZERO
  136:       PARAMETER          ( ZERO = 0.0D0 )
  137: *     ..
  138: *     .. Local Scalars ..
  139:       INTEGER            J4, J4P2
  140:       DOUBLE PRECISION   D, EMIN, SAFMIN, TEMP
  141: *     ..
  142: *     .. External Function ..
  143:       DOUBLE PRECISION   DLAMCH
  144:       EXTERNAL           DLAMCH
  145: *     ..
  146: *     .. Intrinsic Functions ..
  147:       INTRINSIC          MIN
  148: *     ..
  149: *     .. Executable Statements ..
  150: *
  151:       IF( ( N0-I0-1 ).LE.0 )
  152:      $   RETURN
  153: *
  154:       SAFMIN = DLAMCH( 'Safe minimum' )
  155:       J4 = 4*I0 + PP - 3
  156:       EMIN = Z( J4+4 )
  157:       D = Z( J4 )
  158:       DMIN = D
  159: *
  160:       IF( PP.EQ.0 ) THEN
  161:          DO 10 J4 = 4*I0, 4*( N0-3 ), 4
  162:             Z( J4-2 ) = D + Z( J4-1 )
  163:             IF( Z( J4-2 ).EQ.ZERO ) THEN
  164:                Z( J4 ) = ZERO
  165:                D = Z( J4+1 )
  166:                DMIN = D
  167:                EMIN = ZERO
  168:             ELSE IF( SAFMIN*Z( J4+1 ).LT.Z( J4-2 ) .AND.
  169:      $               SAFMIN*Z( J4-2 ).LT.Z( J4+1 ) ) THEN
  170:                TEMP = Z( J4+1 ) / Z( J4-2 )
  171:                Z( J4 ) = Z( J4-1 )*TEMP
  172:                D = D*TEMP
  173:             ELSE
  174:                Z( J4 ) = Z( J4+1 )*( Z( J4-1 ) / Z( J4-2 ) )
  175:                D = Z( J4+1 )*( D / Z( J4-2 ) )
  176:             END IF
  177:             DMIN = MIN( DMIN, D )
  178:             EMIN = MIN( EMIN, Z( J4 ) )
  179:    10    CONTINUE
  180:       ELSE
  181:          DO 20 J4 = 4*I0, 4*( N0-3 ), 4
  182:             Z( J4-3 ) = D + Z( J4 )
  183:             IF( Z( J4-3 ).EQ.ZERO ) THEN
  184:                Z( J4-1 ) = ZERO
  185:                D = Z( J4+2 )
  186:                DMIN = D
  187:                EMIN = ZERO
  188:             ELSE IF( SAFMIN*Z( J4+2 ).LT.Z( J4-3 ) .AND.
  189:      $               SAFMIN*Z( J4-3 ).LT.Z( J4+2 ) ) THEN
  190:                TEMP = Z( J4+2 ) / Z( J4-3 )
  191:                Z( J4-1 ) = Z( J4 )*TEMP
  192:                D = D*TEMP
  193:             ELSE
  194:                Z( J4-1 ) = Z( J4+2 )*( Z( J4 ) / Z( J4-3 ) )
  195:                D = Z( J4+2 )*( D / Z( J4-3 ) )
  196:             END IF
  197:             DMIN = MIN( DMIN, D )
  198:             EMIN = MIN( EMIN, Z( J4-1 ) )
  199:    20    CONTINUE
  200:       END IF
  201: *
  202: *     Unroll last two steps.
  203: *
  204:       DNM2 = D
  205:       DMIN2 = DMIN
  206:       J4 = 4*( N0-2 ) - PP
  207:       J4P2 = J4 + 2*PP - 1
  208:       Z( J4-2 ) = DNM2 + Z( J4P2 )
  209:       IF( Z( J4-2 ).EQ.ZERO ) THEN
  210:          Z( J4 ) = ZERO
  211:          DNM1 = Z( J4P2+2 )
  212:          DMIN = DNM1
  213:          EMIN = ZERO
  214:       ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND.
  215:      $         SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN
  216:          TEMP = Z( J4P2+2 ) / Z( J4-2 )
  217:          Z( J4 ) = Z( J4P2 )*TEMP
  218:          DNM1 = DNM2*TEMP
  219:       ELSE
  220:          Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
  221:          DNM1 = Z( J4P2+2 )*( DNM2 / Z( J4-2 ) )
  222:       END IF
  223:       DMIN = MIN( DMIN, DNM1 )
  224: *
  225:       DMIN1 = DMIN
  226:       J4 = J4 + 4
  227:       J4P2 = J4 + 2*PP - 1
  228:       Z( J4-2 ) = DNM1 + Z( J4P2 )
  229:       IF( Z( J4-2 ).EQ.ZERO ) THEN
  230:          Z( J4 ) = ZERO
  231:          DN = Z( J4P2+2 )
  232:          DMIN = DN
  233:          EMIN = ZERO
  234:       ELSE IF( SAFMIN*Z( J4P2+2 ).LT.Z( J4-2 ) .AND.
  235:      $         SAFMIN*Z( J4-2 ).LT.Z( J4P2+2 ) ) THEN
  236:          TEMP = Z( J4P2+2 ) / Z( J4-2 )
  237:          Z( J4 ) = Z( J4P2 )*TEMP
  238:          DN = DNM1*TEMP
  239:       ELSE
  240:          Z( J4 ) = Z( J4P2+2 )*( Z( J4P2 ) / Z( J4-2 ) )
  241:          DN = Z( J4P2+2 )*( DNM1 / Z( J4-2 ) )
  242:       END IF
  243:       DMIN = MIN( DMIN, DN )
  244: *
  245:       Z( J4+2 ) = DN
  246:       Z( 4*N0-PP ) = EMIN
  247:       RETURN
  248: *
  249: *     End of DLASQ6
  250: *
  251:       END

CVSweb interface <joel.bertrand@systella.fr>