File:  [local] / rpl / lapack / lapack / dlas2.f
Revision 1.6: download - view: text, annotated - select for diffs - revision graph
Fri Aug 13 21:03:52 2010 UTC (13 years, 9 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_0_19, rpl-4_0_18, HEAD
Patches pour OS/2

    1:       SUBROUTINE DLAS2( F, G, H, SSMIN, SSMAX )
    2: *
    3: *  -- LAPACK auxiliary routine (version 3.2) --
    4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
    5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
    6: *     November 2006
    7: *
    8: *     .. Scalar Arguments ..
    9:       DOUBLE PRECISION   F, G, H, SSMAX, SSMIN
   10: *     ..
   11: *
   12: *  Purpose
   13: *  =======
   14: *
   15: *  DLAS2  computes the singular values of the 2-by-2 matrix
   16: *     [  F   G  ]
   17: *     [  0   H  ].
   18: *  On return, SSMIN is the smaller singular value and SSMAX is the
   19: *  larger singular value.
   20: *
   21: *  Arguments
   22: *  =========
   23: *
   24: *  F       (input) DOUBLE PRECISION
   25: *          The (1,1) element of the 2-by-2 matrix.
   26: *
   27: *  G       (input) DOUBLE PRECISION
   28: *          The (1,2) element of the 2-by-2 matrix.
   29: *
   30: *  H       (input) DOUBLE PRECISION
   31: *          The (2,2) element of the 2-by-2 matrix.
   32: *
   33: *  SSMIN   (output) DOUBLE PRECISION
   34: *          The smaller singular value.
   35: *
   36: *  SSMAX   (output) DOUBLE PRECISION
   37: *          The larger singular value.
   38: *
   39: *  Further Details
   40: *  ===============
   41: *
   42: *  Barring over/underflow, all output quantities are correct to within
   43: *  a few units in the last place (ulps), even in the absence of a guard
   44: *  digit in addition/subtraction.
   45: *
   46: *  In IEEE arithmetic, the code works correctly if one matrix element is
   47: *  infinite.
   48: *
   49: *  Overflow will not occur unless the largest singular value itself
   50: *  overflows, or is within a few ulps of overflow. (On machines with
   51: *  partial overflow, like the Cray, overflow may occur if the largest
   52: *  singular value is within a factor of 2 of overflow.)
   53: *
   54: *  Underflow is harmless if underflow is gradual. Otherwise, results
   55: *  may correspond to a matrix modified by perturbations of size near
   56: *  the underflow threshold.
   57: *
   58: *  ====================================================================
   59: *
   60: *     .. Parameters ..
   61:       DOUBLE PRECISION   ZERO
   62:       PARAMETER          ( ZERO = 0.0D0 )
   63:       DOUBLE PRECISION   ONE
   64:       PARAMETER          ( ONE = 1.0D0 )
   65:       DOUBLE PRECISION   TWO
   66:       PARAMETER          ( TWO = 2.0D0 )
   67: *     ..
   68: *     .. Local Scalars ..
   69:       DOUBLE PRECISION   AS, AT, AU, C, FA, FHMN, FHMX, GA, HA
   70: *     ..
   71: *     .. Intrinsic Functions ..
   72:       INTRINSIC          ABS, MAX, MIN, SQRT
   73: *     ..
   74: *     .. Executable Statements ..
   75: *
   76:       FA = ABS( F )
   77:       GA = ABS( G )
   78:       HA = ABS( H )
   79:       FHMN = MIN( FA, HA )
   80:       FHMX = MAX( FA, HA )
   81:       IF( FHMN.EQ.ZERO ) THEN
   82:          SSMIN = ZERO
   83:          IF( FHMX.EQ.ZERO ) THEN
   84:             SSMAX = GA
   85:          ELSE
   86:             SSMAX = MAX( FHMX, GA )*SQRT( ONE+
   87:      $              ( MIN( FHMX, GA ) / MAX( FHMX, GA ) )**2 )
   88:          END IF
   89:       ELSE
   90:          IF( GA.LT.FHMX ) THEN
   91:             AS = ONE + FHMN / FHMX
   92:             AT = ( FHMX-FHMN ) / FHMX
   93:             AU = ( GA / FHMX )**2
   94:             C = TWO / ( SQRT( AS*AS+AU )+SQRT( AT*AT+AU ) )
   95:             SSMIN = FHMN*C
   96:             SSMAX = FHMX / C
   97:          ELSE
   98:             AU = FHMX / GA
   99:             IF( AU.EQ.ZERO ) THEN
  100: *
  101: *              Avoid possible harmful underflow if exponent range
  102: *              asymmetric (true SSMIN may not underflow even if
  103: *              AU underflows)
  104: *
  105:                SSMIN = ( FHMN*FHMX ) / GA
  106:                SSMAX = GA
  107:             ELSE
  108:                AS = ONE + FHMN / FHMX
  109:                AT = ( FHMX-FHMN ) / FHMX
  110:                C = ONE / ( SQRT( ONE+( AS*AU )**2 )+
  111:      $             SQRT( ONE+( AT*AU )**2 ) )
  112:                SSMIN = ( FHMN*C )*AU
  113:                SSMIN = SSMIN + SSMIN
  114:                SSMAX = GA / ( C+C )
  115:             END IF
  116:          END IF
  117:       END IF
  118:       RETURN
  119: *
  120: *     End of DLAS2
  121: *
  122:       END

CVSweb interface <joel.bertrand@systella.fr>