# Solution indexing is very slow

**URL:** <https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896>\
**Category:** Modelling & Simulations\
**Tags:** modelingtoolkit\
**Created:** [August 6, 2024, 2:18pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896 "2024-08-06T14:18:51Z")\
**Posts on this page:** 9\
**Page:** 1

<div class="post-metadata">

**Author:** ![Bart\_van\_de\_Lint](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bart_van_de_lint/32/212161_2.png) [@Bart\_van\_de\_Lint](https://discourse.julialang.org/u/Bart_van_de_Lint)\
**Post date:** [August 6, 2024, 2:18pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/1 "2024-08-06T14:18:51Z")

</div>

When running the following code, the last two lines do the same thing: they get decay2.f at the last point in time. sol[decay2.f][end] is the modelingtoolkit preferred way, while sol.u[end][2] is using normal array indexing.

```julia
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D

function decay(; name)
    @parameters a
    @variables x(t) f(t)
    ODESystem([
            D(x) ~ -a * x + f
        ], t;
        name = name)
end

@named decay1 = decay()
@named decay2 = decay()

connected = compose(
    ODESystem([decay2.f ~ decay1.x
               D(decay1.f) ~ 0], t; name = :connected), decay1, decay2)

equations(connected)
simplified_sys = structural_simplify(connected)
equations(simplified_sys)

x0 = [decay1.x => 1.0
      decay1.f => 0.0
      decay2.x => 1.0]
p = [decay1.a => 0.1
     decay2.a => 0.2]

using OrdinaryDiffEq
prob = ODEProblem(simplified_sys, x0, (0.0, 100.0), p)
sol = solve(prob, Tsit5())
@time sol[decay2.f][end]
@time sol.u[end][2]

```

When timing these two lines of code, I get the following results:  
0.000734 seconds (1.59 k allocations: 114.008 KiB)  
0.000014 seconds (2 allocations: 1024 bytes)

Somehow, sol[decay2.f][end] takes 50 times longer doing the same thing. Are there any possible alternatives or solutions to this? I am developing performance critical simulation code, where this time difference is hurting performance. And using normal array indexing is not safe, as an update of ModelingToolkit could change the order of the solution variables.

---

<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:** [August 6, 2024, 2:53pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/2 "2024-08-06T14:53:27Z")

</div>

You want to cache the symbol lookup function. This is what `getu` is for. [Optimizing through an ODE solve and re-creating MTK Problems · ModelingToolkit.jl](https://docs.sciml.ai/ModelingToolkit/stable/examples/remake/) is somewhat helpful. Try doing it, and post what you get to, and complain, and @cryptic.ax we should write a proper tutorial on this.

---

<div class="post-metadata">

**Author:** ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)\
**Post date:** [August 6, 2024, 4:28pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/3 "2024-08-06T16:28:37Z")

</div>

The following code works for me:

```julia
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D
using SymbolicIndexingInterface

function decay(; name)
    @parameters a
    @variables x(t) f(t)
    ODESystem([
            D(x) ~ -a * x + f
        ], t;
        name = name)
end

@named decay1 = decay()
@named decay2 = decay()

connected = compose(
    ODESystem([decay2.f ~ decay1.x
               D(decay1.f) ~ 0], t; name = :connected), decay1, decay2)

equations(connected)
simplified_sys = structural_simplify(connected)
equations(simplified_sys)

x0 = [decay1.x => 1.0
      decay1.f => 0.0
      decay2.x => 1.0]
p = [decay1.a => 0.1
     decay2.a => 0.2]

using OrdinaryDiffEq
prob = ODEProblem(simplified_sys, x0, (0.0, 100.0), p)
sol = solve(prob, Tsit5())
get_decay2_f = getu(sol, decay2.f)
@time sol[decay2.f][end]
@time sol.u[end][2]
@time get_decay2_f(sol)[end]

```

When you get the first solution you can create a function like `get_decay2_f(sol)` for all variables you are interested in, and later you just use these functions which is fast.

---

<div class="post-metadata">

**Author:** ![Bart\_van\_de\_Lint](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bart_van_de_lint/32/212161_2.png) [@Bart\_van\_de\_Lint](https://discourse.julialang.org/u/Bart_van_de_Lint)\
**Post date:** [August 6, 2024, 4:50pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/4 "2024-08-06T16:50:05Z")

</div>

Thanks, that works!  
I am getting the following results with Fechners code:

```julia
@time sol[decay2.f][end]
@time sol.u[end][2]
@time get_decay2_f(sol)[end]
  0.000141 seconds (276 allocations: 15.031 KiB)
  0.000012 seconds (2 allocations: 1024 bytes)
  0.000028 seconds (147 allocations: 3.453 KiB)

```

---

<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:** [August 6, 2024, 5:36pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/5 "2024-08-06T17:36:47Z")

</div>

`get_decay2_f` should be non-allocating. Is that timing its compilation?

---

<div class="post-metadata">

**Author:** ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)\
**Post date:** [August 6, 2024, 7:05pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/6 "2024-08-06T19:05:13Z")

</div>

> [@ufechner7](#):
>
> ```julia
> get_decay2_f = getu(sol, decay2.f)
> @time sol[decay2.f][end]
> @time sol.u[end][2]
> @time get_decay2_f(sol)[end]
> 
> ```

If I run my code multiple times I get:

```julia
julia> include("mwes/mwe_13b.jl")
  0.000173 seconds (131 allocations: 11.828 KiB)
  0.000017 seconds (2 allocations: 1024 bytes)
  0.000014 seconds (2 allocations: 256 bytes)
4.542765615967835e-5

```

Two allocations, but neither zero and nor 147.

---

<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:** [August 6, 2024, 8:17pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/7 "2024-08-06T20:17:55Z")

</div>

That’s because `get_decay2_f` isn’t const. If you created that function as `const` or used it in a function that global lookup should go away?

---

<div class="post-metadata">

**Author:** ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)\
**Post date:** [August 6, 2024, 8:41pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/8 "2024-08-06T20:41:10Z")

</div>

Well, it doesn’t go away.

```julia
get_decay2_f = getu(sol, decay2.f)

function process_result(sol, get_decay2_f)
    @time sol[decay2.f][end]
    @time sol.u[end][2]
    @time get_decay2_f(sol)[end]
end

process_result(sol)

```

Output:

```julia
julia> process_result(sol)
  0.000019 seconds (132 allocations: 12.812 KiB)
  0.000001 seconds
  0.000001 seconds (3 allocations: 1.234 KiB)
4.542765615967835e-5

```

But it doesn’t matter. Good enough.

---

<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:** [August 6, 2024, 9:23pm UTC](https://discourse.julialang.org/t/solution-indexing-is-very-slow/117896/9 "2024-08-06T21:23:05Z")

</div>

We can open an issue and get that fixed. It’s pretty close though so at least most people will be happy, but embedded will need it fixed.