Julia: implement zeros() by calling calloc

[StefanKarpinski]

There's a clever trick that we could use to create large zero matrices really fast: mmap the file /dev/zero. This is, in fact, exactly what this "file" exists for. The benefit of doing this are:

1. It's nearly instantaneous since no memory actually needs to be allocated or filled with zeros until it's accessed.
2. You can read and write the memory exactly as you normally would: the kernel only allocates memory pages for you when you do something with them.

Since a fair amount of the time, no one actually touches most of the memory in a zeros array, this might be a big win. On the other hand, the drawbacks are:

1. Trade obvious immediate allocation cost for unobvious delayed allocation cost.
2. Can run out of memory on read/write instead of only on allocate.
3. Doesn't work for anything but zeros(), e.g. for ones().

[ViralBShah]

The thing is that we use Array() as the default constructor almost everywhere now, and zeros is not used nearly as much.

Wouldn't it be better if we can handle the pagefaults ourselves, which will take care of all cases that call fill()?

[StefanKarpinski]

That could be possible, but testing with an mmap of /dev/zero might be an easy way to find out what there is to gain.

[ViralBShah]

It will definitely lead to better cache behaviour, and there will be a certain pattern of usage where this will be greatly beneficial for sure.

-viral

On Jul 18, 2011, at 11:23 AM, StefanKarpinski wrote:

That could be possible, but testing with an mmap of /dev/zero might be an easy way to find out what there is to gain.

Reply to this email directly or view it on GitHub:
https://github.com/JuliaLang/julia/issues/130#issuecomment-1593014

[ViralBShah]

See:

http://stackoverflow.com/questions/2688466/why-mallocmemset-slower-than-calloc

[ViralBShah]

It is quite possible that this does not give any observable gains for anything except very large matrices. Is it possible to quickly do an experiment to see if this gives any measurable benefits? Also, we do not use zeros() much in our codebase.

[StefanKarpinski]

Yeah, this is a cool idea, but the way our memory allocation works via a memory pool, it's not very practical. Let's close for now.

[StefanKarpinski]

This is a real performance issue that I've seen in the wild recently and has come up on the mailing list:

https://groups.google.com/forum/#!topic/julia-users/aW4rjUIFq6w

I'm fairly certain that NumPy is using the mmap /dev/zero trick and we should too.

[stevengj]

Why not just use calloc? See also #9147.

[ViralBShah]

calloc does seem like the easiest first thing to try.

[quinnj]

It'd be worth bench-marking Anonymous Mmaps too, since they're now fully supported and "easy" :)

julia> m = Mmap.mmap(Vector{Float64}, 10000)
10000-element Array{Float64,1}:
 0.0
 0.0
 0.0

[StefanKarpinski]

Why not just use calloc?

Yes, excellent point, @stevengj. That's definitely the first thing to do.

[StefanKarpinski]

No reason not to do this. (Also, a really easy three-digit issue.)

[danielmatz]

I used BenchmarkTools to compare Base.zeros against calloc and Mmap.mmap implementations. I've posted the script with my zeros_calloc and zeros_mmap functions and the benchmarks as a gist.

I first compared creating a small vector of zeros. I found that Base.zeros does the best here.

| Function | Median Time |
| --- | --- |
| Base.zeros | 43 ns |
| zeros_calloc | 660 ns |
| zeros_mmap | 1.47 μs |

Then I compared the creation of a somewhat large array. In this case, the Mmap.mmap solution clearly wins.

| Function | Median Time |
| --- | --- |
| Base.zeros | 2.28 ms |
| zeros_calloc | 2.17 ms |
| zeros_mmap | 1.98 μs |

Lastly, I created a somewhat large array of zeros, and then overwrote every element by filling with ones. And now the calloc version wins.

| Function | Median Time |
| --- | --- |
| Base.zeros | 3.50 ms |
| zeros_calloc | 2.83 ms |
| zeros_mmap | 4.74 ms |

Assuming I don't have an error somewhere... Which implementation should we go with?

[quinnj]

I don't think the Mmap.mmap version is viable since it leads to the "shared data" situation and can't be mutated. I guess maybe in a non-vector situation we could consider it. Here's what I'm talking about:

julia> t = Mmap.mmap(Vector{UInt8}, 10)
10-element Array{UInt8,1}:
 0x00
 0x00
 0x00
 0x00
 0x00
 0x00
 0x00
 0x00
 0x00
 0x00

julia> push!(t, 0x00)
ERROR: cannot resize array with shared data
 in push!(::Array{UInt8,1}, ::UInt8) at ./array.jl:480
 in push!(::Array{UInt8,1}, ::UInt8) at /Users/jacobquinn/julia/usr/lib/julia/sys.dylib:?

[yuyichao]

Also note that we can't trivially use calloc for malloc array either since that provide less alignment that what we need.

