What about an `undefs` function in `Base`?

[uniment]

It seems we’re in agreement: we don’t want perfect newbies to be allocating uninitialized arrays, so it makes sense to have a mild hurdle, but Array{Float64}(undef, dims...) is too big a hurdle and is the wrong litmus test—knowledge of Julia’s parametric type system is largely orthogonal to whether someone will fill in their arrays properly.

From here, it seems mostly a debate of allocate or undefs.

The trade-offs I see for allocate vs undefs:

Pros for allocate:

1. Name is consistent as a counterpart to fill
2. Name makes perfectly clear that we are not filling the array with anything, but merely allocating it. This is conceptually pure.

Pros for undefs:

1. Argument signature undefs([T,] dims...) is more consistent with zeros, ones, and rand than with fill
2. If it’s ever decided that zero and one should become array initializers, e.g. Array{Float64}(zero, dims...) then it will have perfect consistency with the use of undef.

Overall, I’m leaning toward allocate. It’s more verbose, and it’s a departure from the style of zeros and ones, but that’s a good thing—it’s a healthy reminder that what we’re doing shouldn’t be taken lightly. It also feels like part of a natural progression for departing from the mathematical expression which called for zeros and ones, toward a programmatically optimized expression using allocate and fill!.

Sidenote: I know Stefan has expressed some dissatisfaction with zeros and ones as not being generic enough, but I disagree with that sentiment. I see it as an artifact of Julia’s strong heritage in numerical computing, which is a heritage to be proud of.

Certainly better than JavaScript’s heritage which causes this type of artifact:

javascript> [typeof([]), typeof({}), typeof([] + {}), typeof([] * {})]
['object', 'object', 'string', 'number']

[John_Gibson]

@DNF

I’m actually in favor of an undefs function, I just think it should be called allocate (or something).

Better to call it undefs to be parallel to ones and zeros. The name fill already a bit off-kilter, as it both allocates and fills. A function named allocate would further muddy the waters.

[mkitti]

If we were to have an allocate method, it really should take an optional, if not required, array type as well. We really want to avoid the imitation situation that NumPy got itself into by not being extensible.

allocate(Array, Int, 4, 5)

The mechanism should be extendable to other array types.

allocate(OffsetArray, Float64, 2:5)
allocate(CuArray, UInt8, (128, 256))

In fact, we have precedence for a such a function in unsafe_wrap.

Why not just call the method Array?

Array(Int32, 16, 32)

Now it’s time for the history lesson. Back in the day in Julia 0.3 you could construct Array via the following syntax:

Array(type, dims)

https://docs.julialang.org/en/v0.3/stdlib/arrays/#Base.Array

The next step were the curly braces in Julia 0.5:

Array{Int}(dims)

https://docs.julialang.org/en/v0.5/stdlib/arrays/#Base.Array

undef and the other array initializers came last.

Array{Int32}(undef, 4, 20)

https://docs.julialang.org/en/v0.7/base/arrays/#Core.Array-Tuple{UndefInitializer,Any}

The other array initializers are nothing and missing, by the way.

postscript:
Also note that we have similar which comes quite close. You can create a undef array without any curly braces.

julia> similar([], Int, 5, 3)
5×3 Matrix{Int64}:
 139988021321736  139988021321736  139988021321736
 139988021321736  139988021321736  139988021321736
 139988021321736  139988021321736  139988021321736
 139988021321736  139988021321736  139988021321736
 139988021321736  139988021321736  139988021321736

julia> similar([], Vector, 2, 3)
2×3 Matrix{Vector}:
 #undef  #undef  #undef
 #undef  #undef  #undef

[DNF]

@John_Gibson

Better to call it undefs to be parallel to ones and zeros.

That was my main point, I think it deliberately should avoid that parallel.

Also, undefs is just a really weird name.

[uniment]

Are we missing anything?

julia> allocate(::Type{A}, ::Type{T}, dims...) where {A,T} = A{T}(undef, dims...)
       allocate(::Type{T}, dims...) where T = allocate(Array, T, dims...)
       allocate(dims...) = allocate(Float64, dims...)
allocate (generic function with 3 methods)

