3-arg `dot` is type-unstable for matrices with eltype `Bool`

[dkarrasch]

This is due to the fact that when we initialize the zero value, we use the eltypes or, alternatively, the first values, before we start to add up further dot products. We also have this inconsistency:

julia> using LinearAlgebra

julia> dot(true, true)
true

julia> dot([true], [true])
1

A minimal fix would be to define

dot(x::Bool, y::Bool) = x*y + false

which changes the return type for scalar arguments, but renders the matrix cases (both 2-arg and 3-arg) type-stable. What do you think?

[araujoms]

I think the correct return type is anyway Int. That is what one gets when one does addition/subtraction with booleans, and that's what dot does in the general case. Even in the unary case that's what Julia does:

julia> +(true)
1

julia> +(false)
0

As for your solution, that's essentially what the vector version is doing, so that is a good way to ensure results will stay consistent. I'm just afraid the extra addition will impact performance. Wouldn't it be better to do

dot(x::Bool, y::Bool) = Base.promote_op(+, Bool)(x * y)

or simply

dot(x::Bool, y::Bool) = Int(x * y)

?

[martinholters]

Or change the Number method to which the Bool case dispatches to include +(...), i.e.

dot(x::Number, y::Number) = +(conj(x) * y)

? I wouldn't expect that to make a difference for any type other than Bool, and if anything, it should hold that +(z) = z + zero(z) which is exactly what is required for type stability.

[dkarrasch]

dot(x::Number, y::Number) = +(conj(x) * y)

I was concerned this might lead to extra work in the complex case, for which

+(z::Complex) = Complex(+real(z), +imag(z))

(I used Float16 underneath to make sure I'm not running BLAS), but I don't see any, so we could actually go with this.

[araujoms]

If you go for this version please write a comment saying it's for fixing the Boolean case, otherwise it will look very mysterious.

[stevengj]

It's not clear to me why the scalar version needs to be changed to give the same type as the vector version? But I guess it doesn't hurt.

Note that a similar issue arises for mapreduce — see JuliaLang/julia#59564, and we could probably use the same unary-+ solution there.