# Inconsistency in \`accumulate\` between \`Array\` and \`CuArray.\`

**URL:** <https://discourse.julialang.org/t/inconsistency-in-accumulate-between-array-and-cuarray/127382>\
**Category:** General Usage\
**Tags:** gpu, gpuarrays, cuda\
**Created:** [March 26, 2025, 2:04pm UTC](https://discourse.julialang.org/t/inconsistency-in-accumulate-between-array-and-cuarray/127382 "2025-03-26T14:04:36Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![0samuraiE](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/0samuraie/32/209825_2.png) [@0samuraiE](https://discourse.julialang.org/u/0samuraiE)\
**Post date:** [March 26, 2025, 2:04pm UTC](https://discourse.julialang.org/t/inconsistency-in-accumulate-between-array-and-cuarray/127382/1 "2025-03-26T14:04:36Z")

</div>

Hi all.  
I’m performing the following computation using `GPUArrays`:

```julia-repl
julia> using GPUArrays

julia> src = CuArray(rand(0:9, 10))
10-element CuArray{Int64, 1, CUDA.DeviceMemory}:
 3
 4
 3
 7
 3
 9
 7
 2
 9
 6

julia> dst = similar(src);

julia> op = (acc, x) -> acc + (1 - ((x >> 0) & 0x1))
#211 (generic function with 1 method)

julia> GPUArrays.neutral_element(::typeof(op), T) = one(T)

julia> accumulate!(op, dst, src; dims=1)
10-element CuArray{Int64, 1, CUDA.DeviceMemory}:
 1
 2
 2
 2
 2
 2
 2
 3
 3
 2

julia> src = collect(src)
10-element Vector{Int64}:
 3
 4
 3
 7
 3
 9
 7
 2
 9
 6

julia> dst = collect(dst);

julia> accumulate!(op, dst, src; init=0)
10-element Vector{Int64}:
 0
 1
 1
 1
 1
 1
 1
 2
 2
 3

```

I noticed that the results differ between CPU and GPU execution.  
Why is that the case? Also, could someone explain what the “neutral element” means in this context?

---

<div class="post-metadata">

**Author:** ![eldee](https://avatars.discourse-cdn.com/v4/letter/e/b5a626/32.png) [@eldee](https://discourse.julialang.org/u/eldee)\
**Post date:** [March 26, 2025, 5:11pm UTC](https://discourse.julialang.org/t/inconsistency-in-accumulate-between-array-and-cuarray/127382/2 "2025-03-26T17:11:39Z")

</div>

My guess is that the issue is that `op` is non-associative (and non-commutative),

```julia-repl
julia> op(10, op(20, 30))
10

julia> op(op(10, 20), 30)
12

julia> op(2, 3)
2

julia> op(3, 2)
4

```

and the GPU version of `accumulate!` will change the order of operations to achieve parallelisation. E.g. it will compute, using infix notation, `((a op b) op c) op d` as `(a op b) op (c op d)`.

* * *

The neutral element for an operation is just some value `n` such that `op(x, n) == x == op(n, x)` for all `x`. In this case `n = 1` indeed is appropriate. The parallel algorithm will presumably always (to avoid branching) combine two elements `x` and `y`, even when no appropriate `y` exists. In that case it will use `y = n`, so that the outcome `op(x, y) == x` is as if you didn’t apply `op` at all.

---

<div class="post-metadata">

**Author:** ![0samuraiE](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/0samuraie/32/209825_2.png) [@0samuraiE](https://discourse.julialang.org/u/0samuraiE)\
**Post date:** [March 26, 2025, 11:58pm UTC](https://discourse.julialang.org/t/inconsistency-in-accumulate-between-array-and-cuarray/127382/3 "2025-03-26T23:58:38Z")

</div>

Thank you for your clear answer.