[StefanKarpinski]

I think we should use manual zeroing for small arrays and calloc (with appropriate alignment) for large arrays. I'm guessing the cutoff should be about a memory page but it would be good to test that experimentally – which you're already doing. Thanks for investigating and pushing this forward!

[stevengj]

https://github.com/stevengj/julia/commit/7c9ab1d88a2b4950a9b889f6b20c522ca50a090c shows how to get 16-byte alignment from calloc.

[yuyichao]

16 bytes is fine. We do 64bytes aligned allocation though.

[stevengj]

The same code can easily be modified to do 64-byte aligned calloc. (Though I thought there was no performance advantage for more than 32-byte alignment, even with AVX?)

[danielmatz]

@stevengj Thanks for the link to your calloc_a16 code, that was very helpful.

I think I'll need to add information to jl_array_flags_t to track whether calloc_a16 was used, so that I can know to call free_a16. I see that that structure is very carefully packed... Can I take a bit away from ndims?

[yuyichao]

You should use the same free.

[danielmatz]

@yuyichao I don't understand how that will work. Steven's calloc_a16 requests 16 more bytes from calloc than the user requested. Then he returns a pointer that is aligned to a 16 byte boundary, but that pointer is never the one that calloc itself returned. Isn't it true that I have to call free on the pointer originally returned by calloc? Steven provides a free_a16 function to do this. That's what led me to believe that I need to use a bit to track whether the memory came from calloc_a16.

[yuyichao]

I mean you should just optimize (for size) and use the existing implementation in jl_{malloc,calloc,free,realloc} when allocating array buffers using malloc. You can use the low bits in the size storage for alignments. This shouldn't cause memory wasting since libc needs to do the over allocation too (for alignment). (Note that you only need to over allocate 64bytes and you are guaranteed to have a sizeof(void*) available to store the size).

You do need to keep track of jl_ptr_to_array(and similar functions) that transfers the ownership of the memory since those should be free'd in using free but that can be done in jl_gc_track_malloc_array and you have at least 4 bits there to mark it.

[danielmatz]

Thank you for the help! A couple of questions...

1. Are you envisioning that I use jl_calloc to allocate the array data only, and then use jl_ptr_to_array to turn it into an array? Or are you thinking I should make a new version of jl_new_array that uses jl_calloc to allocate the the jl_array_t object and the data inline at once?
2. You mention that there are at least 4 bits available in jl_gc_track_malloced_array. I'm not sure what you mean. What structure are you talking about?

[yuyichao]

Are you envisioning that I use jl_calloc to allocate the array data only

Well, also jl_malloc and only for malloc arrays (not the ones we use different array storage) and possibly a different form of it that use 64bytes alignments and leave the 16bytes default alignment unchanged.

and then use jl_ptr_to_array to turn it into an array?

No. jl_ptr_to_array should still accept libc (freeable) pointers directly.

Or are you thinking I should make a new version of jl_new_array that uses jl_calloc to allocate the the jl_array_t object and the data inline at once?

A new version of jl_new_array (and possibly optimized 1d/2d versions) that use jl_calloc when we use (counted wrapper of) malloc now and also change the current usage of (counted wrapper of) malloc to jl_malloc (or the 64bytes aligned version as mentioned above).

You mention that there are at least 4 bits available in jl_gc_track_malloced_array. I'm not sure what you mean. What structure are you talking about?

mallocarray_t::a is guaranteed to be 16 bytes aligned and mallocarray_t::next is guaranteed to be 8 bytes aligned (or at least 4) so you actually have at least 6 bits to store additional informations.

Also for the layout of jl_calloc or jl_malloc return value, we need at least 16bytes alignment for all allocations through those (might even be 64bytes for the one used by array) so the size should be 16bytes aligned too. This means that you have at least 4 bits there to store the offset. Since libc should guarantee 4bytes alignment, the 4 bits here is enough to represent any alignment offset for 64bytes alignement. Even if it doesn't, say if we decide that we want to support 1024 bytes alignment for whatever reason, for the case where you have large offset, you automatically have much more space before this slot to store additional information about the offset. i.e.

<malloc return address>[padding for alignment][additional alignment offset if >=16][sizeof(void*) to store size and alignment offset < 16]<return value 16/64bytes aligned>[data][overallocation for alignment]

[stevengj]

Note that libc guarantees 8-byte alignment from malloc already on all extant 32-bit platforms. On 64-bit systems with Windows, Mac, and Linux/glibc, malloc and calloc are already 16-byte aligned. So, if we just want 16-byte alignment, we can use plain calloc and free on 64-bit systems.

[yuyichao]

For the sizes we want to use malloc/calloc we all want 64bytes alignment.

[danielmatz]

