fix xexpy(0, 1000) (#80)

* fix `xexpy(0, 1000)`

* add xexpy tests

* fixup

* update doc

* Update src/basicfuns.jl

Co-authored-by: David Widmann <devmotion@users.noreply.github.com>

---------

Co-authored-by: David Widmann <devmotion@users.noreply.github.com>

Author: oscardssmith

src/basicfuns.jl

@@ -62,7 +62,7 @@ end
 """
 $(SIGNATURES)
 
-Return `x * exp(y)` for `y > -Inf`, or zero if `y == -Inf`.
+Return `x * exp(y)` for `y > -Inf`, or zero if `y == -Inf` or if `x == 0` and `y` is finite.
 
 ```jldoctest
 julia> xexpy(1.0, -Inf)
@@ -72,7 +72,7 @@ julia> xexpy(1.0, -Inf)
 function xexpy(x::Real, y::Real)
     expy = exp(y)
     result = x * expy
-    return iszero(expy) && !isnan(x) ? zero(result) : result
+    return (iszero(x) && isfinite(y)) || (iszero(expy) && !isnan(x)) ? zero(result) : result
 end
 
 # The following bounds are precomputed versions of the following abstract

test/basicfuns.jl

@@ -73,6 +73,9 @@ end
     for x in (randn(), randn(Float32))
         @test xexpy(x, x) ≈ xexpx(x)
     end
+    @test xexpy(0, 1000) == 0.0
+    @test isnan(xexpy(0, Inf))
+    @test isnan(xexpy(0, NaN))
 end
 
 @testset "logistic & logit" begin