Type stability issue in complex sin?

[malmaud]

Base.return_types(sin,(Complex{Int64},))

1-element Array{Any,1}:
 Union{Complex{Float64},Complex{Int64}}

[jiahao]

There are explicit checks for NaNs in the edge cases:

function sin(z::Complex)
    zr, zi = reim(z)
    if !isfinite(zi) && zr == 0 return Complex(zr, zi) end
    if isnan(zr) && !isfinite(zi) return Complex(zr, zi) end
    if !isfinite(zr) && zi == 0 return Complex(oftype(zr, NaN), zi) end
    if !isfinite(zr) && isfinite(zi) return Complex(oftype(zr, NaN), oftype(zi, NaN)) end
    if !isfinite(zr) && !isfinite(zi) return Complex(zr, oftype(zi, NaN)) end
    Complex(sin(zr)*cosh(zi), cos(zr)*sinh(zi))
end

which is why static analysis thinks it could be type unstable, but really it isn't because those branches can only be triggered for numeric types that support infinite representations. Many of the elementary functions have similar checks.

In principle we could have a separate function for Complex{<:FloatingPoint} but I don't think it is worth the effort to solve a non-problem.

[jiahao]

Actually the current implementation doesn't work exactly right for Complex{<:Rational} because there is no representation of NaN in Rationals, but that is a separate issue and I have a hard time believing that it will show up in real code

julia> sin(1//0 + 2im)
ERROR: ArgumentError: invalid rational: zero(Int64)//zero(Int64)
 in call at rational.jl:8
 in rationalize at rational.jl:95
 in convert at rational.jl:76
 in sin at complex.jl:525

[yuyichao]

But can RationalS be infinite?

[jiahao]

julia> 1//0
1//0

julia> isinf(1//0)
true

[yuyichao]

Ahh, right... (I thought it was not allowed but what's actually not allowed is 0//0...)

[yuyichao]

So on the other hand, for the "separate issue" we can make 0//0 be the NaN for Rational?

[jiahao]

In principle it is possible, but NaN is very tricky to get correct. For example, NaN != NaN in floating point but 0//0 == 0//0 if all you do is compare numerators and denominators. You'd have to ask if NaN semantics are worth the trouble.

[simonbyrne]

So on the other hand, for the "separate issue" we can make 0//0 be the NaN for Rational?

The main argument for the existence of NaNs is that they allow for exception-less code when invalid arguments are used: since Rationals (at least non-BigInt ones) can already throw OverflowError, there is no real reason to avoid another exception.

[StefanKarpinski]

For what it's worth, I originally had 0//0 as a rational NaN but it made the implementations of a lot of basic operations on rationals really crazy. Infinite rationals are both more useful and don't induce nearly so many corner cases in implementations.

[eschnett]

Shouldn't we determine the return type at the beginning of the function, and then ensure that all return statements return this type? This would probably implicitly convert the arguments to Float64 (or Float32 or BigFloat).

Alternatively, we could reject calls with integer or rational types, and require that these be explicitly converted to a floating-point type.

[andreasnoack]

Is this type instability really a problem? I wouldn't expect the input to sin to be Complex{Int} in any performance critical code.

If it is considered a problem, I think the method should be restricted to Complex{T<:FloatingPoint} and we should add a promotion method, i.e. something like sin(z::Complex) = sin(float(z)).

[andreasnoack]

fixed by #13539