Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

UniformScaling should be assignable in setindex!(::AbstractMatrix) #34130

Open
jiahao opened this issue Dec 18, 2019 · 4 comments
Open

UniformScaling should be assignable in setindex!(::AbstractMatrix) #34130

jiahao opened this issue Dec 18, 2019 · 4 comments
Labels

Comments

@jiahao
Copy link
Member

jiahao commented Dec 18, 2019

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 jiahao changed the title UniformScaling should be assignable in setindex(::AbstractMatrix) UniformScaling should be assignable in setindex!(::AbstractMatrix) Dec 18, 2019
@jiahao
Copy link
Member Author

jiahao commented Dec 18, 2019

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
Copy link
Member

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
Copy link

daschw commented Dec 19, 2019

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 mbauman added the domain:linear algebra Linear algebra label Dec 19, 2019
@mbauman
Copy link
Sponsor Member

mbauman commented Dec 19, 2019

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.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

No branches or pull requests

4 participants