Julia: copysign and flipsign aren't type stable with Irrational input

[giordano]

julia> copysign(pi, 1)
π = 3.1415926535897...

julia> copysign(pi, -1)
-3.141592653589793

julia> flipsign(catalan, 1)
catalan = 0.9159655941772...

julia> flipsign(catalan, -1)
-0.915965594177219

I'm not sure what's the best solution. Maybe defining

copysign(x::Irrational, y::Real) = copysign(float(x), y)
flipsign(x::Irrational, y::Real) = flipsign(float(x), y)

?

[TotalVerb]

It seems to me that you're better off just calling copysign(float(x), y) if you want this inferrable. What's the problem with the type instability here?

[giordano]

What's the problem with the type instability here?

I don't have particular problems, but thought it was preferable to have functions in Base type-stable when possible.

[perrutquist]

The fact that -π is Float64, not Irrational may lead to hard-to-detect errors if the programmer assumes that negation will not cause loss of precision. For example:

julia> BigFloat(π) + BigFloat(-π)
1.224646799147353177226065932275001058209749445923078164062861980294536250318213e-16

(I once started writing an IrrationalExpressions package that could resolve this, but it's been on hold for a while.)

[perrutquist]

A simple fix would be to create negative counterparts to each Irrational constant. (There are only about five of them.) If -(::Irrational) returns an Irrational, then the existing copysign and flipsign will become type-stable (albeit in a weak sense as Irrational{:π} and Irrational{:negativeπ} would technically be different types.)

[StefanKarpinski]

We've been quite strict about not making Irrationals an algebraic type by adding more and more arithmetic to them – because that is a very slippery slope – but I suppose negation is simple and restricted enough that this could be defensible.

[giordano]

We've been quite strict about not making Irrationals an algebraic type by adding more and more arithmetic to them – because that is a very slippery slope

I completely agree. I'm often hit by the fact that 2pi * x, with x being a BigFloat, gives an unexpectedly inaccurate result, but I know that and now I'm teaching myself to always use x * 2 * pi. Yes, that may be a bit annoying, but the reason why 2pi * x (and x * 2pi as well) will give inaccurate results are completely understandable.

but I suppose negation is simple and restricted enough that this could be defensible.

That would also complement +pi, which conveniently returns pi itself (but of course that's much easier).

[StefanKarpinski]

I've occasionally thought of extending irrationals to rational multiples of irrational values, but then you have a problem in that some of the computation with big floats may end up happening at runtime, which is extremely undesirable, so it still seems best not to go there.

[nalimilan]

If we add negation, will we be able to resist and not add the inverse later? I'm afraid not...

[StefanKarpinski]

Yes, that's precisely the slippery slope here.

[giordano]

If we add negation, will we be able to resist and not add the inverse later? I'm afraid not...

I think inversion is much more complex than negation. Negation changes only the sign, inversion affects the module. Personally, I wouldn't support even 1x.

[perrutquist]

What if irrational arithmetic were allowed via a macro? So @irrational sqrt(2) would return an object of type Irrational{sqrt(2)} and computation with BigFloat would happen at compile-time via generated functions?
(Although at present it seems type parameters can't be expressions...)

[StefanKarpinski]

That macro could be implemented in a package, which seems like a good way to explore the idea.

[kmsquire]

I'm often hit by the fact that 2pi * x, with x being a BigFloat, gives an unexpectedly inaccurate result

If you use 2pi frequently, you could use the Tau.jl package here.

[eschnett]

It seems that the copysign method called here is

copysign(x::Real, y::Real) = ifelse(signbit(x)!=signbit(y), -x, +x)

I think copysign should be type-stable. To ensure this, we can call promote_type on the two possible results, as in

function copysign(x::Real, y::Real)
    R = promote_type(typeof(+x), typeof(-x))
    R(ifelse(signbit(x)!=signbit(y), -x, +x))
end