# Scalar Indexing Error in PINN Implementation on GPU

**URL:** <https://discourse.julialang.org/t/scalar-indexing-error-in-pinn-implementation-on-gpu/116715>\
**Category:** General Usage\
**Tags:** question, machine-learning, modelingtoolkit\
**Created:** [July 7, 2024, 8:08am UTC](https://discourse.julialang.org/t/scalar-indexing-error-in-pinn-implementation-on-gpu/116715 "2024-07-07T08:08:15Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![Maysam\_Gholampour](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maysam_gholampour/32/218184_2.png) [@Maysam\_Gholampour](https://discourse.julialang.org/u/Maysam_Gholampour)\
**Post date:** [July 7, 2024, 8:08am UTC](https://discourse.julialang.org/t/scalar-indexing-error-in-pinn-implementation-on-gpu/116715/1 "2024-07-07T08:08:15Z")

</div>

Dear All,

I am attempting to run the following Physics-Informed Neural Network (PINN) problem on a GPU. However, I encounter an error related to scalar indexing. Despite my efforts to debug this issue and explore potential solutions, I have been unable to find any relevant threads or guidance.  
Could you please provide some advice or suggestions on how to address and resolve this scalar indexing error?

Thank you for your assistance.

code ===\>

```julia
begin
    using ModelingToolkit, DomainSets
    using Lux
    using NeuralPDE
    using LuxCUDA, ComponentArrays, Random
end

# Parameters, variables, and derivatives
@parameters Xꜛ Yꜛₛₒₗ Yꜛₐ
@variables Θꜛₐ(..) ωꜛₐ(..) 
@variables Θꜛₛₒₗ(..) ξꜛₛₒₗ(..) 

DXꜛ = Differential(Xꜛ)
DYꜛₐ = Differential(Yꜛₐ)
DYꜛₛₒₗ = Differential(Yꜛₛₒₗ)
DYꜛₐ² = Differential(Yꜛₐ)^2
DYꜛₛₒₗ² = Differential(Yꜛₛₒₗ)^2

function Uꜛₛₒₗ(Yꜛₛₒₗ)
    y = Yꜛₛₒₗ * 1.2
    uₛₒₗ = 2 * y * (1.0 - 0.5 * y) / 4.0
    return uₛₒₗ / 1.5
end

function Uꜛₐ(Yꜛₐ)
    y = Yꜛₐ * 1.2
    uₐ = -1.2 - 5.0 *(3.0 ^ 2 - y ^ 2)
    return -uₐ / 2.0
end

i_fg(t_) = 4t_^3 - 3t_^4 + 10.0

function c2(Θ)
    t_ = Θ * 5.0 + 10.0
    return 10.0 * i_fg(t_) * 0.5 / (10.0 * 5000.0)
end
@register_symbolic c2(Θ)

# Equations
eqs = [
    Uꜛₐ(Yꜛₐ) * DXꜛ(Θꜛₐ(Xꜛ,Yꜛₐ)) ~ 5.0 * DYꜛₐ²(Θꜛₐ(Xꜛ,Yꜛₐ)),
    Uꜛₐ(Yꜛₐ) * DXꜛ(ωꜛₐ(Xꜛ,Yꜛₐ)) ~ 6.0 * DYꜛₐ²(ωꜛₐ(Xꜛ,Yꜛₐ)),
    Uꜛₛₒₗ(Yꜛₛₒₗ) * DXꜛ(Θꜛₛₒₗ(Xꜛ,Yꜛₛₒₗ)) ~ 8.0 * DYꜛₛₒₗ²(Θꜛₛₒₗ(Xꜛ,Yꜛₛₒₗ)),
    Uꜛₛₒₗ(Yꜛₛₒₗ) * DXꜛ(ξꜛₛₒₗ(Xꜛ,Yꜛₛₒₗ)) ~ 11.0 * DYꜛₛₒₗ²(ξꜛₛₒₗ(Xꜛ,Yꜛₛₒₗ))
]

bcs = [
    Θꜛₐ(1.0, Yꜛₐ) ~ 1.0,
    ωꜛₐ(1.0, Yꜛₐ) ~ 1.0,
    Θꜛₛₒₗ(0.0, Yꜛₛₒₗ) ~ 1.0,
    ξꜛₛₒₗ(0.0, Yꜛₛₒₗ) ~ 1.0,

    DYꜛₐ(Θꜛₐ(Xꜛ, 0.0)) ~ 0.0,
    DYꜛₐ(ωꜛₐ(Xꜛ, 0.0)) ~ 0.0,
    Θꜛₛₒₗ(Xꜛ, 0.0) ~ 0.0,
    DYꜛₛₒₗ(ξꜛₛₒₗ(Xꜛ, 0.0)) ~ 0.0,

    DYꜛₐ(Θꜛₐ(Xꜛ, 1.0)) + 5.0 * DYꜛₛₒₗ(Θꜛₛₒₗ(Xꜛ, 1.0)) - 5.0 * c2(Θꜛₛₒₗ(Xꜛ, 1.0)) * DYꜛₛₒₗ(ξꜛₛₒₗ(Xꜛ, 1.0)) ~ 0.0,
    DYꜛₐ(ωꜛₐ(Xꜛ, 1.0)) + 1.05 * DYꜛₛₒₗ(ξꜛₛₒₗ(Xꜛ, 1.0)) ~ 0.0
]

domains = [
    Xꜛ ∈ Interval(0.0, 1.0),
    Yꜛₐ ∈ Interval(0.0, 1.0),
    Yꜛₛₒₗ ∈ Interval(0.0, 1.0),
]

@named pdesystem = PDESystem(eqs, bcs, domains,[Xꜛ, Yꜛₐ, Yꜛₛₒₗ], [Θꜛₐ(Xꜛ,Yꜛₐ), ωꜛₐ(Xꜛ,Yꜛₐ), Θꜛₛₒₗ(Xꜛ,Yꜛₛₒₗ), ξꜛₛₒₗ(Xꜛ,Yꜛₛₒₗ)])

# Neural network
dim = 2 # number of dimensions
N_nr = 8 # number of neurons
chain = [Lux.Chain(Dense(dim, N_nr, Lux.σ), Dense(N_nr, N_nr, Lux.σ), Dense(N_nr, 1)) for _ in 1:4]

strategy = QuasiRandomTraining(100)
gpu_dev = gpu_device()
# Seeding
rng = Random.default_rng()
Random.seed!(rng, 0)

ps = Lux.setup(Random.default_rng(), chain)[1] 
ps = ps |> gpu_dev
ps = f64(ps)

discretization = PhysicsInformedNN(chain, strategy,init_params = ps)
prob = discretize(pdesystem, discretization)

```

error ===\>

```julia
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. 

```

---

<div class="post-metadata">

**Author:** ![Maysam\_Gholampour](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maysam_gholampour/32/218184_2.png) [@Maysam\_Gholampour](https://discourse.julialang.org/u/Maysam_Gholampour)\
**Post date:** [July 8, 2024, 10:51am UTC](https://discourse.julialang.org/t/scalar-indexing-error-in-pinn-implementation-on-gpu/116715/2 "2024-07-08T10:51:36Z")

</div>

I think this issue arises from c2(Θꜛₛₒₗ(Xꜛ, 1.0)) in Bcs but I don’t know how to handle it. It might stem form scalar indexing in discretization.

thanks in advance for your assistance and taking the time to help.

@ChrisRackauckas

error ===\>

```julia
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 C:\Users\lab30501\.julia\packages\GPUArraysCore\GMsgk\src\GPUArraysCore.jl:155        
  [3] _assertscalar(op::String, behavior::GPUArraysCore.ScalarIndexing)
    @ GPUArraysCore C:\Users\lab30501\.julia\packages\GPUArraysCore\GMsgk\src\GPUArraysCore.jl:128        
  [4] assertscalar(op::String)
    @ GPUArraysCore C:\Users\lab30501\.julia\packages\GPUArraysCore\GMsgk\src\GPUArraysCore.jl:116        
  [5] getindex(A::CuArray{Float64, 2, CUDA.DeviceMemory}, I::Int64)
    @ GPUArrays C:\Users\lab30501\.julia\packages\GPUArrays\8Y80U\src\host\indexing.jl:48
  [6] iterate
    @ .\abstractarray.jl:1217 [inlined]
  [7] iterate
    @ .\abstractarray.jl:1215 [inlined]
  [8] _zip_iterate_some
    @ .\iterators.jl:428 [inlined]
  [9] _zip_iterate_some
    @ .\iterators.jl:430 [inlined]
 [10] _zip_iterate_all
    @ .\iterators.jl:420 [inlined]
 [11] iterate
    @ .\iterators.jl:410 [inlined]
 [12] _append!(a::Vector{Float64}, ::Base.HasShape{2}, iter::CuArray{Float64, 2, CUDA.DeviceMemory})      
    @ Base .\array.jl:1197
 [13] append!
    @ .\array.jl:1187 [inlined]
 [14] make_idx(data::Vector{Float64}, x::CuArray{Float64, 2, CUDA.DeviceMemory}, last_val::Int64)
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:186   
 [15] #28
    @ C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:161 [inlined]
 [16] (::ComponentArrays.var"#28#29"{Vector{…}, Base.RefValue{…}})(::Pair{Symbol, CuArray{…}})
    @ ComponentArrays .\none:0
 [17] iterate
    @ .\generator.jl:47 [inlined]
 [18] merge(a::@NamedTuple{}, itr::Base.Generator{@Kwargs{…}, ComponentArrays.var"#28#29"{…}})
    @ Base .\namedtuple.jl:361
 [19] make_idx(data::Vector{Float64}, nt::@NamedTuple{weight::CuArray{…}, bias::CuArray{…}}, last_val::Int64)
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:159   
 [20] #28
    @ C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:161 [inlined]
 [21] (::ComponentArrays.var"#28#29"{Vector{…}, Base.RefValue{…}})(::Pair{Symbol, @NamedTuple{…}})
    @ ComponentArrays .\none:0
 [22] iterate
    @ .\generator.jl:47 [inlined]
 [23] merge(a::@NamedTuple{}, itr::Base.Generator{@Kwargs{…}, ComponentArrays.var"#28#29"{…}})
    @ Base .\namedtuple.jl:361
 [24] make_idx(data::Vector{…}, nt::@NamedTuple{…}, last_val::Int64)
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:159   
 [25] #28
    @ C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:161 [inlined]
 [26] (::ComponentArrays.var"#28#29"{Vector{…}, Base.RefValue{…}})(::Pair{Symbol, @NamedTuple{…}})        
    @ ComponentArrays .\none:0
 [27] iterate
    @ .\generator.jl:47 [inlined]
 [28] merge(a::@NamedTuple{}, itr::Base.Generator{@Kwargs{…}, ComponentArrays.var"#28#29"{…}})
    @ Base .\namedtuple.jl:361
 [29] make_idx(data::Vector{…}, nt::@NamedTuple{…}, last_val::Int64)
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:159   
 [30] make_carray_args(A::Type{…}, nt::@NamedTuple{…})
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:151   
 [31] make_carray_args(nt::@NamedTuple{…})
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:144   
 [32] ComponentArray(nt::@NamedTuple{Θꜛₐ::@NamedTuple{…}, ωꜛₐ::@NamedTuple{…}, Θꜛₛₒₗ::@NamedTuple{…}, ξꜛₛₒₗ::@NamedTuple{…}})
    @ ComponentArrays C:\Users\lab30501\.julia\packages\ComponentArrays\joxVV\src\componentarray.jl:64    
 [33] symbolic_discretize(pde_system::PDESystem, discretization::PhysicsInformedNN{…})
    @ NeuralPDE C:\Users\lab30501\.julia\packages\NeuralPDE\z18Qg\src\discretize.jl:415
 [34] discretize(pde_system::PDESystem, discretization::PhysicsInformedNN{…})
    @ NeuralPDE C:\Users\lab30501\.julia\packages\NeuralPDE\z18Qg\src\discretize.jl:708
 [35] top-level scope
    @ d:\Maysam\01 - Projects\05_Ionic_Liquid_Dsign_Simulation\_Julia_code\_NPDE_2d\MWE.jl:88
Some type information was truncated. Use `show(err)` to see complete types.

```

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [July 9, 2024, 5:10am UTC](https://discourse.julialang.org/t/scalar-indexing-error-in-pinn-implementation-on-gpu/116715/3 "2024-07-09T05:10:15Z")

</div>

@SathvikBhagavan is the one that has been working on related issues. Is this known?

---

<div class="post-metadata">

**Author:** ![Maysam\_Gholampour](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maysam_gholampour/32/218184_2.png) [@Maysam\_Gholampour](https://discourse.julialang.org/u/Maysam_Gholampour)\
**Post date:** [July 10, 2024, 12:50am UTC](https://discourse.julialang.org/t/scalar-indexing-error-in-pinn-implementation-on-gpu/116715/4 "2024-07-10T00:50:32Z")

</div>

Thanks for your response.

I changed

> [@Maysam\_Gholampour](#):
>
> `ps = ps |> gpu_dev`

to  
`ps = ps .|> ComponentArray |> gpu_dev`.

It works.