# How to use OffsetArray with CUDA

**URL:** <https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400>\
**Category:** General Usage\
**Tags:** cuda, offsetarrays\
**Created:** [August 20, 2024, 2:07am UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400 "2024-08-20T02:07:25Z")\
**Posts on this page:** 9\
**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:** [August 20, 2024, 2:07am UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/1 "2024-08-20T02:07:25Z")

</div>

Hello all.  
In the following code, sum returns `ERROR: Scalar indexing is disallowed`.  
How should `sum` be dispatched or How should I do?

```julia
A = OffsetArray(CUDA.rand(10), -4:5);
a = @views A[-1:2];
sum(a)

```

---

<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:** [August 20, 2024, 3:01am UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/2 "2024-08-20T03:01:30Z")

</div>

A simple implementation is shown. But I think there must be a better implementation.

```julia
Base.sum(A::OffsetArray{T,N,CuArray{T,N,M}}) where {T,N,M} = sum(parent(A))
function Base.sum(A::SubArray{T,N,OffsetArray{T,N,CuArray{T,N,M}}}) where {T,N,M}
    indices = A.indices
    offsets = parent(A).offsets
    sum(@view parent(parent(A))[CartesianIndices(ntuple(n -> indices[n] .- offsets[n], N))])
end

```

---

<div class="post-metadata">

