dladiv()

subroutine dladiv	(	double precision	a,
		double precision	b,
		double precision	c,
		double precision	d,
		double precision	p,
		double precision	q )

DLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.

Download DLADIV + dependencies [TGZ] [ZIP] [TXT]

Purpose:

!>
!> 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
!> 
!> 

Parameters

[in]	A	!> A is DOUBLE PRECISION !>
[in]	B	!> B is DOUBLE PRECISION !>
[in]	C	!> C is DOUBLE PRECISION !>
[in]	D	!> D is DOUBLE PRECISION !> The scalars a, b, c, and d in the above expression. !>
[out]	P	!> P is DOUBLE PRECISION !>
[out]	Q	!> Q is DOUBLE PRECISION !> The scalars p and q in the above expression. !>

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Definition at line 88 of file dladiv.f.

*
*  -- LAPACK auxiliary routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      DOUBLE PRECISION   A, B, C, D, P, Q
*     ..
*
*  =====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   BS
      parameter( bs = 2.0d0 )
      DOUBLE PRECISION   HALF
      parameter( half = 0.5d0 )
      DOUBLE PRECISION   TWO
      parameter( two = 2.0d0 )
*
*     .. Local Scalars ..
      DOUBLE PRECISION   AA, BB, CC, DD, AB, CD, S, OV, UN, BE, EPS
*     ..
*     .. External Functions ..
      DOUBLE PRECISION   DLAMCH
      EXTERNAL           dlamch
*     ..
*     .. External Subroutines ..
      EXTERNAL           dladiv1
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, max
*     ..
*     .. Executable Statements ..
*
      aa = a
      bb = b
      cc = c
      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
*

Here is the call graph for this function:

Here is the caller graph for this function: