# Solving ODE with inheritance of initial condition with parallelization by threads

**URL:** <https://discourse.julialang.org/t/solving-ode-with-inheritance-of-initial-condition-with-parallelization-by-threads/105418>\
**Category:** General Usage\
**Created:** [October 26, 2023, 5:55am UTC](https://discourse.julialang.org/t/solving-ode-with-inheritance-of-initial-condition-with-parallelization-by-threads/105418 "2023-10-26T05:55:22Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![Sergey\_Novak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergey_novak/32/37716_2.png) [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)\
**Post date:** [October 26, 2023, 5:55am UTC](https://discourse.julialang.org/t/solving-ode-with-inheritance-of-initial-condition-with-parallelization-by-threads/105418/1 "2023-10-26T05:55:23Z")

</div>

Hello!  
I have an ODE and I need to perform inheritance of initial conditions when changing two parameters. I can do this without parallelizing threads, and everything works correctly. But when I start parallelizing, I get the wrong result. The variable `u0_start_2d` takes incorrect values. I assume that one thread interferes with the work of another thread. But I don’t understand how this happens.

A fragment of code that I have paralleled:

```julia
Threads.@threads for p2_index_range in eachindex(p2_range)
        for p1_index_range in eachindex(p1_range)

            if p1_index_range == 1
                u0_start_2d = SVector{6}(matrix_u0_end[p1_index_range, p2_index_range, :]);
                continue;
            end
            
            loc_params = copy(params);
            loc_params[index_p1] = p1_range[p1_index_range];
            loc_params[index_p2] = p2_range[p2_index_range];

            prob = ODEProblem(sys, u0_start_2d, tspan, loc_params);
            sol = solver(prob, integ_set);

            if sol.retcode == ReturnCode.MaxIters
                exitcode = 1;
            end

            u0_end_2d = sol[end];

            matrix_u0_start[p1_index_range, p2_index_range, :] = u0_start_2d;
            matrix_u0_end[p1_index_range, p2_index_range, :] = u0_end_2d;
            matrix_params[p1_index_range, p2_index_range, :] = [loc_params[index_p1], loc_params[index_p2]]
            
            u0_start_2d = u0_end_2d;

            if exitcode == 1
                exit();
            end;
        end
    end

```

Link on full project: [source code - Google Drive](https://drive.google.com/drive/folders/17XH5n2rACxlsXiZoCKnUz6pmJobILLEr?usp=sharing)

---

<div class="post-metadata">

**Author:** ![abraemer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/abraemer/32/51403_2.png) [@abraemer](https://discourse.julialang.org/u/abraemer)\
**Post date:** [October 26, 2023, 6:08am UTC](https://discourse.julialang.org/t/solving-ode-with-inheritance-of-initial-condition-with-parallelization-by-threads/105418/2 "2023-10-26T06:08:59Z")

</div>

So something from before the loop might be interesting but I guess you do something like:

```julia
u0_start_2d = ....
Threads.@threads for p2_index_range in eachindex(p2_range)
        for p1_index_range in eachindex(p1_range)

            if p1_index_range == 1
                u0_start_2d = SVector{6}(matrix_u0_end[p1_index_range, p2_index_range, :]);
                continue;
            end
            ...
            u0_start_2d = u0_end_2d;

```

Which means each thread sees the same variable `u0_start_2d` and overwrite each others values!  
A fix can simply be moving the variable inside the outer loop:

```julia
### u0_start_2d = .... # from here
Threads.@threads for p2_index_range in eachindex(p2_range)
        u0_start_2d = .... # to here
        for p1_index_range in eachindex(p1_range)

```

then the variable will be local in each thread and is safe to use and overwrite.

You could also refactor this a bit to:

```julia
Threads.@threads for p2_index_range in eachindex(p2_range)
        u0_start_2d = SVector{6}(matrix_u0_end[1, p2_index_range, :])
        for p1_index_range in eachindex(p1_range)
           p1_index_range == 1 && continue

           loc_params = copy(params)
           # rest of code...

```

Also note you don’t need “;” at the end of each line - we are not in C/Java land anymore and thus free of the oppression of semicolons 😄

---

<div class="post-metadata">

**Author:** ![Sergey\_Novak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergey_novak/32/37716_2.png) [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)\
**Post date:** [October 26, 2023, 8:17am UTC](https://discourse.julialang.org/t/solving-ode-with-inheritance-of-initial-condition-with-parallelization-by-threads/105418/3 "2023-10-26T08:17:58Z")

</div>

Thank you very much  
Everything is working now

---

<div class="post-metadata">

**Author:** ![DrPapa](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/drpapa/32/6835_2.png) [@DrPapa](https://discourse.julialang.org/u/DrPapa)\
**Post date:** [October 26, 2023, 1:29pm UTC](https://discourse.julialang.org/t/solving-ode-with-inheritance-of-initial-condition-with-parallelization-by-threads/105418/4 "2023-10-26T13:29:05Z")

</div>

@Sergey_Novak If you are not aware of if, you may be interested in the `EnsembleProblem` interface. [https://docs.sciml.ai/DiffEqDocs/dev/features/ensemble/](https://docs.sciml.ai/DiffEqDocs/dev/features/ensemble/)