julia> allocate(2, 2)
2×2 Matrix{Float64}:
 5.0e-324  7.13625e-312
 5.0e-324  7.13625e-312

julia> allocate(Bool, 2, 2)
2×2 Matrix{Bool}:
 0  0
 1  1

julia> using CUDA

julia> allocate(CuArray, Float64, 2, 2)
2×2 CuArray{Float64, 2, CUDA.Mem.DeviceBuffer}:
 0.0  0.0
 0.0  0.0

[mkitti]

You should be able to accept dims as a single tuple.

julia> A = Array{Int8}(undef, (2,3,4));

julia> B = Array{Int16}(undef, size(A));

See this issue for other thoughts:

design of array constructors

[jar1]

In this API there’s no room to say where the memory comes from, we just pull it out of the sky. If I have a block of memory, I might want to split it up and make arrays and dictionaries out of it, e.g. using allocate!(block, Array, Float64, (2,2)).

[uniment]

As I’ve written it here, allocate just forwards all responsibility to Array, which handles the rest.

julia> A = allocate(Int8, (2, 3, 4))
2×3×4 Array{Int8, 3}:
[:, :, 1] =
 0  0  0
 0  0  0

[:, :, 2] =
 0   16  -105
 0  -63   -84

[:, :, 3] =
  -6  0  16
 127  0  81

[:, :, 4] =
 -105   -6  0
  -84  127  0

julia> B = allocate(eltype(A), size(A))
       (typeof(A), size(A)) == (typeof(B), size(B))
true

[mkitti]

@jar1

In this API there’s no room to say where the memory comes from, we just pull it out of the sky

That’s one direction I thought of going with ArrayAllocators.jl. One could create a BumpAllocator. I still need to think about how to handle garbage collection with that though.

This is a good point though that the above is conflating memory allocation with array construction. ArrayAllocators.allocate is primarily meant for memory allocation. This is invoked behind the scenes after the number of bytes needed is computed.

julia> using ArrayAllocators

julia> ArrayAllocators.allocate(calloc, 1024)
Ptr{Nothing} @0x0000000004270990

julia> ArrayAllocators.allocate(malloc, 1024)
Ptr{Nothing} @0x0000000004271a20

julia> ArrayAllocators.allocate(MemAlign(1024), 1024)
Ptr{Nothing} @0x0000000004277400

[mkitti]

@uniment

julia> allocate(::Type{A}, ::Type{T}, dims...) where {A,T} = A{T}(undef, dims...)
       allocate(::Type{T}, dims...) where T = allocate(Array, T, dims...)
       allocate(dims...) = allocate(Float64, dims...)

One thing I’m not a fan of here is that we are peeling off optional arguments from the beginning of the argument list whereas I usually want to think about optional arguments occurring at the end or as keywords. There is also an ambiguity of what allocate(Array, 5, 3) might do. Does it make an Array{Array} or Array{Float64}? Maybe those would be better as keyword arguments.

allocate(dims...; array_type = Array, element_type = Float64, init = undef) =
    array_type{element_type}(init, dims...)

[Sukera]

@mkitti

As Carsten Bauer mentioned, there is a better solution than undefs for the novice user where they can obtain nearly the same performance. This is calloc, a standard C function. On operating systems, there also more specialized methods in the same vein.

“Better” is a very subjective term here, as has already been explained in the very thread you linked. There is more to the story there than “just” dropping in calloc as a replacement to malloc and be done with it.

@mkitti

I have since packaged this into ArrayAllocators.jl. Let me offer a brief demonstration.

This benchmark is misleading because it does not zero the memory until sum is called. Please don’t sweep the big caveat of using calloc under the rug and present it as a pure win in performance - it just isn’t. I’ve already benchmarked this in the other thread extensively.

@DNF

Would that work for all kinds of data types? zeros really only makes sense for numerical data. What about allocations of unitialized data in general?

No - this only “works” for types where a zero bitpattern also happens to coincide with the zero of that type, as I mentioned in the linked thread:

@Sukera

In fact, for user defined types zeros has to use malloc with initialization because a zero-filled block of memory may not be valid for zero of a user defined type! As you note, it may be boxed or a number of other things, like not being isbits or simply not having a block of zeros as its zero representation.

