JuliaDiff / DualNumbers.jl

Julia package for representing dual numbers and for performing dual algebra
Other
80 stars 30 forks source link

Using a strong zero in Dual(0.0, 1)^0 to avoid NaN #84

Closed jwscook closed 3 years ago

jwscook commented 3 years ago

I found another minor bug-ette:

julia> Dual(0.0, 1)^0
Dual{Float64}(1.0,NaN)

The behaviour from this PR gives:

julia> Dual(0.0, 1)^0
1.0 - 0.0ɛ

What do you reckon?

dlfivefifty commented 3 years ago

Shouldn't this also work for (::Dual)^0.0 based on:

julia> 0.0^0.0
1.0

julia> h = 0.000001; ((0+h)^0.0 - 0.0^0.0)/h
0.0

Also make sure there are tests for dual(0,0)^0 (which should return dual(1,0)) as well as 0.0^dual(0,1) (which should return dual(1,NaN)) and finally dual(0,1)^dual(0,1) should return dual(1,0).

jwscook commented 3 years ago

Pleases and thank yous go a long way. Dual(0, 0)^0 should return a DomainError - see the tests.

mlubin commented 3 years ago

For reference, this is a longstanding and known issue, see https://github.com/JuliaDiff/ForwardDiff.jl/blob/3fb16ac21ac06b7c42c807e898395ad3bd9f54ae/docs/src/user/advanced.md#fixing-naninf-issues. I guess DualNumbers isn't especially focused on performance, so the 5%-10% overhead may be acceptable.

dlfivefifty commented 3 years ago

Dual(0, 0)^0 should return a DomainError

Why?

julia> 0^0
1

julia> (0+eps())^0
1.0
jwscook commented 3 years ago

In the previous implementation Dual(0, 1)^0 performed a 0^-1 which threw a DomainError. Dual{T} where {T<:Integer} still throws rather than promote to float.

I've also fixed the type instability of Dual^Dual to promote from Dual{T} where {T<:Integer} to Dual{T} where {T<:AbstractFloat} because it has to use a log(::Integer) anyway.

jwscook commented 3 years ago

Is there anything remaining that needs to be done before this is merged (other than semver bump)?

dlfivefifty commented 3 years ago

Pleas do the semvar bump