[BACK]Return to dladiv.f CVS log [TXT] [DIR] Up to [local] / rpl / lapack / lapack

Diff for /rpl/lapack/lapack/dladiv.f between versions 1.4 and 1.18

version 1.4, 2010/08/06 15:32:25 version 1.18, 2017/06/17 11:06:21
Line 1 Line 1
   *> \brief \b DLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.
   *
   *  =========== DOCUMENTATION ===========
   *
   * Online html documentation available at
   *            http://www.netlib.org/lapack/explore-html/
   *
   *> \htmlonly
   *> Download DLADIV + dependencies
   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dladiv.f">
   *> [TGZ]</a>
   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dladiv.f">
   *> [ZIP]</a>
   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dladiv.f">
   *> [TXT]</a>
   *> \endhtmlonly
   *
   *  Definition:
   *  ===========
   *
   *       SUBROUTINE DLADIV( A, B, C, D, P, Q )
   *
   *       .. Scalar Arguments ..
   *       DOUBLE PRECISION   A, B, C, D, P, Q
   *       ..
   *
   *
   *> \par Purpose:
   *  =============
   *>
   *> \verbatim
   *>
   *> DLADIV performs complex division in  real arithmetic
   *>
   *>                       a + i*b
   *>            p + i*q = ---------
   *>                       c + i*d
   *>
   *> The algorithm is due to Michael Baudin and Robert L. Smith
   *> and can be found in the paper
   *> "A Robust Complex Division in Scilab"
   *> \endverbatim
   *
   *  Arguments:
   *  ==========
   *
   *> \param[in] A
   *> \verbatim
   *>          A is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[in] B
   *> \verbatim
   *>          B is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[in] C
   *> \verbatim
   *>          C is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[in] D
   *> \verbatim
   *>          D is DOUBLE PRECISION
   *>          The scalars a, b, c, and d in the above expression.
   *> \endverbatim
   *>
   *> \param[out] P
   *> \verbatim
   *>          P is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[out] Q
   *> \verbatim
   *>          Q is DOUBLE PRECISION
   *>          The scalars p and q in the above expression.
   *> \endverbatim
   *
   *  Authors:
   *  ========
   *
   *> \author Univ. of Tennessee
   *> \author Univ. of California Berkeley
   *> \author Univ. of Colorado Denver
   *> \author NAG Ltd.
   *
   *> \date January 2013
   *
   *> \ingroup doubleOTHERauxiliary
   *
   *  =====================================================================
       SUBROUTINE DLADIV( A, B, C, D, P, Q )        SUBROUTINE DLADIV( A, B, C, D, P, Q )
 *  *
 *  -- LAPACK auxiliary routine (version 3.2) --  *  -- LAPACK auxiliary routine (version 3.7.0) --
 *  -- 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 2006  *     January 2013
 *  *
 *     .. Scalar Arguments ..  *     .. Scalar Arguments ..
       DOUBLE PRECISION   A, B, C, D, P, Q        DOUBLE PRECISION   A, B, C, D, P, Q
 *     ..  *     ..
 *  *
 *  Purpose  *  =====================================================================
 *  =======  
 *  *
 *  DLADIV performs complex division in  real arithmetic  *     .. Parameters ..
         DOUBLE PRECISION   BS
         PARAMETER          ( BS = 2.0D0 )
         DOUBLE PRECISION   HALF
         PARAMETER          ( HALF = 0.5D0 )
         DOUBLE PRECISION   TWO
         PARAMETER          ( TWO = 2.0D0 )
 *  *
 *                        a + i*b  *     .. Local Scalars ..
 *             p + i*q = ---------        DOUBLE PRECISION   AA, BB, CC, DD, AB, CD, S, OV, UN, BE, EPS
 *                        c + i*d  *     ..
 *  *     .. External Functions ..
 *  The algorithm is due to Robert L. Smith and can be found        DOUBLE PRECISION   DLAMCH
 *  in D. Knuth, The art of Computer Programming, Vol.2, p.195        EXTERNAL           DLAMCH
 *  *     ..
 *  Arguments  *     .. External Subroutines ..
 *  =========        EXTERNAL           DLADIV1
 *  *     ..
 *  A       (input) DOUBLE PRECISION  *     .. Intrinsic Functions ..
 *  B       (input) DOUBLE PRECISION        INTRINSIC          ABS, MAX
 *  C       (input) DOUBLE PRECISION  *     ..
 *  D       (input) DOUBLE PRECISION  *     .. Executable Statements ..
 *          The scalars a, b, c, and d in the above expression.  *
 *        AA = A
 *  P       (output) DOUBLE PRECISION        BB = B
 *  Q       (output) DOUBLE PRECISION        CC = C
 *          The scalars p and q in the above expression.        DD = D
         AB = MAX( ABS(A), ABS(B) )
         CD = MAX( ABS(C), ABS(D) )
         S = 1.0D0
   
         OV = DLAMCH( 'Overflow threshold' )
         UN = DLAMCH( 'Safe minimum' )
         EPS = DLAMCH( 'Epsilon' )
         BE = BS / (EPS*EPS)
   
         IF( AB >= HALF*OV ) THEN
            AA = HALF * AA
            BB = HALF * BB
            S  = TWO * S
         END IF
         IF( CD >= HALF*OV ) THEN
            CC = HALF * CC
            DD = HALF * DD
            S  = HALF * S
         END IF
         IF( AB <= UN*BS/EPS ) THEN
            AA = AA * BE
            BB = BB * BE
            S  = S / BE
         END IF
         IF( CD <= UN*BS/EPS ) THEN
            CC = CC * BE
            DD = DD * BE
            S  = S * BE
         END IF
         IF( ABS( D ).LE.ABS( C ) ) THEN
            CALL DLADIV1(AA, BB, CC, DD, P, Q)
         ELSE
            CALL DLADIV1(BB, AA, DD, CC, P, Q)
            Q = -Q
         END IF
         P = P * S
         Q = Q * S
   *
         RETURN
   *
   *     End of DLADIV
   *
         END
   
   *> \ingroup doubleOTHERauxiliary
   
   
         SUBROUTINE DLADIV1( A, B, C, D, P, Q )
   *
   *  -- LAPACK auxiliary routine (version 3.7.0) --
   *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
   *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
   *     January 2013
   *
   *     .. Scalar Arguments ..
         DOUBLE PRECISION   A, B, C, D, P, Q
   *     ..
 *  *
 *  =====================================================================  *  =====================================================================
 *  *
   *     .. Parameters ..
         DOUBLE PRECISION   ONE
         PARAMETER          ( ONE = 1.0D0 )
   *
 *     .. Local Scalars ..  *     .. Local Scalars ..
       DOUBLE PRECISION   E, F        DOUBLE PRECISION   R, T
 *     ..  *     ..
 *     .. Intrinsic Functions ..  *     .. External Functions ..
       INTRINSIC          ABS        DOUBLE PRECISION   DLADIV2
         EXTERNAL           DLADIV2
   *     ..
   *     .. Executable Statements ..
   *
         R = D / C
         T = ONE / (C + D * R)
         P = DLADIV2(A, B, C, D, R, T)
         A = -A
         Q = DLADIV2(B, A, C, D, R, T)
   *
         RETURN
   *
   *     End of DLADIV1
   *
         END
   
   *> \ingroup doubleOTHERauxiliary
   
         DOUBLE PRECISION FUNCTION DLADIV2( A, B, C, D, R, T )
   *
   *  -- LAPACK auxiliary routine (version 3.7.0) --
   *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
   *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
   *     January 2013
   *
   *     .. Scalar Arguments ..
         DOUBLE PRECISION   A, B, C, D, R, T
   *     ..
   *
   *  =====================================================================
   *
   *     .. Parameters ..
         DOUBLE PRECISION   ZERO
         PARAMETER          ( ZERO = 0.0D0 )
   *
   *     .. Local Scalars ..
         DOUBLE PRECISION   BR
 *     ..  *     ..
 *     .. Executable Statements ..  *     .. Executable Statements ..
 *  *
       IF( ABS( D ).LT.ABS( C ) ) THEN        IF( R.NE.ZERO ) THEN
          E = D / C           BR = B * R
          F = C + D*E           IF( BR.NE.ZERO ) THEN
          P = ( A+B*E ) / F              DLADIV2 = (A + BR) * T
          Q = ( B-A*E ) / F           ELSE
               DLADIV2 = A * T + (B * T) * R
            END IF
       ELSE        ELSE
          E = C / D           DLADIV2 = (A + D * (B / C)) * T
          F = D + C*E  
          P = ( B+A*E ) / F  
          Q = ( -A+B*E ) / F  
       END IF        END IF
 *  *
       RETURN        RETURN
 *  *
 *     End of DLADIV  *     End of DLADIV12
 *  *
       END        END
Removed from v.1.4  
changed lines
  Added in v.1.18

CVSweb interface <joel.bertrand@systella.fr>