# Double for loop using filtered ranges

**URL:** <https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722>\
**Category:** Performance\
**Created:** [March 11, 2019, 3:06am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722 "2019-03-11T03:06:04Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![MrRobot](https://avatars.discourse-cdn.com/v4/letter/m/b38774/32.png) [@MrRobot](https://discourse.julialang.org/u/MrRobot)\
**Post date:** [March 11, 2019, 3:06am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722/1 "2019-03-11T03:06:04Z")

</div>

I have some MWE below. What I want is to have a subsection of a range, interact with the rest of the range, but not itself.

Also, I posted this on [stackoverflow](https://stackoverflow.com/questions/55094612/filtering-a-range-in-julia) … I should have put it here first I think…

For instance if the range is `1:100` , I want to have a for loop that will have each index in `4:6` , interact with all values of `1:100` BUT NOT `4:6` .

I want to do this using ranges/filters to avoid generating temporary arrays.

In my case the total range is the number of atoms in the system. The sub-range, is the atoms in a specific molecule. I need to do calculations where each atom in a molecule interacts with all other atoms, but not the atoms in the same molecule.

## Further

I am trying to avoid using if statements because that can potentially complicate parallelizing the code later (such is my impression, correct me if I am wrong). Doing this with an if statement would be

```julia-auto
for i=4:6
    for j = 1:100
        if j == 4 || j==5 || j==6
            continue
        end
        println(i, " ", j) #roughly 0.1-100 microseconds to run the real call here
    end
end

```

I have actual indexing in my code, I would never hardcode values like the above… But I want to avoid that `if` statement.

## Trials

The following does what I want, but I now realize that using `filter` is bad when it comes to memory and the amount used scales linearly with `b` .

```julia-auto
a = 4:6
b = 1:100
for i in a
    for j in filter((b) -> !(b in a),b) 
        print(i, " ", j) 
    end
end

```

Is there a way to get the double for loop I want where the outer is a sub-range of the inner, but the inner does not include the outer sub-range and most importantly is fast and does not create alot of memory usage like filter does?

---

<div class="post-metadata">

**Author:** ![MrRobot](https://avatars.discourse-cdn.com/v4/letter/m/b38774/32.png) [@MrRobot](https://discourse.julialang.org/u/MrRobot)\
**Post date:** [March 11, 2019, 3:24am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722/2 "2019-03-11T03:24:04Z")

</div>

A user on stackoverflow recommended the following package

> **[GitHub - JuliaData/InvertedIndices.jl: A simple index type that allows for...](https://github.com/JuliaData/InvertedIndices.jl)**
>
> A simple index type that allows for inverted selections - GitHub - JuliaData/InvertedIndices.jl: A simple index type that allows for inverted selections

It looks like it does what I am after, but it doesn’t seem to be registered :S

---

<div class="post-metadata">

**Author:** ![bennedich](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bennedich/32/4894_2.png) [@bennedich](https://discourse.julialang.org/u/bennedich)\
**Post date:** [March 11, 2019, 4:56am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722/3 "2019-03-11T04:56:34Z")

</div>

You can use a generator to do the filtering. Below I also changed your inner function (`print` is not so nice for benchmarking) and made it return a summed value to be able to verify correctness:

```julia
func(i, j) = sin(i) * cos(j)

function loop(range, n)
    s = 0.0
    first, last = range[1], range[end]
    for i = range
        for j = (k for k = 1:n if k < first || k > last)
            s += func(i, j)
        end
    end
    s
end

```

Usage and timing:

```julia
julia> @btime loop(4:6, 100)
  5.140 μs (0 allocations: 0 bytes)
2.239502224477359

```

Could you explain a bit more how you want to parallelize this, and why using `if` statements or similar is not an option? Do you intend to parallelize the inner loop directly? What will the ranges be for your real code, and for how many molecules will you run this?

---

<div class="post-metadata">

**Author:** ![MrRobot](https://avatars.discourse-cdn.com/v4/letter/m/b38774/32.png) [@MrRobot](https://discourse.julialang.org/u/MrRobot)\
**Post date:** [March 11, 2019, 5:21am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722/4 "2019-03-11T05:21:37Z")

</div>

Great stuff… I should get your business card one of these days.

> [@bennedich](#):
>
> Could you explain a bit more how you want to parallelize this, and why using `if` statements or similar is not an option? Do you intend to parallelize the inner loop directly? What will the ranges be for your real code, and for how many molecules will you run this?

I must start by saying I do not know much about parallel coding. I have used openMP before. My reason for not wanting if statements is that I perhaps naively think the early returns will interfere with the management of the threads. Edit: in this case, each pass has the same amount of if statements, so running things in parallel should not affect timing at all.

I would like to make my code GPU usable, rather than cpu parallel if I can. I don’t know what dragons come along with that choice. However, I think the easiest form of parallel coding is cpu openMP, which I can’t seem to find for Julia. I have prime algorithms for openMP but it must be able to handle some sort of reduction. `@threads` does not seem to do this, so it gives me wrong summations. Using `mpi` is an alright option, but I think gpu would be preferable.

Molecules… well, thousands, which can lead to upwards of 30,000 atoms. I don’t plan on making my code enter the millions of atoms range, but I think it should be able to do 30,000. To be honest though, I will get a great deal done with it using only 5000 atoms.

So long as it is serial I have to stick to around 3000 atoms or less. My monte carlo code is inevitable N^2, but the molecular dynamics can be made N \log N

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [March 11, 2019, 8:33am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722/5 "2019-03-11T08:33:28Z")

</div>

> [@bennedich](#):
>
> `first, last = range[1], range[end]`

Shadowing both `range`, `first` and `last`?🤨 That seems, at least, like a dangerous habit 😉

It could be more generic and shorter by using `!(k in myrange)`.

---

<div class="post-metadata">

**Author:** ![bennedich](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bennedich/32/4894_2.png) [@bennedich](https://discourse.julialang.org/u/bennedich)\
**Post date:** [March 11, 2019, 9:31am UTC](https://discourse.julialang.org/t/double-for-loop-using-filtered-ranges/21722/6 "2019-03-11T09:31:21Z")

</div>

> [@MrRobot](#):
>
> I have prime algorithms for openMP but it must be able to handle some sort of reduction. `@threads` does not seem to do this, so it gives me wrong summations.

What do you mean by `@threads` not supporting reductions?

Do you intend to parallelize this on a single computer only (using all cores), or a cluster?

When I asked about the number of molecules I meant if you’ll be calling this `loop` function multiple times, or just a single time? Since if you’ll be calling it multiple times, it might be wiser to parallelize those calls, rather than parallelizing the inner-most loop.

> [@DNF](#):
>
> Shadowing both `range` , `first` and `last` ?😉 That seems, at least, like a dangerous habit 😉

Haha, yes good point. The only thing missing was calling it `sum` instead of `s` 🙂

> [@DNF](#):
>
> It could be more generic and shorter by using `!(k in myrange)` .

Again good point, I was initially thinking an implementation like this:

```julia
function loop2(molecule_range, n)
    s = 0.0
    first_idx, last_idx = molecule_range[1], molecule_range[end]
    for i = molecule_range
        for j = 1:first_idx-1
            s += func(i, j)
        end
        for j = last_idx+1:n
            s += func(i, j)
        end
    end
    s
end

```

(but probably cleaned up to remove the code duplication), which would be a lot faster for large ranges:

```julia
julia> @btime loop(2:99, 100)
  12.864 μs (0 allocations: 0 bytes)
-0.6483986089918431

julia> @btime loop2(2:99, 100)
  2.191 μs (0 allocations: 0 bytes)
-0.6483986089918431

```

But for small ranges it’s likely not worth it, in which case I agree that your suggestion is cleaner.