Julia: `*(::Number, ::AbstractArray)` calls `.*`

[mauro3]

This example

type A
    a::Float64
end

*(a::A, b::Real) = A(a.a*b)
*(b::Real, a::A) = a*b

[A(1)] * 5

fails with

ERROR: LoadError: MethodError: `.*` has no method matching .*(::A, ::Int64)
Closest candidates are:
  .*(!Matched::FloatRange{T<:FloatingPoint}, ::Real)
  .*(!Matched::Range{T}, ::Real)
  .*(!Matched::Range{T}, ::Number)
  ...
 in .* at array.jl:811
 in * at abstractarray.jl:446
 in include at ./boot.jl:250
 in include_from_node1 at loading.jl:129
 in process_options at ./client.jl:308
 in _start at ./client.jl:407

Of course it could be fixed by adding .*(a::A, b::Real) = a*b but seems to me that that breaks the abstraction of treating A and Real as scalar like entities. Same goes for other operators.

I had a look at the code, relevant definitions feature both in abstracarray.jl https://github.com/JuliaLang/julia/blob/master/base/abstractarray.jl#L445, in array.jl https://github.com/JuliaLang/julia/blob/master/base/array.jl#L799 and https://github.com/JuliaLang/julia/blob/master/base/operators.jl#L92.

[jiahao]

seems to me that that breaks the abstraction of treating A and Real as scalar like entities

Given only the information A<:Any, it is incorrect to deduce that A is scalar or even scalar-like (what does that even mean?). (Counterexample: A = AbstractArray.) Such erroneous assumptions have led to confusion over issues like in #9170, where @andreasnoack and I attempted to raise awareness of this fallacy. Others are apparently happy with special-casing AbstractArrays as the only non-scalars, but this is simply untenable, as we already have quite a few objects in base which are non-scalar non-AbstractArrays. As a result, we have an unhappy coexistence of scalar and non-scalar types under the umbrella of Any, and we don't have a place in the type hierarchy that comfortably says scalar but not necessarily number.

Without redesigning the type hierarchy, I think the only thing you can do is define A <: Number in your example and define a promotion rule for (A, Int).

[mauro3]

Thanks! I'll need to ponder this for a bit.

[mauro3]

Conceptually, it would make sense to me that if *(a::A,b::B)=... is defined then all of these should work:

a = A(...); b = B(...)
a*b
[a]*b
a*[b]
[a].*[b]

and it shouldn't matters what A and B actually are, it should just work for A<:Any, B<:Any. However, there are several problems:
(1) say, if b is actually a matrix, then the following logical inconsistency would arise:

julia> a=1.; b = rand(2,2);

julia> a*b; # works

julia> [a]*b; # this should work by my above premise
ERROR: DimensionMismatch("*")
 in gemm_wrapper! at linalg/matmul.jl:272
 in * at linalg/matmul.jl:76

julia> [a a;a a]*b # and this should be equivalent to element wise multiplication!
2x2 Array{Float64,2}:
 1.01438  1.04651
 1.01438  1.04651

I think this is (part of) what @jiahao referred to in his post. To make this work, it would somehow be necessary to tag b to be a scalar, so Julia would know that the * is just a scaling and not a matrix multiplication (not possible at the moment).

(2) the problem with implementing it: what is the Array-type of [a]*b? (This is the problem map faces as well and cannot be solved in generality until we get function types. Edit: and also popped up in #11119) Thus the current implementation *(::AbstractArray, ::T)=... for T<:Number works because then the array-type can be determined through the promotion system.

But I think, even with above limitations, the current system could be improved on. When B<:Number then [a]*b for A<:Any should be well-defined as ==[a*b] as only scaling-multiplication is an option (my understanding of Number is that it is always a scalar).

To fix this, I think just the definition of the .* and * needs to be interchanged (along with +, etc), which
I find backwards anyway. Currently it is:

*(A::Number, B::AbstractArray) = A .* B # abstractarray.jl#L445

function (.*){T}(A::Number, B::AbstractArray{T}) # array.jl#L799
    F = similar(B, promote_array_type(typeof(A),T))
    for i in eachindex(B)
        @inbounds F[i] = (.*)(A, B[i])
    end
    return F
end

But * should be the main function and .* should be tacked-on for convenience. The caveat is that the return type for my initial example would be Array{Any,1} as the promotion is not defined. But I think it would still be an improvement. Or am I missing something? If not I'll try and make a PR.

[mauro3]

Another point, if instead of * define .*:

type A
    a::Float64
end

.*(a::A, b::Real) = A(a.a*b)
.*(b::Real, a::A) = a.*b

now this works:

julia> A(1) .* 5
A(5.0)

julia> [A(1)] * 5
1-element Array{Any,1}:
 A(5.0)

This also seems backwards.

[jiahao]

Indeed, there is something unsatisfactory about how the dotted operators are defined. The broadcasting operators only make sense on arrays, and they should automatically inherit the meanings of the undotted operators acting on individual elements.

[mauro3]

This should get fixed with #17623, I think.

[stevengj]

Fixed by #17623:

julia> import Base.*

julia> type A
           a::Float64
       end

julia> *(a::A, b::Real) = A(a.a*b)
* (generic function with 165 methods)

julia> *(b::Real, a::A) = a*b
* (generic function with 166 methods)

julia> [A(1)] * 5
1-element Array{A,1}:
 A(5.0)

julia> [A(1)] .* 5
1-element Array{A,1}:
 A(5.0)

[tkelman]

worth adding that as a test case, if you haven't added anything equivalent yet?

[oscardssmith]

Is there any reason why we can't add Scalar as a type that numbers, chars, and bools etc are subtypes of?

[stevengj]

No need for a type. Broadcast defaults to treating arguments as "scalar" if they are not arrays etc

[stevengj]

Will add a test case