Note that the layout of Array is about to change (although not in ways that will get rid of undef) Add `Memory` type by oscardssmith · Pull Request #51319 · JuliaLang/julia · GitHub
3 Likes
It’s not crazy, but we just have a variety of methods to construct it that aren’t attached to the type itself.
I’d rather not initialize like how C++ does, users seem to agree that it’s complicated and doesn’t even stop indeterminate values generally. The rule of thumb, very much emphasized in that blog, is to never rely on automatic things for initialization.
That second line would make a vector of 10 empty vectors, did you mean > b(10, a); to make a 10x10? That’s called the fill constructor, but unlike Julia’s fill, it copies the value to each element. Either way, can be done with a comprehension in Julia.
Thanks for pointing at sizehint!
Yeah, since the array is dynamic (can grow and shrink), it must distinguish capacity and size internally. Otherwise, for example, the code
a = Int[]
for i = 1:n
push!(a, i)
end
would cause reallocation on every iteration, turning this into a quadratic time operation.
What I meant was for this to be part of the array interface, and user controllable. sizehint! does seem to do this, although it appears to make few promises.
Yeah, I don’t mean that a = Vector{T}(n) should give a::Vector{Union{T,Nothing}}.
I mean that if you can’t use some mechanism to initialize at creation, then you can 1) construct a vector a::Vector{Union{T,Nothing}} initialized with nothing, 2) populate it with elements of type T, and 3) construct from a the vector b::Vector{T}. If the vector is never fully populated with elements of type T, or if this is done in a non-local manner, then maybe Vector{Union{T,Nothing}} (or some other type, like a dictionary or sparse vector) is a better type than Vector{T}.
A functional programming/type theory view is that code that produces data of some type is akin to a constructive proof that an element of that type exists (see Curry–Howard correspondence - Wikipedia). To me, this does not match Julia’s Vector{T} as every element in such a vector should be of type T - I feel Julia “fakes” this type.
It’s like
- Simple function f: Julia, I ready to do work! I can handle any vector of
T’s as input! - Julia: Excellent news f! We happen to have a lot to work today, could you please compute on this data for me?
- Simple function f: Sure can! I promise to do my best to complete as soon as possible!
- Simple function f: It is a vector of
T’s, right? - Julia: Yes, its totally a vector of
T’s, I promise! - Simple function f: *goes to work*
… - Simple function f: Julia, I found a nullpointer in the data!
- Julia: ****, no one will ever know, *blows bomb*
Julia doesn’t necessarily break Curry-Howard here. Constructing the Vector{T} itself works perfectly fine; it’s retrieving an element in that which ultimately throws. That is, having a container object is not equivalent with having that container actually be filled with sensible/retrievable/initialized data.
If julia would not throw here, then yes, it’d break constructor assumptions, which is exactly what happens with isbits types that only allow specific bitpatterns to exist, and is why I brought up that as an example above/in the other thread. This isn’t really related to what happens in the case of mutable data though.
Please don’t assign a gender to a programming language.
In particular, do not sexualize the term “Julia” or any other aspects of the project. While “Julia” is a female name in many parts of the world, the programming language is not a person and does not have a gender.
Sorry, edited.
1 Like
Then couldn’t we just pass the initialization mechanism to the C level?
As I noted above, _new_array_ indeed does perform zero-initialization, so I don’t see why more sophisticated initialization would be impossible.
I don’t understand this. To me this appears inconsistent with
and with
I would appreciate any reference to documentation or code that explains the memory safety design you speak of.
Vector{Int}'s elements don’t have pointers to check, Int is an isbits type. No pointers, no dereferencing problems like segfaults or corruption. Granted, you still shouldn’t index it before proper initialization because it’s indeterminate, just interprets preexisting bits in memory.
I mean f uses a::Vector{Int} without check, and Vector{Int} is not an isbitstype type.
Vector{Int} and Vector{Vector{Int}} have different element types. See the difference between your first post
and this post:
Int is an isbits type so indexing an element with its type doesn’t require pointer checks.
A method itself doesn’t need to check arguments for nonexisting instances, the error would be thrown before even reaching the method.
julia> global b::Vector{Int}
julia> f(a::Vector{Int}) = @inbounds a[1]
f (generic function with 1 method)
julia> f(b)
ERROR: UndefVarError: b not defined
...
Same goes for isbits types, you don’t get random memory bits.
julia> let
local c::Int
c+1
end
ERROR: UndefVarError: c not defined
...
Yes, that’s my point!
Recall, I was replying to
In my example, the reference to a is known to be valid, so no check is needed.
If we always had fully initialized arrays, then wouldn’t the same reasoning apply for the vector elements? When accessing an element of a::Vector{T}, for reference type T, why would Julia check the reference?
Correct logic, but we don’t have those. Variables can be determined to have defined or undefined references at compile-time, so the compiler may elide those checks (e.g. observe foo() = @inbounds a[1] when the global a is defined versus undefined). Container elements and struct fields cannot because methods that initialize the overall instance are written directly in Julia. The only way to accomplish what you’re aiming for is to sacrifice that feature by doing everything like:
, which goes against the design of Julia.
1 Like
Everyone agrees that it would be preferable if arrays and objects could be guaranteed to always come into existence fully constructed. But we also need to have optimally efficient ways to do all kinds of array or object constructions in the language. If there’s some construction that simply cannot be done efficiently, that’s not really acceptable. So far, the two proposed approaches to avoiding partial initialization are:
- Require providing a default object for all reference types that are used in arrays.
- Make the array construction API take an iterator to provide initial values.
I had a big comment written about this, but then I got to an example that made my points better than my wall of text, so here’s that example. Suppose you want to fill an array with Pascal’s triangle values. In Julia now you can do this:
function pascal(::Type{T}, n::Integer) where {T<:Real}
A = Matrix{T}(undef, n, n)
A[:, 1] .= one(T)
A[1, :] .= one(T)
for j = 2:n, i = 2:n
A[i, j] = A[i-1, j] + A[i, j-1]
end
return A
end
When T is a value type like Int then there’s no issue. But there can be number types that aren’t value types, like BigInt. This example could also be generalized to allow vectors as elements and allow the caller to provide initial values for the first row/column that aren’t all ones, in which case the element type would be something more complex—and even more importantly something mutable (more later).
Let’s focus on the T = BigInt case. This is what the allocation looks like currently:
julia> A = Matrix{BigInt}(undef, 3, 3)
3×3 Matrix{BigInt}:
#undef #undef #undef
#undef #undef #undef
#undef #undef #undef
Yes, those pesky undefs are there, but we didn’t have to allocate any BigInts. Suppose we did default initialization. For BigInts, zero is a good default value, so we can emulate this by replacing Matix{T}(undef, n, n) with zeros(T, n, n). Let’s time them for T = BigInt, n = 3:
julia> @belapsed Matrix{BigInt}(undef, 3, 3)
1.5405811623246494e-8
julia> @belapsed zeros(BigInt, 3, 3)
3.876814516129032e-8
That’s a 2.5x slowdown for default initialization. Not great. Keep in mind that in the pascal function we’re just going to replace all these values anyway, so this 2.5x slowdown just to fill the initial array with something valid is wasted.
We should also note that in this case default initialization is using the same instance of BigInt(0) to fill the entire matrix. If we consider BigInts to be immutable that’s fine, but they do actually have a mutable API, not exposed by Base, but available through the MutableArithmetics package, which can be used to write faster code and which lets us expose this fact:
julia> using MutableArithmetics
julia> A = zeros(BigInt, 3, 3)
3×3 Matrix{BigInt}:
0 0 0
0 0 0
0 0 0
julia> add!!(A[1], 123)
123
julia> A
3×3 Matrix{BigInt}:
123 123 123
123 123 123
123 123 123
Oops! All the entires changed. Now, one may say that we shouldn’t be surprised here and that the issue is that BigInts should be treated as immutable unless you allocate one yourself and can be sure that it’s safe to mutate it. But the proposal is to use default filling for all reference types, not just ones like BigInt which are conceptually immutable. In fact, the original example was to use Int[] as the default value for Vector{Int}. In that case we’d have this:
julia> A = fill(Int[], 3, 3)
3×3 Matrix{Vector{Int64}}:
[] [] []
[] [] []
[] [] []
julia> push!(A[1], 123)
1-element Vector{Int64}:
123
julia> A
3×3 Matrix{Vector{Int64}}:
[123] [123] [123]
[123] [123] [123]
[123] [123] [123]
And in this case you can’t use the defense that the elements aren’t supposed to be mutable because the whole point of arrays is that you can change their contents. So if we’re going to do default filling, then if we use the same default instance for the whole array we’ve got this pretty bad footgun.
What if we instead we used a newly constructed default value for each array slot? That would be equivalent to doing [ zero(T) for i=1:n, j=1:n ] in place of Matrix{T}(undef, n, n). Let’s time that as well for T = BigInt, n = 3:
julia> @belapsed [ zero(BigInt) for i=1:3, j=1:3 ]
1.7863372859025033e-7
Ah, now we get a 11x slowdown relative to Matrix{T}(undef, n, n) instead of just 2.5x. On the other hand, we don’t have the issue with contagious mutation:
julia> A = [ zero(BigInt) for i=1:3, j=1:3 ]
3×3 Matrix{BigInt}:
0 0 0
0 0 0
0 0 0
julia> add!!(A[1], 123)
123
julia> A
3×3 Matrix{BigInt}:
123 0 0
0 0 0
0 0 0
What’s worse about this one is that, whereas the relative slowdown of zeros(BigInt, n, n) versus Matrix{T}(undef, n, n) stays around 2.5x slower even for larger matrices, [ zero(T) for i=1:n, j=1:n ] is gets relatively worse and worse the bigger n gets—for n = 5 it’s already 23x slower; for n = 100 it’s 89x slower.
So filling arrays with a default value is pretty bad from a performance perspective, at least in the case where you just need an array that you’re immediately going to replace the contents of. It’s bad enough that it really can’t be the only option—there needs to be something else.
The proposed alternative is using an iterator to generate initial values for an array. That’s fine when it works, but how would you implement the pascal function as an iterator? This is a serious question. There isn’t an elegant way to do it that I can tell because the initialization wants to look at the previously initialized values, which iteration won’t let you do. Assuming the iterator goes in column major order you could remember what you assigned in the previous column, but where do you store those values? You can’t look them up in the array you’re initializing because that array object doesn’t exist yet—you’re still constructing it. So you’d need to allocate an O(n) temporary storage buffer to remember the values you generated in the last column of the array.
Keep in mind that the Pascal triangle initialization is hardly the gnarliest initialization possible. It’s not too hard to come up with reasonable cases where either the best initialization order for an array is all over the place or each next element depends on the previously initialized elements in some arbitrarily complex way so that you’d need to keep a second copy of the entire array so far.
Regarding C++: you’ll note that people have asked how to get around this exact limitation of std::vector for performance reasons. In C++, of course, you always have other options that allow uninitialiezed memory, which is what some answers there suggest. In Julia we wouldn’t have that option because unlike std::vector, which builds upon lower level C++ primitives like raw pointers, Julia’s Array is the primitive mutable collection time and there’s nothing lower to resort to. In other words, we can’t not have a way to do this efficiently.
14 Likes
This is sufficiently tricky that it took me a while to figure out how to do it. Here’s the code:
function pascal′(::Type{T}, n::Integer) where {T<:Real}
x = one(T)
v = zeros(T, n)
[ begin
if i == 1 && j != 1
x = zero(T)
end
x = v[i] += x
end for i = 1:n, j = 1:n ]
end
This works but it seems significantly less clear than the original, not to mention less efficeint because it needs a temporary vector v to keep track of what was generated in the previous column.
Oh, I should also mention that we’ve been focused here on arrays but partially initialized objects are also a place where undef happens. The canonical motivating example is a recursive type like this:
mutable struct X
f::X
function X()
x = new() # x.f is undef here
x.f = x
end
end
If all objects have to be constructed fully initialized, how can you do this? Note that this isn’t entirely academic since this is basically the structure of a circular linked list. I should also note, however, that unlike arrays, the compiler can and does analyze on an individual field basis, whether each field can be undef or not, so if a type doesn’t ever construct partially initialized objects then the compiler will know that there’s no need for undef checks.
4 Likes
Wow, thanks for this very detailed post! Lots to unpack and to digest! (for a future day with more time)
I have a few direct reactions:
I like how the example illustrates that code written for, e.g. some numerical computation, should work well with both value types and reference types. I had not really though about this.
I fully agree that a default instance/value is a bad idea.
You compare the cost of initialization when going from “no initialization” to “initialization”, of course this can look bad (yes, I know “no initialization” requires “null initialization”, but still). Luckily, I think, initialization is usually not the dominant cost in a computation - and if it is, then that’s a pretty good sign you are using a dense structure when you should be using a sparse.
The Pascal example connects back to earlier discussion in this thread about incremental construction and sizehint!
function pascal(::Type{T}, n::Integer) where {T<:Real}
A = T[]
sizehint!(A, n*n)
idx(i, j) = n * (j - 1) + i
for j = 1:n
for i = 1:n
if i == 1 || j == 1
push!(A, one(T))
else
push!(A, A[idx(i-1, j)] + A[idx(i, j-1)])
end
end
end
return reshape(A, (n,n))
end
Yes, your index-style implementation looks much better than this one.
I understand your point, and I very much appreciate the emphasis on performance in the design. However, this doesn’t appear to be an absolute. Even Julia have tradeoffs here (I believe).
For example, going back to the Pascal calculation. Assume you wanted to use high precision floats and opted for DoubeFloats.jl. Now, since Double64 is not an isbitstype type, you have your “number objects” heap allocated [this is false, see below] – which is probably not ideal for performance in this example. And AFAICT, there is no way to tell Julia that a vector should store elements in-place, or to annotate a mutable struct (or other non-isbitstype struct) to be a value type instead of a reference type (which would cause it to be stored in-place in a vector).
This change, allocating the objects in-place, I think, would be a standard way to optimize the code in C++ (without changing algorithm). It could then be further tweaked to use a tiling iteration pattern to optimize the cache usage. I imagine (without having checked this) that this would give a good speedup.
This is not saying that “C++ is better than Julia”: C++ is complicated, Julia is simpler. Simplicity is also valuable. (IMHO, of course)
EDIT: My statement above about Double64 being heap allocated is false. Double64 actually an isbitstype type. I should have checked. Apologies. (However, the argument still holds for types that are not isbitstype, I believe.)
EDIT:
Is this a typo? isbitstype(Float64) == true.
Are there any isbitstype elements that aren’t stored in-place?
Well, mutable struct instances are referenced to implement mutability, so they’re all reference types. As I understand it, isbitstypes are equivalent to value types, but the term formally refers to Val in Julia instead.
1 Like
Somewhat relevant, escape analysis, that’s being worked on for Julia, could enable some optimizations that you’d like:
https://docs.julialang.org/en/v1.11-dev/devdocs/EscapeAnalysis/
stack allocation of mutable objects
1 Like
This change is not trivial to do. Consider this, which must continue to work:
mutable struct A
el::Int
end
a = A(1)
b = [ a ]
b[1].el = 2
a.el == 2 # true
Evidently, trying to declare b to have its elements stored in-line doesn’t work with the non-moving GC we have today; a already exists before b is ever allocated, so either the in-line semantics of b breaks, the use of a in b copies (thus breaking the == test) or an error is thrown. Neither of these seem good.
That Julia stores mutable values as references is by design; Julia doesn’t have the reference semantics exposed to be manipulated/worked with by the user. There are undoubtly situations where treating a collection and the objects it contains as one entity for the purposes of memory management would be useful, but that’s not the semantics Vector (or mutable objects for that matter) has today.
That being said, if you do want to have those kinds of semantics, you can somewhat get to them using Ref (pointing it to an array like Ref([ a ], 1)), which must explicitly be indexed in order to retrieve the stored (mutable) object.
2 Likes