optimization for associative reductions

[Alexander-Barth]

This PR addresses #188

[rafaqz]

I think we also need the version with a function as first argument.

And also is count really in that list? Especially with the function? I thought it would need some tweaking.

[meggart]

I think the root problem is that some internal dispatch changed in Julia Base. All of these were working in previous Julia versions because they used to fall back to reduce or foldl for which we have an efficient implementation. So I am wondering if we could rather make this work again instead of hand-picking a few reduction functions that can be fixed and maybe missing some on the way.

[rafaqz]

Oh right the dispatch on internals may have broken, sorry I should have realised that. We also just switched this to the main method in DimensionalData.jl

https://github.com/meggart/DiskArrays.jl/blob/ccd3092ca6f32c2663ded1d74176aa9dd222a5ef/src/mapreduce.jl#L7-L12

Sorry @Alexander-Barth for the misdirection! I guess we just need to edit the mapreduce implementation instead. Its basically the generalisation of the methods in this PR.

But probably we can keep the tests to catch these things in future, especially if we explicitly test the access count of those arrays after the reductions, to make sure we are not hitting fallbacks. We clearely are not doing that at the moment.

[meggart]

But probably we can keep the tests to catch these things in future, especially if we explicitly test the access count of those arrays after the reductions, to make sure we are not hitting fallbacks. We clearely are not doing that at the moment.

I am not sure about this we have these tests: https://github.com/meggart/DiskArrays.jl/blob/ccd3092ca6f32c2663ded1d74176aa9dd222a5ef/test/runtests.jl#L114 which also test the access count. Also the test that was added in this PR already passes the access_count test on the current main:

using DiskArrays
using DiskArrays.TestTypes
using Test
A = rand(1:10,30,30)
DA = AccessCountDiskArray(A,chunksize=(2,2));

r = sum(DA);

@test getindex_count(DA)==length(eachchunk(DA))

Here we test that every chunk is accessed exactly once, please try this on main, so I would guess the slow performance in #188 must come from somewhere else or some operation we are not yet hitting with our tests.

[meggart]

Ok, I did some digging and realized that currently we do not hit our mapreduce methods but rather the iterator fallback, which does read the data chunk by chunk but is more inefficient and numerically unstable compared to this implementation. Will launch Cthulhu to find out where this goes wrong

[meggart]

Ok, I just opened #191 which also solves the performance issues described in #188. My suggestion for this PR would be to merge anyway, because of the following benefits:

1. The PR adds some more safety in case Julia Base switches their mapreducedim internal signatures again in future Julia versions.
2. In the case of sum for very large or very small chunks this PR should also be more accurate than the fallback mapreduce, i.e. when already hitting floating point issues within a chunk or when adding too many chunk in Float32 and Float16 precision (would be great to have some examples where this happens in the unit tests, I could also do this later if you don't have time)
3. WIth this PR any and all can be much faster because of early-loop-exit when the first true resp false value is encountered which we would not have in the generic fallback

So I would definitely want to add the specialized functions from this PR for sum, prod, any and all, ideally with some accompanying test that check these advantages. For count, maximum, and minimum I am quite neutral as the only benefit would be protection from future Base interface changes, would be great to hear your opinions.

[rafaqz]

Seems this has been pretty productive! good to have all those things working and tested.

I think youre right - sum, prod, any and all are the real beneficiaries of not using mapreduce. maximum and minimum are probably worth keeping too for reliability.

I don't quite understand how count is working over chunks here because it returns an Int but does not work on an iterator of Int. I thought it would need to be sum as the outer function and count as the inner.

[meggart]

Superseded by #196