**Author:** ![roflmaostc](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/roflmaostc/32/30123_2.png) [@roflmaostc](https://discourse.julialang.org/u/roflmaostc)\
**Post date:** [August 20, 2024, 11:10am UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/3 "2024-08-20T11:10:12Z")

</div>

I think OffsetArray is not quite compatible with CUDA naively, even broadcasting fails:

```julia
julia> A = OffsetArray(CUDA.rand(10), -4:5);

julia> sum(A)
ERROR: Scalar indexing is disallowed.
Invocation of getindex resulted in scalar indexing of a GPU array.
This is typically caused by calling an iterating implementation of a method.
Such implementations *do not* execute on the GPU, but very slowly on the CPU,
and therefore should be avoided.

If you want to allow scalar iteration, use `allowscalar` or `@allowscalar`
to enable scalar iteration globally or for the operations in question.
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] errorscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:155
  [3] _assertscalar(op::String, behavior::GPUArraysCore.ScalarIndexing)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:128
  [4] assertscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:116
  [5] getindex
    @ ~/.julia/packages/GPUArrays/bbZD0/src/host/indexing.jl:50 [inlined]
  [6] getindex
    @ ~/.julia/packages/OffsetArrays/hwmnB/src/OffsetArrays.jl:438 [inlined]
  [7] _mapreduce(f::typeof(identity), op::typeof(Base.add_sum), ::IndexLinear, A::OffsetVector{Float32, CuArray{…}})
    @ Base ./reduce.jl:438
  [8] _mapreduce_dim
    @ ./reducedim.jl:365 [inlined]
  [9] mapreduce
    @ ./reducedim.jl:357 [inlined]
 [10] _sum
    @ ./reducedim.jl:1015 [inlined]
 [11] _sum
    @ ./reducedim.jl:1014 [inlined]
 [12] sum(a::OffsetVector{Float32, CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}})
    @ Base ./reducedim.jl:1010
 [13] top-level scope
    @ REPL[15]:1
 [14] top-level scope
    @ ~/.julia/packages/CUDA/htRwP/src/initialization.jl:206
Some type information was truncated. Use `show(err)` to see complete types.

julia> A .+ A
ERROR: Scalar indexing is disallowed.
Invocation of getindex resulted in scalar indexing of a GPU array.
This is typically caused by calling an iterating implementation of a method.
Such implementations *do not* execute on the GPU, but very slowly on the CPU,
and therefore should be avoided.

If you want to allow scalar iteration, use `allowscalar` or `@allowscalar`
to enable scalar iteration globally or for the operations in question.
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] errorscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:155
  [3] _assertscalar(op::String, behavior::GPUArraysCore.ScalarIndexing)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:128
  [4] assertscalar(op::String)
    @ GPUArraysCore ~/.julia/packages/GPUArraysCore/GMsgk/src/GPUArraysCore.jl:116
  [5] getindex
    @ ~/.julia/packages/GPUArrays/bbZD0/src/host/indexing.jl:50 [inlined]
  [6] getindex
    @ ~/.julia/packages/OffsetArrays/hwmnB/src/OffsetArrays.jl:438 [inlined]
  [7] _broadcast_getindex
    @ ./broadcast.jl:675 [inlined]
  [8] _getindex
    @ ./broadcast.jl:705 [inlined]
  [9] _broadcast_getindex
    @ ./broadcast.jl:681 [inlined]
 [10] getindex
    @ ./broadcast.jl:636 [inlined]
 [11] macro expansion
    @ ./broadcast.jl:1004 [inlined]
 [12] macro expansion
    @ ./simdloop.jl:77 [inlined]
 [13] copyto!
    @ ./broadcast.jl:1003 [inlined]
 [14] copyto!
    @ ./broadcast.jl:956 [inlined]
 [15] copy
    @ ./broadcast.jl:928 [inlined]
 [16] materialize(bc::Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{…}, Nothing, typeof(+), Tuple{…}})
    @ Base.Broadcast ./broadcast.jl:903
 [17] top-level scope
    @ REPL[16]:1
 [18] top-level scope
    @ ~/.julia/packages/CUDA/htRwP/src/initialization.jl:206
Some type information was truncated. Use `show(err)` to see complete types.

julia> A
10-element OffsetArray(::CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, -4:5) with eltype Float32 with indices -4:5:
 0.19389287
 0.71536046
 0.8175033
 0.57097757
 0.058391646
 0.72871023
 0.6042697
 0.88648033
 0.76349247
 0.9775343

```

There is some hints here if it might help:

> <https://github.com/JuliaArrays/OffsetArrays.jl/pull/57>
>
> Adapt is a lightweight dependency that allows wrapper packages like OffsetArray …to be
> converted between different underlying arrays. The primary use-case is to be able to
> convert a CPU OffsetArray to a OffsetArray whose memory is a \`CuArray\`.

---

<div class="post-metadata">

**Author:** ![maleadt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maleadt/32/10097_2.png) [@maleadt](https://discourse.julialang.org/u/maleadt)\
**Post date:** [August 20, 2024, 11:49am UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/4 "2024-08-20T11:49:53Z")

</div>

Yeah, Julia isn’t currently great wrt. wrapped arrays and preserving functionality from the contained array type where needed. I typically link to [Use with multiple wrappers · Issue #21 · JuliaGPU/Adapt.jl · GitHub](https://github.com/JuliaGPU/Adapt.jl/issues/21) for this, and this would need some work in Base to resolve (e.g., [AbstractWrappedArray, or another approach for wrapped array identification · Issue #51910 · JuliaLang/julia · GitHub](https://github.com/JuliaLang/julia/issues/51910)). We try to support Base’s array wrappers as much as possible, and for other types like OffsetArray a package extension that fixes or overrides dispatch where needed could be added.

If you simply want compatibility (i.e., without triggering scalar indexing errors, but also without executing on the GPU) you can use unified memory, see [CUDA.jl 5.4: Memory management mayhem ⋅ JuliaGPU](https://juliagpu.org/post/2024-05-28-cuda_5.4/#unified_memory_iteration)

---

<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:** [August 20, 2024, 3:01pm UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/5 "2024-08-20T15:01:55Z")

</div>

Thank you both. So, as it stands, we need to define this roundabout wrapper for ourselves.

- unified memory seems not to support float64.
- adapt seems not to support broadcast or reducemap.

---

<div class="post-metadata">

**Author:** ![maleadt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maleadt/32/10097_2.png) [@maleadt](https://discourse.julialang.org/u/maleadt)\
**Post date:** [August 20, 2024, 3:36pm UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/6 "2024-08-20T15:36:48Z")

</div>

> [@0samuraiE](#):
>
> unified memory seems not to support float64.

That is not the case. Can you share what you are running into?

---

<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:** [August 20, 2024, 11:18pm UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/7 "2024-08-20T23:18:14Z")

</div>

I tried

```julia
A = CUDA.zeros(10)
A = cu(A, unified=true)

```

Maybe CuArray{Float64,1,CUDA.UnifiedMemory} works? but I cannot try yet.

Additionally it seems that unified memory is allocated on cpu memory. Is this cause any performance issues?

---

<div class="post-metadata">

**Author:** ![maleadt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maleadt/32/10097_2.png) [@maleadt](https://discourse.julialang.org/u/maleadt)\
**Post date:** [August 21, 2024, 7:17am UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/8 "2024-08-21T07:17:48Z")

</div>

The `cu` function is to be used with CPU inputs; It’s a user-friendly constructor.

---

<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:** [August 21, 2024, 11:58pm UTC](https://discourse.julialang.org/t/how-to-use-offsetarray-with-cuda/118400/9 "2024-08-21T23:58:41Z")

</div>

Thank you. I understand.