diff --git a/src/math/cmplx/asin.go b/src/math/cmplx/asin.go index 61880a257d..062f324ce2 100644 --- a/src/math/cmplx/asin.go +++ b/src/math/cmplx/asin.go @@ -49,11 +49,8 @@ import "math" // Asin returns the inverse sine of x. func Asin(x complex128) complex128 { - if imag(x) == 0 { - if math.Abs(real(x)) > 1 { - return complex(math.Pi/2, 0) // DOMAIN error - } - return complex(math.Asin(real(x)), 0) + if imag(x) == 0 && math.Abs(real(x)) <= 1 { + return complex(math.Asin(real(x)), imag(x)) } ct := complex(-imag(x), real(x)) // i * x xx := x * x @@ -65,12 +62,8 @@ func Asin(x complex128) complex128 { // Asinh returns the inverse hyperbolic sine of x. func Asinh(x complex128) complex128 { - // TODO check range - if imag(x) == 0 { - if math.Abs(real(x)) > 1 { - return complex(math.Pi/2, 0) // DOMAIN error - } - return complex(math.Asinh(real(x)), 0) + if imag(x) == 0 && math.Abs(real(x)) <= 1 { + return complex(math.Asinh(real(x)), imag(x)) } xx := x * x x1 := complex(1+real(xx), imag(xx)) // 1 + x*x @@ -140,10 +133,6 @@ func Acosh(x complex128) complex128 { // Atan returns the inverse tangent of x. func Atan(x complex128) complex128 { - if real(x) == 0 && imag(x) > 1 { - return NaN() - } - x2 := real(x) * real(x) a := 1 - x2 - imag(x)*imag(x) if a == 0 { diff --git a/src/math/cmplx/cmath_test.go b/src/math/cmplx/cmath_test.go index 7a5c485a31..8d705622fd 100644 --- a/src/math/cmplx/cmath_test.go +++ b/src/math/cmplx/cmath_test.go @@ -435,6 +435,24 @@ var tanhSC = []complex128{ NaN(), } +// branch cut continuity checks +// points on each axis at |z| > 1 are checked for one-sided continuity from both the positive and negative side +// all possible branch cuts for the elementary functions are at one of these points + +var zero = 0.0 +var eps = 1.0 / (1 << 53) + +var branchPoints = [][2]complex128{ + {complex(2.0, zero), complex(2.0, eps)}, + {complex(2.0, -zero), complex(2.0, -eps)}, + {complex(-2.0, zero), complex(-2.0, eps)}, + {complex(-2.0, -zero), complex(-2.0, -eps)}, + {complex(zero, 2.0), complex(eps, 2.0)}, + {complex(-zero, 2.0), complex(-eps, 2.0)}, + {complex(zero, -2.0), complex(eps, -2.0)}, + {complex(-zero, -2.0), complex(-eps, -2.0)}, +} + // functions borrowed from pkg/math/all_test.go func tolerance(a, b, e float64) bool { d := a - b @@ -508,6 +526,11 @@ func TestAcos(t *testing.T) { t.Errorf("Acos(%g) = %g, want %g", vcAcosSC[i], f, acosSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Acos(pt[0]), Acos(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Acos(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestAcosh(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -520,6 +543,11 @@ func TestAcosh(t *testing.T) { t.Errorf("Acosh(%g) = %g, want %g", vcAcoshSC[i], f, acoshSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Acosh(pt[0]), Acosh(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Acosh(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestAsin(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -532,6 +560,11 @@ func TestAsin(t *testing.T) { t.Errorf("Asin(%g) = %g, want %g", vcAsinSC[i], f, asinSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Asin(pt[0]), Asin(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Asin(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestAsinh(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -544,6 +577,11 @@ func TestAsinh(t *testing.T) { t.Errorf("Asinh(%g) = %g, want %g", vcAsinhSC[i], f, asinhSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Asinh(pt[0]), Asinh(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Asinh(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestAtan(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -556,6 +594,11 @@ func TestAtan(t *testing.T) { t.Errorf("Atan(%g) = %g, want %g", vcAtanSC[i], f, atanSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Atan(pt[0]), Atan(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Atan(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestAtanh(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -568,6 +611,11 @@ func TestAtanh(t *testing.T) { t.Errorf("Atanh(%g) = %g, want %g", vcAtanhSC[i], f, atanhSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Atanh(pt[0]), Atanh(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Atanh(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestConj(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -635,6 +683,11 @@ func TestLog(t *testing.T) { t.Errorf("Log(%g) = %g, want %g", vcLogSC[i], f, logSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Log(pt[0]), Log(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Log(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestLog10(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -685,6 +738,11 @@ func TestPow(t *testing.T) { t.Errorf("Pow(%g, %g) = %g, want %g", vcPowSC[i][0], vcPowSC[i][0], f, powSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Pow(pt[0], 0.1), Pow(pt[1], 0.1); !cVeryclose(f0, f1) { + t.Errorf("Pow(%g, 0.1) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestRect(t *testing.T) { for i := 0; i < len(vc); i++ { @@ -733,6 +791,11 @@ func TestSqrt(t *testing.T) { t.Errorf("Sqrt(%g) = %g, want %g", vcSqrtSC[i], f, sqrtSC[i]) } } + for _, pt := range branchPoints { + if f0, f1 := Sqrt(pt[0]), Sqrt(pt[1]); !cVeryclose(f0, f1) { + t.Errorf("Sqrt(%g) not continuous, got %g want %g", pt[0], f0, f1) + } + } } func TestTan(t *testing.T) { for i := 0; i < len(vc); i++ { diff --git a/src/math/cmplx/sqrt.go b/src/math/cmplx/sqrt.go index 72f81e907c..0fbdcdedd3 100644 --- a/src/math/cmplx/sqrt.go +++ b/src/math/cmplx/sqrt.go @@ -57,13 +57,14 @@ import "math" // The result r is chosen so that real(r) ≥ 0 and imag(r) has the same sign as imag(x). func Sqrt(x complex128) complex128 { if imag(x) == 0 { + // Ensure that imag(r) has the same sign as imag(x) for imag(x) == signed zero. if real(x) == 0 { - return complex(0, 0) + return complex(0, imag(x)) } if real(x) < 0 { - return complex(0, math.Sqrt(-real(x))) + return complex(0, math.Copysign(math.Sqrt(-real(x)), imag(x))) } - return complex(math.Sqrt(real(x)), 0) + return complex(math.Sqrt(real(x)), imag(x)) } if real(x) == 0 { if imag(x) < 0 {