Permalink
Browse files

Fix complex atanh/acos/asin

Fix infinity?(complex)
  • Loading branch information...
Jeremy Clifford
Jeremy Clifford committed Apr 21, 2012
1 parent 291aa4b commit 6a1fb7132cf363e2418848fc78a5f72141dcbe73
@@ -313,7 +313,7 @@ inline overload nan?(x:T) = nan?(real(x)) or nan?(imag(x));
inline overload finite?(x:T) = floatExponent(abs(x)) != exponentMask((ComplexRealType(T)));
[T when Complex?(T)]
-inline overload infinity?(x:T) = abs(x) == infinity(T);
+inline overload infinity?(x:T) = abs(x) == infinity(ComplexBaseType(T));
[T when Complex?(T)]
inline overload negativeInfinity?(x:T) = negativeInfinity?(real(x));
@@ -67,7 +67,9 @@ acos(x:Double) {
// Return cosine of complex double value.
-overload acos(x:ComplexDouble) --> res:ComplexDouble {
+[C when Complex?(C)]
+overload acos(x:C) --> res:C {
+ alias T = ComplexBaseType(C);
var y = asin (x);
- res <-- M_PI_2 - real(y) - imag(y);
+ res <-- Complex(T(M_PI_2) - real(y), -imagValue(y));
}
@@ -77,16 +77,16 @@ overload asin(x:C) --> res:C {
alias T = ComplexBaseType(C);
- if (not finite?(x) or nan?(x))
+ if (nan?(real(x)) or nan?(imagValue(x)))
if (real(x) == T(0))
res <-- x;
else if (infinity?(real(x)) or infinity?(imagValue(x)))
- res <-- Complex(nan(T), copysign(infinity(T),imagValue(x)));
+ res <-- Complex(nan(T), copysign(infinity(T),imagValue(x)));
else
- res <-- nan(C);
+ res <-- nan(C);
else {
- var y = asinh(x);
- res <-- Complex(imagValue(y),real(y));
+ var y = asinh(Complex(-imagValue(x), real(x)));
+ res <-- Complex(imagValue(y), -real(y));
}
}
@@ -3,8 +3,7 @@
import math.native.kernel.ieee754.*;
import math.native.log1p.(log1p);
import math.native.fpclassify.(fpclassify,FP_ZERO,FP_INFINITE);
-import math.native.log.(log);
-import math.native.atan2.(atan2);
+import math.native.atan.(atan);
import math.native.copysign.(copysign);
import math.native.constants.(M_PI_2);
import numbers.floats.(nan);
@@ -28,7 +27,8 @@ overload atanh(a:Double) {
if(ix<CV3) { // x < 0.5
t = x+x;
t = 0.5*log1p(t+t*x/(ONE-x));
- } else
+ }
+ else
t = 0.5*log1p((x+x)/(ONE-x));
if(GET_IEEE_HIGH(xdb)>=0) return t;
else return -t;
@@ -55,12 +55,7 @@ overload atanh(x:ComplexDouble) --> res:ComplexDouble {
res <-- x;
}
else {
- var i2 = imagValue(x) * imagValue(x);
- var num = 1.0 + real(x);
- num = i2 + num * num;
- var den = 1.0 - real(x);
- den = i2 + den * den;
- var den2 = 1.0 - real(x) * real(x) - i2;
- res <-- Complex(0.25 * (log(num) - log(den)), 0.5 * atan2(2.0 * imagValue(x), den2));
- }
+ var y = atan(Complex(-imagValue(x), real(x)));
+ res <-- Complex(imagValue(y), -real(y));
+ }
}

0 comments on commit 6a1fb71

Please sign in to comment.