Merge pull request #15649 from JuliaLang/anj/hypot

Some changes to hypot

Author: andreasnoack

base/docs/helpdb/Base.jl

@@ -314,13 +314,6 @@ Compute a "2d histogram" with respect to the bins delimited by the edges given i
 """
 hist2d!
 
-"""
- hypot(x, y)
-
-Compute the ``\\sqrt{x^2+y^2}`` avoiding overflow and underflow.
-"""
-hypot
-
 """
  airybi(x)
 

base/math.jl

@@ -146,36 +146,51 @@ sqrt(x::Float32) = box(Float32,sqrt_llvm(unbox(Float32,x)))
 sqrt(x::Real) = sqrt(float(x))
 @vectorize_1arg Number sqrt
 
-hypot(x::Real, y::Real) = hypot(promote(float(x), float(y))...)
-function hypot{T<:AbstractFloat}(x::T, y::T)
- x = abs(x)
- y = abs(y)
- if x < y
- x, y = y, x
+"""
+ hypot(x, y)
+
+Compute the hypotenuse ``\\sqrt{x^2+y^2}`` avoiding overflow and underflow.
+"""
+hypot(x::Number, y::Number) = hypot(promote(x, y)...)
+function hypot{T<:Number}(x::T, y::T)
+ ax = abs(x)
+ ay = abs(y)
+ if ax < ay
+ ax, ay = ay, ax
  end
- if x == 0
- r = y/one(x)
+ if ax == 0
+ r = ay / one(ax)
  else
- r = y/x
- if isnan(r)
- isinf(x) && return x
- isinf(y) && return y
- return r
- end
+ r = ay / ax
+ end
+
+ rr = ax * sqrt(1 + r * r)
+
+ # Use type of rr to make sure that return type is the same for
+ # all branches
+ if isnan(r)
+ isinf(ax) && return oftype(rr, Inf)
+ isinf(ay) && return oftype(rr, Inf)
+ return oftype(rr, r)
+ else
+ return rr
  end
- x * sqrt(one(r)+r*r)
 end
+@vectorize_2arg Number hypot
+
+"""
+ hypot(x...)
+
+Compute the hypotenuse ``\\sqrt{\\sum x_i}`` avoiding overflow and underflow.
+"""
+hypot(x::Number...) = vecnorm(x)
 
 atan2(y::Real, x::Real) = atan2(promote(float(y),float(x))...)
 atan2{T<:AbstractFloat}(y::T, x::T) = Base.no_op_err("atan2", T)
 
-for f in (:atan2, :hypot)
- @eval begin
- ($f)(y::Float64, x::Float64) = ccall(($(string(f)),libm), Float64, (Float64, Float64,), y, x)
- ($f)(y::Float32, x::Float32) = ccall(($(string(f,"f")),libm), Float32, (Float32, Float32), y, x)
- @vectorize_2arg Number $f
- end
-end
+atan2(y::Float64, x::Float64) = ccall((:atan2,libm), Float64, (Float64, Float64,), y, x)
+atan2(y::Float32, x::Float32) = ccall((:atan2f,libm), Float32, (Float32, Float32), y, x)
+@vectorize_2arg Number atan2
 
 max{T<:AbstractFloat}(x::T, y::T) = ifelse((y > x) | (signbit(y) < signbit(x)),
  ifelse(isnan(y), x, y), ifelse(isnan(x), y, x))

doc/stdlib/math.rst

@@ -884,7 +884,13 @@ Mathematical Functions
 
  .. Docstring generated from Julia source
 
- Compute the :math:`\sqrt{x^2+y^2}` avoiding overflow and underflow.
+ Compute the hypotenuse :math:`\sqrt{x^2+y^2}` avoiding overflow and underflow.
+
+.. function:: hypot(x...)
+
+ .. Docstring generated from Julia source
+
+ Compute the hypotenuse :math:`\sqrt{\sum x_i}` avoiding overflow and underflow.
 
 .. function:: log(x)
 

test/math.jl

@@ -89,6 +89,8 @@ for T in (Float32, Float64)
  @test_approx_eq exp2(x) exp2(big(x))
  @test_approx_eq expm1(x) expm1(big(x))
  @test_approx_eq hypot(x,y) hypot(big(x),big(y))
+ @test_approx_eq hypot(x,x,y) hypot(hypot(big(x),big(x)),big(y))
+ @test_approx_eq hypot(x,x,y,y) hypot(hypot(big(x),big(x)),hypot(big(y),big(y)))
  @test_approx_eq log(x) log(big(x))
  @test_approx_eq log10(x) log10(big(x))
  @test_approx_eq log1p(x) log1p(big(x))
@@ -177,6 +179,11 @@ for T in (Float32, Float64)
  @test isinf(log1p(-one(T)))
  @test isnan(log1p(convert(T,NaN)))
  @test_throws DomainError log1p(convert(T,-2.0))
+ @test hypot(T(0), T(0)) === T(0)
+ @test hypot(T(Inf), T(Inf)) === T(Inf)
+ @test hypot(T(Inf), T(x)) === T(Inf)
+ @test hypot(T(Inf), T(NaN)) === T(Inf)
+ @test isnan(hypot(T(x), T(NaN)))
 end
 @test_approx_eq exp10(5) exp10(5.0)
 @test_approx_eq exp2(Float16(2.)) exp2(2.)