Complex Differentiation

[juliadiff]

Complex Differentiation

[dlfivefifty]

Just use a Complex{<:Dual}.

Note DualNumbers.jl takes the other approach you advocate but is planned to be redesigned as ForwardDiff.Dual

[dlfivefifty]

What could be added is support for derivative(f, im) but that doesn't require a redesign of Dual.

[akriegman]

Thanks for the response. I tried that solution and there were some problems.

Complex{Dual} worked for some functions like polynomials, but didn't work for functions like sin. If we look at the source code for sin(::Complex), which gets called when we try sin(::Complex{Dual}):

function sin(z::Complex{T}) where T
    F = float(T)
    zr, zi = reim(z)
    if zr == 0
        Complex(F(zr), sinh(zi))
    elseif !isfinite(zr)
        if zi == 0 || isinf(zi)
            Complex(F(NaN), F(zi))
        else
            Complex(F(NaN), F(NaN))
        end
    else
        s, c = sincos(zr)
        Complex(s * cosh(zi), c * sinh(zi))
    end
end

So if we call derivative(sin, im) (after overloading derivative and extract_derivative to use Complex{Dual}) we end up in the first code path because zr === Dual(0, 1) == 0. However, this gives the wrong answer because zr != sin(zr) * cosh(zi). These expressions have the same real part, but different infinitesimal parts.

Simply passing a Dual into the implementation of sin(::Real) from Base wouldn't work for similarish reasons, which is why ForwardDiff provides a method for sin(::Dual)

function Base.sin(d::Dual{T}) where T
    s, c = sincos(value(d))
    return Dual{T}(s, c * partials(d))
end

If we wanted to use Complex{Dual}, we would have to make this method and methods like it take both Dual and Complex{Dual}. Alternatively, we could use Dual{Complex}. It seems to me that the trade off is between replacing Dual with Union{Dual, Complex{Dual}} throughout the codebase or replacing Real with Number. The second option seems preferable, but the first might work too with some trickery.

As for DualNumbers.jl, that is an option, but it doesn't have a derivative interface to hide the dual numbers, and it isn't actively maintained. Maybe I could look into using TaylorSeries.jl. In any case I think ForwardDiff would benefit from this feature.

What could be added is support for derivative(f, im) but that doesn't require a redesign of Dual.

To be clear, I tried this when I was experimenting with Complex{Dual}, and ran into the problems described above.

[dlfivefifty]

Changing Duals supertype to Number will be incredibly disruptive and break a lot of code. I think there is no chance a PR doing this will be merged so I would recommend not putting time into that.

I don't think you would have to add Unions everywhere, only special cases where the default breaks (such as sin). For example, exp already works:

julia> exp(Complex(Dual(0,1),Dual(1,0)))
Dual{Nothing}(0.5403023058681398,0.5403023058681398) + Dual{Nothing}(0.8414709848078965,0.8414709848078965)*im

[akriegman]

Changing Duals supertype to Number will be incredibly disruptive and break a lot of code. I think there is no chance a PR doing this will be merged so I would recommend not putting time into that.

Okay, fair enough. Also, I just noticed that I posted this issue on DualNumbers and not ForwardDiff, my bad 😬

Taking a closer look at ForwardDiff, we wouldn't actually have to change that many places to take Complex{<:Dual} in addition to Dual. Most of the special function definitions for Duals are done with unary_dual_definition and binary_dual_definition, so we just have to change those and a couple other places. I might give that a go.