You say you want me to use jl_malloc. In the current array code, I see calls to jl_gc_alloc and jl_gc_managed_malloc. Can you give me more context on why you want me to switch to using jl_malloc? What are the differences between all of these functions? It seems like jl_malloc allocates memory with an extra 16 bytes at the front but without any special alignment. And I see that jl_gc_alloc allocates memory with a tag at the front and uses either a pool or allocates a "big" chunk of 64 byte aligned memory. And then jl_gc_managed_malloc allocates 64 byte aligned memory, but with no extra bytes at the front.

Also, looking at the current implementation of _new_array_, assuming I understand these allocation functions a little bit, it sure seems like it is already doing 16 byte alignment for small arrays and 64 byte alignment for larger arrays. This makes me think I'm not understanding what you were saying about doing 64 byte alignment for arrays. Maybe after I get some context on the switch to jl_malloc this part will make more sense, too...

Thank you again for the help!

[yuyichao]

You say you want me to use jl_malloc.

To replace the current use of malloc wrappers, which is jl_gc_managed_malloc

Can you give me more context on why you want me to switch to using jl_malloc?

Because that means you don't have to distinguish between malloc and calloc array when freeing/resizing them (ref https://github.com/JuliaLang/julia/issues/130#issuecomment-257057036) and it shouldn't waste more memory then what we currently do (ref https://github.com/JuliaLang/julia/issues/130#issuecomment-257064226)

What are the differences between all of these functions?
It seems like jl_malloc allocates memory with an extra 16 bytes at the front but without any special alignment.
And then jl_gc_managed_malloc allocates 64 byte aligned memory, but with no extra bytes at the front.

jl_malloc and jl_gc_managed_malloc are malloc wrappers. jl_malloc allocates a sizeof(void*) to store the size so that jl_realloc can be implemented. jl_gc_managed_malloc allocates 64bytes aligned memory for array.

And I see that jl_gc_alloc allocates memory with a tag at the front and uses either a pool or allocates a "big" chunk of 64 byte aligned memory.

jl_gc_alloc allocates managed (tagged) memory and is irrelevant here.

Also, looking at the current implementation of _new_array_, assuming I understand these allocation functions a little bit, it sure seems like it is already doing 16 byte alignment for small arrays and 64 byte alignment for larger arrays.

Correct, and you should not change the small array to use calloc, only the big ones that are using the malloc wrappers, which should all be 64bytes aligned IIRC.

[stevengj]

@yuyichao, for large-size arrays, don't malloc/calloc ordinarily give page-aligned data?

My understanding was that we only need 16-byte alignment for correctness (since SIMD operations may fail on data with less alignment), whereas 64-byte alignment is merely an optimization (some SIMD operations are faster). In which case we can rely on the OS to give us 64-byte alignment most of the time for large arrays, without worrying whether it is a strict guarantee, no?

[yuyichao]

for large-size arrays, don't malloc/calloc ordinarily give page-aligned data?

Not AFAICT. libc still need to store metadata somewhere.

My understanding was that we only need 16-byte alignment for correctness

We aren't specify that anywhere and it is not needed for correctness on 32bit. (edit: and it is only needed for correctness on 64bits platforms for Int128 or similar 128bits C types)

In which case we can rely on the OS to give us 64-byte alignment most of the time for large arrays, without worrying whether it is a strict guarantee, no?

The issue is that OS does not give 64bytes aligned pointer most of the time last time I checked. Not even 32bytes aligned and this causes >40% timing fluctuation last time I checked.

[yuyichao]

In which case we can rely on the OS to give us 64-byte alignment most of the time for large arrays, without worrying whether it is a strict guarantee, no?

In fact, I've just checked using allocations of 1MB, which should be roughly the most relevant size since the alignment doesn't seem to matter too much when the data doesn't fit in the L3 cache anymore. My glibc gives a pretty even distribution of alignments (256 ones for 0, 16, 32, 48 mod 64 each.)

[stevengj]

On MacOS, large memory allocations are guaranteed to be page-aligned. (I just tried it, and allocations as small as 1024 bytes seem to be consistently at least 256-byte aligned.)

[yuyichao]

I'm certainly fine with a different implementation that relies on malloc alignments on mac.

[danielmatz]

Thank you both for your patience and for taking the time to work through this with me. I think I'm well on my way.

I understand now that we want to have wrappers around jl_malloc and jl_calloc that do alignment manually, so that we don't have to remember which one was used.

But this leads me to a few questions:

1. Do we want to even keep jl_malloc_aligned around (the version that calls posix_memalign or something similar)? If we keep it, we'll have two flavors of aligned malloc. I only see it used in the array code and a few times in the threading code. I'll of course be replacing its usage in the array code, but can I just replace its usage in the threading code, too?
2. Can I go ahead and remove jl_gc_managed_malloc and jl_gc_managed_realloc, then? They don't seem to be used anywhere else.
3. The isaligned field in jl_array_flags_t won't have a purpose once I am no longer calling jl_malloc_aligned from the array code. Can I remove that? Any thoughts on what to do with that bit?

[yuyichao]

Do we want to even keep jl_malloc_aligned around.
I'll of course be replacing its usage in the array code, but can I just replace its usage in the threading code, too?

We should keep it.

Can I go ahead and remove jl_gc_managed_malloc and jl_gc_managed_realloc, then? They don't seem to be used anywhere else.

Yes. Note that you might want to replace the realloc in the array code with something that pass in the old size explicitly so that you can do less memcpy

The isaligned field in jl_array_flags_t won't have a purpose once I am no longer calling jl_malloc_aligned from the array code. Can I remove that? Any thoughts on what to do with that bit?

It can be removed. Leave the bit empty (name it _dummy or sth similar).

[yuyichao]

Note that you might want to replace the realloc in the array code with something that pass in the old size explicitly so that you can do less memcpy

i.e. You might be able to reuse the name and the signature and keep the jl_gc_managed_* functions if the new functions is compatible in semantics. (not strictly necessary).

[danielmatz]

Out of curiosity, why should we keep jl_malloc_aligned?

Yes. Note that you might want to replace the realloc in the array code with something that pass in the old size explicitly so that you can do less memcpy

I'm not following. How would passing in the size allow me to do less memcpy? Other than the copying realloc does itself, I don't see any memcpy going on...

[yuyichao]

Out of curiosity, why should we keep jl_malloc_aligned?

Useful for thing that mustn't trigger the GC.

How would passing in the size allow me to do less memcpy?

Because the array code knows the exact size it needs, which is in general smaller than the allocation size.

[danielmatz]

Useful for thing that mustn't trigger the GC.

Oh, yes, I see that now.

And looking at the isaligned field more closely, I see that I can't get rid of it, because of jl_ptr_to_array.

[yuyichao]

And looking at the isaligned field more closely, I see that I can't get rid of it, because of jl_ptr_to_array.

jl_ptr_to_array can't be resized AFAICT.

[yuyichao]

jl_ptr_to_array can't be resized AFAICT.

Which is also why I said you can move that info into mallocarray_t instead.

[danielmatz]

jl_ptr_to_array can't be resized AFAICT.

Even if jl_ptr_to_array is called with the own_buffer option? I don't see how we can tell the difference once flags.how is set to 2...

[yuyichao]

Even if jl_ptr_to_array is called with the own_buffer option?

Correct.

I don't see how we can tell the difference once flags.how is set to 2...

isshared. Array that owns shared data can't be resized.

[danielmatz]

isshared. Array that owns shared data can't be resized.

Ah, yes, thanks.

[danielmatz]

Because the array code knows the exact size it needs, which is in general smaller than the allocation size.

I guess I was only considering the case where an array fills its allocated space, and then has to increase its allocation. In that case, realloc would already only be copying the minimum amount of data. But I guess you are saying that it is common to resize the array's allocation while still only using a subset of the allocation?

[danielmatz]

I think my changes are working properly, but when I make the final changes to the zeros function and then run make to build the sysimg, it gets hung up. Any tips on how to debug that?

[yuyichao]

Attach gdb and see why it hangs.

[danielmatz]

The backtrace is tens of thousands of frames deep. The top of the stack is thousands of calls to inst_tuple_w_. Below that are thousands of calls to jl_apply_generic. Is that normal?

[yuyichao]

It's normal if you messed sth up like memory allocations. There're many ways I'd try to debug it including running in rr, dump local variables, print in allocation etc. It's impossible to tell without seeing the code.

[danielmatz]

My apologies for setting this issue aside for so long. I should have time again now to continue looking at it.

I've been testing my changes again. It seems that they work fine everywhere except for const GLOBAL_RNG = MersenneTwister() in random.jl. That is where things hang when I run make. If I comment out that line, the build completes (though the new Julia executable doesn't work, since this global variable is missing).

Is there some interaction between memory allocation, globals, and the pre-compilation steps that I am missing?

If I can't figure this out soon, would it be OK to submit a WIP pull request? Maybe someone else would quickly spot my mistake.

[yuyichao]

Submitting a wip pr or at least have a pointer to the wip code would be useful.

[danielmatz]

I think I figured out my issue. I had replaced the zeros methods for all eltypes, but there are some where zeroing the memory is not appropriate. I've now restored the fill! implementation as the fallback, and will only opt in to using calloc when appropriate for the eltype.

So I should be back on track now. I have the WIP PR up, too.

[antoine-levitt]

Could someone comment on the status of this? There were a couple of followup PRs but it doesn't seem that this functionality is in. Was it decided it was not beneficial? Or is it but nobody got around to it?

[StefanKarpinski]

I think nobody got around to it.

[ViralBShah]

I think the latest on the topic was https://github.com/JuliaLang/julia/pull/22953