In which case, using calloc by default means the memory is initialized twice, once by the kernel when it gives you the memory that’s filled with 0, and once when julia has to inevitably call zero to initialize the data properly. This will lead to users being confused about why zeros is slower than it needs to be and should thus be avoided.

@mkitti

The knowledge I’m using above is:

1. The underlying memory representation of Vector{Vector{Int}} is an array of pointers.
2. #undef actually means those pointers are C_NULL (in other words, 0)
3. When set, the pointers point to a jl_array structure which has 48 bytes of meta data before the array contents.

All of these are internal implementation details you cannot rely on. Do not assume these to be true.

---

The interface for allocation I would prefer to all options presented here looks like this:

"""
Allocate a single object of type `T`, using memory managed by `Allocator`.
Return a `T`, throws an OutOfMemory exception when it fails to allocate memory.
"""
allocate(::Allocator, ::T)

"""
Allocate `n` objects of type `T`, using memory managed by `Allocator`.
Return a collection of `T`, throws an OutOfMemory exception when it fails to allocate memory.
"""
allocate(::Allocator, ::T, n)

This is how Zig does it, and for good reason - you don’t need more. Rust has some more fanciness for deallocation. Both of these can be used in the compiler to hoist allocations to the stack if need be. The second one semantically returns an Array (or a fixed-size equivalent) when passed a default GC.

@jar1

In this API there’s no room to say where the memory comes from, we just pull it out of the sky. If I have a block of memory, I might want to split it up and make arrays and dictionaries out of it, e.g. using allocate!(block, Array, Float64, (2,2)).

That’s because “splitting the memory up” to create other objects out of is not a safe operation - if you allocate a large array and use that as the “backing memory” for multiple objects, you’re well on the way of reimplementing a memory manager yourself. For the regular julia GC, it would have to keep that big chunk around as long as even a single other object (that may use only a tiny portion of that memory) is still accessible. You really need the allocator passed in the API to handle that part of the operation for you.

[mkitti]

@Sukera

No - this only “works” for types where a zero bitpattern also happens to coincide with the zero of that type, as I mentioned in the linked thread:

Does Base allocation handle this any differently?

julia> struct Foo
           x::Int
           Foo() = new(5)
       end

julia> Foo()
Foo(5)

julia> Array{Foo}(undef, 5)
5-element Vector{Foo}:
 Foo(0)
 Foo(0)
 Foo(0)
 Foo(0)
 Foo(0)

[Sukera]

I really don’t want to rehash the same conversation again, but undef just happens to give you whatever is in memory. That’s precisely what undef means - there are no constructors involved, ever.

julia> struct Foo
           bar::Int
       end
julia> Base.zero(::Type{Foo}) = Foo(1)

julia> zeros(Foo, 5)
10-element Vector{Foo}:
 Foo(1)
 Foo(1)
 Foo(1)
 Foo(1)
 Foo(1)

julia> while true
           x = Array{Foo}(undef, 5)
           if any(!=(Foo(0)), x)
               display(x)
               break
           end
       end
5-element Vector{Foo}:
 Foo(140600521258544)
 Foo(140600521258608)
 Foo(140600521258672)
 Foo(140600521257344)
 Foo(0)

In contrast, zeros does call zero:

julia> @edit zeros(Foo, (5,))

Opens editor here:

for (fname, felt) in ((:zeros, :zero), (:ones, :one))
    @eval begin
        $fname(dims::DimOrInd...) = $fname(dims)
        $fname(::Type{T}, dims::DimOrInd...) where {T} = $fname(T, dims)
        $fname(dims::Tuple{Vararg{DimOrInd}}) = $fname(Float64, dims)
        $fname(::Type{T}, dims::NTuple{N, Union{Integer, OneTo}}) where {T,N} = $fname(T, map(to_dim, dims))
        function $fname(::Type{T}, dims::NTuple{N, Integer}) where {T,N}
            a = Array{T,N}(undef, dims)
            fill!(a, $felt(T))
            return a
        end
        function $fname(::Type{T}, dims::Tuple{}) where {T}
            a = Array{T}(undef)
            fill!(a, $felt(T))
            return a
        end
    end
