Why is `UniformScaling` not a subtype of `AbstractArray`?


#1

Say I have the following line of code in Julia

C = A*B

where A and B are matrices, i.e.

typeof(A) == typeof(B) <: AbstractMatrix

If A is a multiple of the appropriately sized identity matrix then it is convenient to represent it as A = UniformScaling(x) for some scalar x. That works great except when you need to create A in one function then pass it into another. UniformScaling is not a subtype of AbstractArray (this is the case in Julia 0.5 and on master) so the following code doesn’t work

f(A::AbstractMatrix,B::AbstractMatrix) = A*B

A = UniformScaling(1.0);
B = sprandn(100,100,0.1);
C = f(A,B)

which you can quickly confirm at the REPL:

julia> f(A::AbstractMatrix,B::AbstractMatrix) = A*B
f (generic function with 1 method)

julia> A = UniformScaling(1.0);

julia> B = sprandn(100,100,0.1);

julia> C = f(A,B)
ERROR: MethodError: no method matching f(::UniformScaling{Float64}, ::SparseMatrixCSC{Float64,Int64})
Closest candidates are:
  f(::AbstractArray{T,2}, ::AbstractArray{T,2}) at REPL[1]:1

julia> 

I’m wondering if that is the intended behaviour or a bug? The one reason I can think of off the top of my head for why UniformScaling shouldn’t be a subtype of AbstractArray is that methods that should be supported for any array type (e.g. size, ndims) don’t exist for UniformScaling. I don’t know enough to know if that’s reason to nix the idea or not.

If the consensus is that making UniformScaling a subtype of AbstractArray is a good idea, I’ll make a PR. If not, does anyone know of another way to generate a placeholder for a very large identity matrix with a negligible memory footprint that can be passed to functions that take AbstractMatrix inputs?


#2

I think this is the case: it doesn’t support the full Array interface. But I would like to see a UniformScaling with size that does satisfy the full interface.

BTW, this is usually fixed by just not strictly typing the functions so much.


#3

Thanks for the reply Chris. Requiring something to support the full array interface in order to be an array is fair enough. I agree that adding something akin to UniformScaling with a size would be nice. I’ll have to think a bit more about that but in the meantime I guess I’ll just have to loosen my function signature typing.


#4

Was an issue or PR ever created? Almost a year later, I ran into essentially the same problem for a built-in matrix function:

julia> kron(rand(3,2),UniformScaling(3))
ERROR: MethodError: no method matching kron(::Array{Float64,2}, ::UniformScaling{Int64})
Closest candidates are:
  kron(::Any, ::Any, ::Any, ::Any...) at operators.jl:466
  kron(::Union{Array{T,1}, Array{T,2}} where T, ::SparseArrays.SparseMatrixCSC) at /Users/vagrant/worker/worker/package_osx64/build/usr/share/julia/site/v0.7/SparseArrays/src/linalg.jl:823
  kron(::AbstractArray{T,2}, ::AbstractArray{S,2}) where {T, S} at /Users/vagrant/worker/worker/package_osx64/build/usr/share/julia/site/v0.7/LinearAlgebra/src/dense.jl:379
  ...

julia> versioninfo()
Julia Version 0.7.0-DEV.4698
Commit a9eaa4f (2018-03-26 18:49 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin14.5.0)
  CPU: Intel(R) Core(TM) i7-4790K CPU @ 4.00GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-3.9.1 (ORCJIT, haswell)
Environment:

I’m also puzzled why the deprecation warning for eye tells me to use Matrix(1.0I,3,3) rather than the sparser and more precise Diagonal(ones(3)). Is there a downside to the latter, or was this a case of “just give them an array”? (Of course, I suspect that ones probably isn’t long for this world, either.)

I know there have been a lot of electrons spilled on the problems with convenience functions such as eye, and on the meaning of UniformScaling. From outside the debate, it doesn’t feel fully sorted out yet. As a long-time Matlab fanboy, I hope the Julia powers that be don’t lose sight of the incredible value Matlab created (and numpy piggybacked on) through unity and simplicity for interactive computing.


#5

eye returned a dense array and not a Diagonal, so it’s not the same (but is better).


#6

As far as I know (and can tell by searching the main Julia repo) there haven’t been any PRs or issues directly asking for UniformScaling to become a subtype of AbstractArray. From this issue
It seems that at least Andreas Noack would prefer that they not be sized, which would be a prerequisite to AbstractArray status.


#7

What if something like the following is defined:

(J::UniformScaling{T})(d::Int) where {T} = Diagonal(fill(J.λ, d))

Then you could construct dimensioned matrices …

Julia-0.6.2> I(3)
3×3 Diagonal{Int64}:
 1  ⋅  ⋅
 ⋅  1  ⋅
 ⋅  ⋅  1

Julia-0.6.2> 1.0I(3)
3×3 Diagonal{Float64}:
 1.0   ⋅    ⋅
  ⋅   1.0   ⋅
  ⋅    ⋅   1.0

#8

Hmm, that might be a great idea! Could combine this with something like the Eye type from FillArrays.jl to get zero allocations.

Edit: or rather

struct ScalarMatrix{T} <: AbstractMatrix{T}
    val::T
    n::Int
end
Base.size(A::ScalarMatrix) = (A.n, A.n)
Base.getindex(A::ScalarMatrix, i::Int, j::Int) = ifelse(i == j, A.val, zero(A.val))
Base.:*(x::Number, A::ScalarMatrix) = ScalarMatrix(x * A.val, A.n)
Base.:*(A::ScalarMatrix, x::Number) = x * A

(J::UniformScaling{T})(d::Int) where {T} = ScalarMatrix(J.λ, d)

#9

I guess this was proposed in https://github.com/JuliaLang/julia/issues/23156. Unfortunately the proposal was not implemented, despite some initial support.