Mercurial > projects > ldc
diff druntime/src/compiler/dmd/complex.c @ 759:d3eb054172f9
Added copy of druntime from DMD 2.020 modified for LDC.
author | Tomas Lindquist Olsen <tomas.l.olsen@gmail.com> |
---|---|
date | Tue, 11 Nov 2008 01:52:37 +0100 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/druntime/src/compiler/dmd/complex.c Tue Nov 11 01:52:37 2008 +0100 @@ -0,0 +1,107 @@ +/* + * Placed into the public domain. + * Written by Walter Bright + * www.digitalmars.com + */ + + +#include <math.h> + +typedef struct Complex +{ + long double re; + long double im; +} Complex; + +Complex _complex_div(Complex x, Complex y) +{ + Complex q; + long double r; + long double den; + + if (fabs(y.re) < fabs(y.im)) + { + r = y.re / y.im; + den = y.im + r * y.re; + q.re = (x.re * r + x.im) / den; + q.im = (x.im * r - x.re) / den; + } + else + { + r = y.im / y.re; + den = y.re + r * y.im; + q.re = (x.re + r * x.im) / den; + q.im = (x.im - r * x.re) / den; + } + return q; +} + +Complex _complex_mul(Complex x, Complex y) +{ + Complex p; + + p.re = x.re * y.re - x.im * y.im; + p.im = x.im * y.re + x.re * y.im; + return p; +} + +long double _complex_abs(Complex z) +{ + long double x,y,ans,temp; + + x = fabs(z.re); + y = fabs(z.im); + if (x == 0) + ans = y; + else if (y == 0) + ans = x; + else if (x > y) + { + temp = y / x; + ans = x * sqrt(1 + temp * temp); + } + else + { + temp = x / y; + ans = y * sqrt(1 + temp * temp); + } + return ans; +} + +Complex _complex_sqrt(Complex z) +{ + Complex c; + long double x,y,w,r; + + if (z.re == 0 && z.im == 0) + { + c.re = 0; + c.im = 0; + } + else + { + x = fabs(z.re); + y = fabs(z.im); + if (x >= y) + { + r = y / x; + w = sqrt(x) * sqrt(0.5 * (1 + sqrt(1 + r * r))); + } + else + { + r = x / y; + w = sqrt(y) * sqrt(0.5 * (r + sqrt(1 + r * r))); + } + if (z.re >= 0) + { + c.re = w; + c.im = z.im / (w + w); + } + else + { + c.im = (z.im >= 0) ? w : -w; + c.re = z.im / (c.im + c.im); + } + } + return c; +}