Functions that return arrays with eltype as input should use container type instead?

[tkelman]

ref #11280 (comment) and a few other places I've mentioned this. Functions like zeros, ones, eye, rand, mmap (maybe), etc which currently take an eltype as inputs and return some subtype of AbstractArray, would be much more general if they instead took the container type (with element type as a parameter) as input, so the current case would instead look like zeros(Array{Float64}, dims...). This would allow extending these functions to other types of AbstractArray in a much nicer way than right now. We could potentially get rid of a bunch of ugly-named functions like bitrand and spzeros and speye and sprandbool with this change.

I think it would by necessity be breaking, though maybe clever enough subtyping of zeros{T<:AbstractArray}(::Type{T}, dims...) would allow putting in deprecations (or defaults? might be fragile) for "scalar" types, if we have enough of a concept of interfaces to determine what makes a type a scalar (that's a different issue).

This would have to be a breaking change for arrays-of-arrays (if any of these functions currently work for arrays of arrays, that is?).

[ScottPJones]

Not sure anybody gives a hoot about my opinion, but I think this is a very good idea, and is exactly the sort of thing that should only happen once in the language, and the coming 0.5 "Arraymeggedon" is exactly the perfect time to do so. 👍

[ScottPJones]

Should it be called zeros() though? Could this better be more general? inittype(initclass, type, dims...)? Just a thought off the top of my head... zeros(T, dims...) = inittype(INITZERO, T, dims...),
ones(T,dims...) = inittype(INITONE, T, dims...), initclass might be from an enum... [ignore my names for the function and classes... I'm not great at naming so early in the morning])

[tkelman]

That doesn't work for things like rand and eye though, it's a bit narrow-minded and probably well-served already by fill!. Completely eliminating zeros and ones would probably not go over all that well with the linear algebraists among us, though I suppose nothing is guaranteed to be sacred from banishment to MatlabCompat.jl, at least not before Julia hits 1.0.

[ScottPJones]

Who said anything about eliminating them? That's why I showed how they could be defined using a new general function... I think you are the one who will get the 🍅s for daring to suggest taking away spzeros, etc. 😀

[tkelman]

That is actually already done in terms of fill! at

julia/base/array.jl

Lines 234 to 240 in d198e90

	for (fname, felt) in ((:zeros,:zero), (:ones,:one))
	@eval begin
	($fname)(T::Type, dims...) = fill!(Array(T, dims...), ($felt)(T))
	($fname)(dims...) = fill!(Array(Float64, dims...), ($felt)(Float64))
	($fname){T}(A::AbstractArray{T}) = fill!(similar(A), ($felt)(T))
	end
	end

- I guess the proposal here is to identify the rest of the places where we do this inconsistently and come up with something more uniform, apply it more broadly. Currently you can do zeros(A::AbstractArray), but not zeros{T<:AbstractArray}(::Type{T}, dims...).

The decision item is which is more general, needing to take an already-constructed instance of the AbstractArray type as input, or taking the array type itself as input to allow these functions to define their own direct construction?

[ScottPJones]

Well, I already stated that I liked the idea of passing the array type (actually, couldn't it be any collection type?)

[jiahao]

For matrices of matrices, the type information alone is insufficient to determine the appropriate zero, eye or one. The shape information of the inner nested matrix is also required. None of these functions currently work on matrices of matrices because, e.g. zero(::AbstractArray{T}) calls zero(T), but does not propagate the shape information of the current instance of T to the second call.

Thinking about zero(::SparseMatrix{<:SparseMatrix{T}}) and zero(::SparseMatrix{Matrix{T}}) makes my head hurt.

[tkelman]

Yeah, sparse arrays of sparse arrays, while godawfully cool in concept, would be pretty tough to support generically. A nerd can dream.

To give a concrete example where I ran into this as a major problem, I was making a separate sparse matrix type in a package using a different representation than CSC, and immediately found that there was no signature for spzeros (or zeros) that I could extend to construct an object of my new type. I wound up calling the function foozeros and being very displeased with that workaround.

[jiahao]

Arrays of arrays are not pleasant to work with in Julia. Recall #8759? There's definitely room for a variant of fill that does recursive copies.

[ScottPJones]

@tkelman, why are sparse arrays of sparse arrays so difficult? Is it just to handle the mathematical functions? (I've been using and implementing structures [multidimensional associative arrays], for the past 29 years, and they could be viewed/used as sparse arrays of sparse arrays)... sorry in advance for the potentially dumb question

[tkelman]

There are assumptions sprinkled around in various places that you can always call zero on the element type of an array once you try to do non-trivial mathematical operations on them. Constructing the data structure in the first place should work okay as long as the type and constructor are general enough to support that, but doing non-trivial math on arrays of arrays, or even just on sparse arrays while preserving and exploiting sparsity in as many operations as possible (ref JuliaLang/LinearAlgebra.jl#136), runs into trouble with the current state of base. It's gradually improving but there are still some missing pieces in terms of formalizing interfaces and making nontrivial operations generic over both structure and element types.

edit: Also please please refrain from conflating missing data with sparsity in the linear algebraic sense. They are rather different in implementation aspects and semantics. You're very familiar with data structures that are optimized for missingness, those exist outside of base in the Julia ecosystem, primarily for statistical purposes right now. Most of what I'm discussing is about linear algebra and data structures optimized for that purpose.