end

which defines

function zeros(::Type{T}, dims::NTuple{N, Integer}) where {T,N}
    a = Array{T}(undef, dims)
    fill!(a, zero(T))
    return a
end

[Eben60]

That similar can be used to create arrays having non-similar dimensions is not exactly intuitive. Still, the form

similar(Float64[], 2,3)

may be easier to remember than

Array{Float64}(undef, 2, 3)

So, thank you

[aplavin]

@Eben60

That similar can be used to create arrays having non-similar dimensions is not exactly intuitive.

similar doesn’t refer to similar dimensions specifically. Its meaning is basically “as similar to the passed array as possible”, optionally specifying the target eltype or dimensions.

[sijo]

@John_Gibson

The name fill already a bit off-kilter, as it both allocates and fills. A function named allocate would further muddy the waters.

We could add a fill method that takes an additional parameter to specify a pre-allocated array that should be filled. The allocate and fill names would make sense then, and it would not be surprising that fill can also allocate itself as a convenience, if no array is given (i.e. the current behavior).

[mcabbott]

@aplavin

Its meaning is basically “as similar to the passed array as possible”, optionally specifying the target eltype or dimensions.

Note that in addition to methods accepting an instance, similar also has some accepting a type. These are closer to the undefs idea, but they aren’t very flexible:

x = [1 2]  # Array, could be SArray, or CuArray, etc.
similar(x)               # instance
similar(x, 1, 3)         # instance -- change size
similar(x, Float32, 3)   # instance -- change eltype & ndims
similar(typeof(x), 1, 3) # type -- must specify size (unless static)

# These all give an error:
similar(typeof(x), Float32, 1, 3)  # type -- cannot change eltype
similar(typeof(x), 3)  # type -- cannot change ndims
similar(Array, Int, 3)  # type -- cannot use partial type

Defining the last three here would, as a side-effect, give another spelling for Array{Int}(undef, 3). But it would also be useful elsewhere, for handling stranger arrays, where you cannot assume that the first type parameter is the element type. This would allow e.g. stack(Vector{Int}[]) to work (and be type-stable) without seeing an instance.

@sijo

fill method that takes an additional parameter to specify a pre-allocated array

This is fill!(A, val).

[DNF]

This is what fill! does currently. The trailing ! indicates that it mutates one of its input arguments.

[uniment]

Is there a syntax for array comprehensions of different array types, such as CuArray? Maybe such a thing would diminish the need for pre-allocation and for loops to begin with. Something like this, perhaps?

Type{CuArray{Int}}[i*j for i=1:2, j=1:2]

It’s not consistent with the ElType[...] syntax we’re accustomed to, but StaticArrays already broke that expectation with SA[...]. It’s nice that Type{...} doesn’t provide any constructors, so it wouldn’t break anything.

[mkitti]

I created a package ArrayInitializers.jl that provides a set of initializer types that can be passed as the first argument to an array constructor. Since the initializer often contains information about the element type, one can avoid the use of curly braces.

julia> using Pkg

julia> pkg"add https://github.com/mkitti/ArrayInitializers.jl"
...

julia> using ArrayInitializers

julia> Array(zeroinit(Int), 3, 9)
3×9 Matrix{Int64}:
 0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0

julia> Array(zeroinit(Float64), 3, 9)
3×9 Matrix{Float64}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0

julia> Array(oneinit(Float64), 3, 9)
3×9 Matrix{Float64}:
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0
 1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0  1.0

julia> const fives = init(5)
ArrayInitializers.FillArrayInitializer{Int64}(5)

julia> Array(fives, 3)
3-element Vector{Int64}:
 5
 5
 5

julia> Array(undeftype(Int), 3)
3-element Vector{Int64}:
       1
       3
 7105903

This also works with other kind of arrays such as OffsetArrays or BitArrays:

julia> using OffsetArrays

julia> OffsetArray(fives, 2:5)
4-element OffsetArray(::Vector{Int64}, 2:5) with eltype Int64 with indices 2:5:
 5
 5
 5
 5

julia> BitArray(oneinit, 5)
5-element BitVector:
 1
 1
 1
 1
 1