diff options
Diffstat (limited to 'src/fortran/lapack/dladiv.f')
| -rw-r--r-- | src/fortran/lapack/dladiv.f | 62 |
1 files changed, 62 insertions, 0 deletions
diff --git a/src/fortran/lapack/dladiv.f b/src/fortran/lapack/dladiv.f new file mode 100644 index 0000000..b6a74d1 --- /dev/null +++ b/src/fortran/lapack/dladiv.f @@ -0,0 +1,62 @@ + SUBROUTINE DLADIV( A, B, C, D, P, Q ) +* +* -- LAPACK auxiliary routine (version 3.1) -- +* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. +* November 2006 +* +* .. Scalar Arguments .. + DOUBLE PRECISION A, B, C, D, P, Q +* .. +* +* Purpose +* ======= +* +* DLADIV performs complex division in real arithmetic +* +* a + i*b +* p + i*q = --------- +* c + i*d +* +* The algorithm is due to Robert L. Smith and can be found +* in D. Knuth, The art of Computer Programming, Vol.2, p.195 +* +* Arguments +* ========= +* +* A (input) DOUBLE PRECISION +* B (input) DOUBLE PRECISION +* C (input) DOUBLE PRECISION +* D (input) DOUBLE PRECISION +* The scalars a, b, c, and d in the above expression. +* +* P (output) DOUBLE PRECISION +* Q (output) DOUBLE PRECISION +* The scalars p and q in the above expression. +* +* ===================================================================== +* +* .. Local Scalars .. + DOUBLE PRECISION E, F +* .. +* .. Intrinsic Functions .. + INTRINSIC ABS +* .. +* .. Executable Statements .. +* + IF( ABS( D ).LT.ABS( C ) ) THEN + E = D / C + F = C + D*E + P = ( A+B*E ) / F + Q = ( B-A*E ) / F + ELSE + E = C / D + F = D + C*E + P = ( B+A*E ) / F + Q = ( -A+B*E ) / F + END IF +* + RETURN +* +* End of DLADIV +* + END |