Annotation of rpl/lapack/lapack/dladiv.f, revision 1.13
1.11 bertrand 1: *> \brief \b DLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.
1.8 bertrand 2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download DLADIV + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dladiv.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dladiv.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dladiv.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DLADIV( A, B, C, D, P, Q )
22: *
23: * .. Scalar Arguments ..
24: * DOUBLE PRECISION A, B, C, D, P, Q
25: * ..
26: *
27: *
28: *> \par Purpose:
29: * =============
30: *>
31: *> \verbatim
32: *>
33: *> DLADIV performs complex division in real arithmetic
34: *>
35: *> a + i*b
36: *> p + i*q = ---------
37: *> c + i*d
38: *>
1.13 ! bertrand 39: *> The algorithm is due to Michael Baudin and Robert L. Smith
! 40: *> and can be found in the paper
! 41: *> "A Robust Complex Division in Scilab"
1.8 bertrand 42: *> \endverbatim
43: *
44: * Arguments:
45: * ==========
46: *
47: *> \param[in] A
48: *> \verbatim
49: *> A is DOUBLE PRECISION
50: *> \endverbatim
51: *>
52: *> \param[in] B
53: *> \verbatim
54: *> B is DOUBLE PRECISION
55: *> \endverbatim
56: *>
57: *> \param[in] C
58: *> \verbatim
59: *> C is DOUBLE PRECISION
60: *> \endverbatim
61: *>
62: *> \param[in] D
63: *> \verbatim
64: *> D is DOUBLE PRECISION
65: *> The scalars a, b, c, and d in the above expression.
66: *> \endverbatim
67: *>
68: *> \param[out] P
69: *> \verbatim
70: *> P is DOUBLE PRECISION
71: *> \endverbatim
72: *>
73: *> \param[out] Q
74: *> \verbatim
75: *> Q is DOUBLE PRECISION
76: *> The scalars p and q in the above expression.
77: *> \endverbatim
78: *
79: * Authors:
80: * ========
81: *
82: *> \author Univ. of Tennessee
83: *> \author Univ. of California Berkeley
84: *> \author Univ. of Colorado Denver
85: *> \author NAG Ltd.
86: *
1.13 ! bertrand 87: *> \date January 2013
1.8 bertrand 88: *
89: *> \ingroup auxOTHERauxiliary
90: *
91: * =====================================================================
1.1 bertrand 92: SUBROUTINE DLADIV( A, B, C, D, P, Q )
93: *
1.13 ! bertrand 94: * -- LAPACK auxiliary routine (version 3.5.0) --
1.1 bertrand 95: * -- LAPACK is a software package provided by Univ. of Tennessee, --
96: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.13 ! bertrand 97: * January 2013
1.1 bertrand 98: *
99: * .. Scalar Arguments ..
100: DOUBLE PRECISION A, B, C, D, P, Q
101: * ..
102: *
103: * =====================================================================
104: *
1.13 ! bertrand 105: * .. Parameters ..
! 106: DOUBLE PRECISION BS
! 107: PARAMETER ( BS = 2.0D0 )
! 108: DOUBLE PRECISION HALF
! 109: PARAMETER ( HALF = 0.5D0 )
! 110: DOUBLE PRECISION TWO
! 111: PARAMETER ( TWO = 2.0D0 )
! 112: *
1.1 bertrand 113: * .. Local Scalars ..
1.13 ! bertrand 114: DOUBLE PRECISION AA, BB, CC, DD, AB, CD, S, OV, UN, BE, EPS
! 115: * ..
! 116: * .. External Functions ..
! 117: DOUBLE PRECISION DLAMCH
! 118: EXTERNAL DLAMCH
! 119: * ..
! 120: * .. External Subroutines ..
! 121: EXTERNAL DLADIV1
1.1 bertrand 122: * ..
123: * .. Intrinsic Functions ..
1.13 ! bertrand 124: INTRINSIC ABS, MAX
1.1 bertrand 125: * ..
126: * .. Executable Statements ..
127: *
1.13 ! bertrand 128: AA = A
! 129: BB = B
! 130: CC = C
! 131: DD = D
! 132: AB = MAX( ABS(A), ABS(B) )
! 133: CD = MAX( ABS(C), ABS(D) )
! 134: S = 1.0D0
! 135:
! 136: OV = DLAMCH( 'Overflow threshold' )
! 137: UN = DLAMCH( 'Safe minimum' )
! 138: EPS = DLAMCH( 'Epsilon' )
! 139: BE = BS / (EPS*EPS)
! 140:
! 141: IF( AB >= HALF*OV ) THEN
! 142: AA = HALF * AA
! 143: BB = HALF * BB
! 144: S = TWO * S
! 145: END IF
! 146: IF( CD >= HALF*OV ) THEN
! 147: CC = HALF * CC
! 148: DD = HALF * DD
! 149: S = HALF * S
! 150: END IF
! 151: IF( AB <= UN*BS/EPS ) THEN
! 152: AA = AA * BE
! 153: BB = BB * BE
! 154: S = S / BE
! 155: END IF
! 156: IF( CD <= UN*BS/EPS ) THEN
! 157: CC = CC * BE
! 158: DD = DD * BE
! 159: S = S * BE
! 160: END IF
! 161: IF( ABS( D ).LE.ABS( C ) ) THEN
! 162: CALL DLADIV1(AA, BB, CC, DD, P, Q)
1.1 bertrand 163: ELSE
1.13 ! bertrand 164: CALL DLADIV1(BB, AA, DD, CC, P, Q)
! 165: Q = -Q
1.1 bertrand 166: END IF
1.13 ! bertrand 167: P = P * S
! 168: Q = Q * S
1.1 bertrand 169: *
170: RETURN
171: *
172: * End of DLADIV
173: *
174: END
1.13 ! bertrand 175:
! 176:
! 177:
! 178: SUBROUTINE DLADIV1( A, B, C, D, P, Q )
! 179: *
! 180: * -- LAPACK auxiliary routine (version 3.5.0) --
! 181: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 182: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 183: * January 2013
! 184: *
! 185: * .. Scalar Arguments ..
! 186: DOUBLE PRECISION A, B, C, D, P, Q
! 187: * ..
! 188: *
! 189: * =====================================================================
! 190: *
! 191: * .. Parameters ..
! 192: DOUBLE PRECISION ONE
! 193: PARAMETER ( ONE = 1.0D0 )
! 194: *
! 195: * .. Local Scalars ..
! 196: DOUBLE PRECISION R, T
! 197: * ..
! 198: * .. External Functions ..
! 199: DOUBLE PRECISION DLADIV2
! 200: EXTERNAL DLADIV2
! 201: * ..
! 202: * .. Executable Statements ..
! 203: *
! 204: R = D / C
! 205: T = ONE / (C + D * R)
! 206: P = DLADIV2(A, B, C, D, R, T)
! 207: A = -A
! 208: Q = DLADIV2(B, A, C, D, R, T)
! 209: *
! 210: RETURN
! 211: *
! 212: * End of DLADIV1
! 213: *
! 214: END
! 215:
! 216: DOUBLE PRECISION FUNCTION DLADIV2( A, B, C, D, R, T )
! 217: *
! 218: * -- LAPACK auxiliary routine (version 3.5.0) --
! 219: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 220: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 221: * January 2013
! 222: *
! 223: * .. Scalar Arguments ..
! 224: DOUBLE PRECISION A, B, C, D, R, T
! 225: * ..
! 226: *
! 227: * =====================================================================
! 228: *
! 229: * .. Parameters ..
! 230: DOUBLE PRECISION ZERO
! 231: PARAMETER ( ZERO = 0.0D0 )
! 232: *
! 233: * .. Local Scalars ..
! 234: DOUBLE PRECISION BR
! 235: * ..
! 236: * .. Executable Statements ..
! 237: *
! 238: IF( R.NE.ZERO ) THEN
! 239: BR = B * R
! 240: if( BR.NE.ZERO ) THEN
! 241: DLADIV2 = (A + BR) * T
! 242: ELSE
! 243: DLADIV2 = A * T + (B * T) * R
! 244: END IF
! 245: ELSE
! 246: DLADIV2 = (A + D * (B / C)) * T
! 247: END IF
! 248: *
! 249: RETURN
! 250: *
! 251: * End of DLADIV12
! 252: *
! 253: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dladiv.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack