Julia: UniformScaling should be assignable in setindex!(::AbstractMatrix)

[jiahao]

I was surprised that assignment of a UniformScaling into a matrix does not work. For example:

A = rand(4,6); A[1:4, 1:4] = I #Make this part of A the identity matrix

produces this error on Julia 1.3:

ERROR: ArgumentError: indexed assignment with a single value to many locations is not supported; perhaps use broadcasting `.=` instead?

The issue here is distinct from #23197 in that broadcasting semantics are not relevant; the context
of the setindex! call defines precisely what the implicit dimensions of UniformScaling ought to be in the assignment operation.

Here is a simple implementation that would work as intended in the original statement:

function Base.setindex!(A::AbstractMatrix, I::UniformScaling, is, js)
    n = min(maximum(is), maximum(js))
    imin = minimum(is)
    jmin = minimum(js)
    for i in is, j in js
        if i-imin == j-jmin
            A[i,j] = I.λ
        else
            A[i,j] = 0
        end
    end
end

Some use cases that work:

julia> A=rand(4,6); A[1:4,1:4] = I;A # Arguably quite common to want
4×6 Array{Float64,2}:
 1.0  0.0  0.0  0.0  0.147456  0.745478
 0.0  1.0  0.0  0.0  0.763384  0.445349
 0.0  0.0  1.0  0.0  0.886189  0.978133
 0.0  0.0  0.0  1.0  0.37166   0.272733

julia> A=rand(4,6); A[1:4,1:3] = I;A # Somewhat less common to want
4×6 Array{Float64,2}:
 1.0  0.0  0.0  0.58555   0.230785  0.197464
 0.0  1.0  0.0  0.364411  0.589343  0.83179
 0.0  0.0  1.0  0.758986  0.300437  0.190048
 0.0  0.0  0.0  0.208139  0.514587  0.782728

julia> A=rand(4,6); A[2:3,2:4] = I;A #A bit unusual to do this but it would work nonetheless
4×6 Array{Float64,2}:
 0.586206  0.219209  0.768404  0.0562488  0.403744  0.402825
 0.14575   1.0       0.0       0.0        0.284117  0.860951
 0.239251  0.0       1.0       0.0        0.688153  0.0935465
 0.171968  0.989361  0.99173   0.327323   0.439517  0.290657

[jiahao]

With this proposed implementation, a modification of the example given by @StefanKarpinski here would work (at least algebraically; not avoiding the intermediate construction of A[1:4,1:4] + I):

julia> A=zeros(4,6); A[1:4,1:4] += I; A
4×6 Array{Float64,2}:
 1.0  0.0  0.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0  0.0  0.0
 0.0  0.0  0.0  1.0  0.0  0.0

[fredrikekre]

julia> A = rand(3,3); A[[2,1],[1,3]] = I; A
3×3 Array{Float64,2}:
 1.0       0.840775   0.0     
 0.0       0.0859019  0.0     
 0.138221  0.30412    0.267277

[daschw]

Wouldn't

function Base.setindex!(A::AbstractMatrix, I::UniformScaling, is, js)
    for i in eachindex(is), j in eachindex(js)
        A[is[i], js[j]] = i == j ? I.λ : 0
    end
end

provide the expected results?

[mbauman]

Specializing setindex!(::AbstractMatrix, ::UniformScaling, is, js) would be ambiguity-inducing for all custom matrix types. I suppose we could do it on the internal Base._setindex!, or we could just implement it for ::Matrix.