# Elegant way to construct a set of multi-index

**URL:** <https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081>\
**Category:** Numerics\
**Created:** [November 27, 2018, 11:18pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081 "2018-11-27T23:18:18Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 27, 2018, 11:18pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/1 "2018-11-27T23:18:18Z")

</div>

Good evening,

I wish to build an array containing specific multi-index. More precisely, in the 2-D case the multi-index are given by `(i, j)` where ` 0 <= i , j <= N+1` such that `i+j <= N+1`. The way I do it is make a (N+1)x(N+1) array and cut it in two along its anti-diagonal, namely

```julia
βtab = CartesianIndices((N+1, N+1))
βtab = Tuple.( vec(rotl90( UpperTriangular( rotr90(βtab) ) )) )
filter!(x -> sum(x) ≠ 0, βtab)

```

Does Julia have features that could be used to implement this in a more elegant way ?

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [November 27, 2018, 11:49pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/3 "2018-11-27T23:49:34Z")

</div>

I’d just use logical indexing and construct the logical array just like you explained it:

```julia
julia> N = 5; A = reshape(1:N^2, N, N);

julia> B = [i+j <= N+1 for i in 1:N, j in 1:N]
5×5 Array{Bool,2}:
 true true true true true
 true true true true false
 true true true false false
 true true false false false
 true false false false false

julia> A[B]
15-element Array{Int64,1}:
  1
  2
  3
  4
  5
  6
  7
  8
  9
 11
 12
 13
 16
 17
 21

```

---

<div class="post-metadata">

**Author:** ![jandehaan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jandehaan/32/6805_2.png) [@jandehaan](https://discourse.julialang.org/u/jandehaan)\
**Post date:** [November 28, 2018, 12:04am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/4 "2018-11-28T00:04:40Z")

</div>

Second attempt.

```julia
julia> N = 5
julia> a=[]
julia> for j in 1:N+1
           a=[a;[(i,j) for i in j:N+1]]
       end
julia> a
1-element Array{Tuple{Int64,Int64},1}:
 (1, 1)
 (2, 1)
 (3, 1)
 (4, 1)
 (5, 1)
 (6, 1)
 (1, 2)
 (2, 2)
 (3, 2)
 (4, 2)
 (5, 2)
 (1, 3)
 (2, 3)
 (3, 3)
 (4, 3)
 (1, 4)
 (2, 4)
 (3, 4)
 (1, 5)
 (2, 5)
 (1, 6)

```

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 28, 2018, 12:24am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/5 "2018-11-28T00:24:13Z")

</div>

Thank you for the suggestions! I should have been more precise, my goal is to implement it so that Julia can automatically treat the N-D case.  
From both suggestions, it would mean that Julia can automatically increase the number of loop (since they rely on `[something for i in I]`).  
I cannot seem to exploit the geometry of the problem. The code should link M+1 nodes with straight lines and fill the interior. That is, in 2D we link (0,0), (N,0) and (0,N) [triangle]. In 3D we link (0,0,0), (N,0,0), (0,N,0) and (0,0,N) [tetrahedron].

---

<div class="post-metadata">

**Author:** ![jandehaan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jandehaan/32/6805_2.png) [@jandehaan](https://discourse.julialang.org/u/jandehaan)\
**Post date:** [November 28, 2018, 2:41am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/6 "2018-11-28T02:41:38Z")

</div>

Attempt 3 with caveats

- it is not elegant but it gets the job done.
- Your explanation of the geometry of the problem got me confused. I’m not sure if my solution addresses the correct problem.

```julia
macro makenodes(N,D)
    quote
        a = []
        @makenodes($N,$D,1)
        a
    end
end

```

```julia
macro makenodes(N, D, i)
    varname = Symbol(:d, i)
    if i > D
        # generate the statement inside the inner loop
        vars = ()
        for j in 1:D
            vars = (vars..., Symbol(:d,j))
        end
        expr1 = Expr(:call, :+, vars...)
        expr2 = Expr(:call, :push!, :a, Expr(:tuple, vars...))
        return Expr(:if, Expr(:call, :(<=), expr1, N+1), expr2)
    else
        # generate as many loops as there are dimensions
        return quote
            for $varname in 0:$(N+1)
                @makenodes($N, $D, $(i+1))
            end
        end
    end
 end

```

```julia
julia> @makenodes(5,3)
84-element Array{Any,1}:
 (0, 0, 0)
 (0, 0, 1)
 (0, 0, 2)
 (0, 0, 3)
 (0, 0, 4)
 (0, 0, 5)
 (0, 0, 6)
 (0, 1, 0)
 (0, 1, 1)
 (0, 1, 2)
 (0, 1, 3)
 (0, 1, 4)
 (0, 1, 5)
 ⋮        
 (3, 2, 1)
 (3, 3, 0)
 (4, 0, 0)
 (4, 0, 1)
 (4, 0, 2)
 (4, 1, 0)
 (4, 1, 1)
 (4, 2, 0)
 (5, 0, 0)
 (5, 0, 1)
 (5, 1, 0)
 (6, 0, 0)

```

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [November 28, 2018, 3:26am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/7 "2018-11-28T03:26:12Z")

</div>

> [@OlivierHnt](#):
>
> my goal is to implement it so that Julia can automatically treat the N-D case.  
> From both suggestions, it would mean that Julia can automatically increase the number of loop

You could do something with the [`@nloops` macro in `Base.Cartesian`](https://docs.julialang.org/en/v1/devdocs/cartesian/#Base.Cartesian.@nloops), for example:

```julia
@generated function gentups(::Val{d}, N::Integer) where {d}
    quote
        L = NTuple{d,Int}[]
        $(Symbol("i_",d+1)) = 1
        @nloops $d i (j->i_{j+1}:N) begin
            push!(L, @ntuple $d j->i_j)
        end
        return L
    end
end
gentups(d::Integer, N::Integer) = gentups(Val(Int(d)), Int(N))

```

so that we get, for example:

```julia
julia> gentups(3, 5)
35-element Array{Tuple{Int64,Int64,Int64},1}:
 (1, 1, 1)
 (2, 1, 1)
 (3, 1, 1)
 (4, 1, 1)
 (5, 1, 1)
 (2, 2, 1)
 (3, 2, 1)
 (4, 2, 1)
 (5, 2, 1)
 (3, 3, 1)
 (4, 3, 1)
 (5, 3, 1)
 (4, 4, 1)
 ⋮        
 (5, 4, 2)
 (5, 5, 2)
 (3, 3, 3)
 (4, 3, 3)
 (5, 3, 3)
 (4, 4, 3)
 (5, 4, 3)
 (5, 5, 3)
 (4, 4, 4)
 (5, 4, 4)
 (5, 5, 4)
 (5, 5, 5)

```

(I’m not sure if this is the output you want, or if you also want all permutations of each tuple?)

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 28, 2018, 3:37am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/8 "2018-11-28T03:37:36Z")

</div>

@jandehaan Oh nice. That seems to be it. We build all the m-tuple = (x\_1, …, x\_m) (where m = dimension) whose entry are the integer from 0 to N (you chose 0 to N+1 in your function) satisfying x\_1 + … + x\_m \<= N.  
The geometry approach did not seem to give anything pretty anyways.

@stevengj I have to look more precisely, but it does not seem to give the right tuples. If I set N= 2 and d = 3 I got an empty array. However, thank you for showing me the @nloops macro! That will definitely come in handy.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [November 28, 2018, 3:40am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/9 "2018-11-28T03:40:33Z")

</div>

> [@OlivierHnt](#):
>
> but it does not seem to give the right tuples. If I set N= 2 and d = 3 I got an empty array.

I had a bug, which is now corrected above. It now gives:

```julia
julia> gentups(3,2)
4-element Array{Tuple{Int64,Int64,Int64},1}:
 (1, 1, 1)
 (2, 1, 1)
 (2, 2, 1)
 (2, 2, 2)

```

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 28, 2018, 3:47am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/10 "2018-11-28T03:47:34Z")

</div>

Indeed, as you said, it is just missing the permutations to get the full structure! I am looking into adding this into your code.

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 28, 2018, 5:01am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/11 "2018-11-28T05:01:23Z")

</div>

Ok so there was a slight issue with the code. From my understanding, your code does

```julia
for i_3 in 1:N
    for i_2 in i_3:N
        for i_1 in i_2:N

```

But, I think the loops needed are

```julia
for i_3 in 1:N
    for i_2 in 1:N-i_3+1
        for i_1 in 1:N-i_2-i_3+2
            println((i_1, i_2, i_3))
        end
    end
end

```

Would you know how to pass this to `@nloops` ?

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [November 28, 2018, 11:21am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/12 "2018-11-28T11:21:43Z")

</div>

I think by far the easiest is

```julia
tf = [sum(I) <= N+1 for I in Iterators.product(ntuple(d->0:N+1, D)...)]

```

where `D` is the dimensionality.

---

<div class="post-metadata">

**Author:** ![jandehaan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jandehaan/32/6805_2.png) [@jandehaan](https://discourse.julialang.org/u/jandehaan)\
**Post date:** [November 28, 2018, 11:54am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/13 "2018-11-28T11:54:52Z")

</div>

> [@tim.holy](#):
>
> tf = [sum(I) \<= N+1 for I in Iterators.product(ntuple(d-\>0:N+1, D)…)]

Now that is elegant!  
To get the individual tuples:

```julia
julia> makenodes(N,D) = [I for I in Iterators.product(ntuple(d->0:N, D)...) if sum(I) <= N]
julia>makenodes(5,3)
84-element Array{Tuple{Int64,Int64,Int64},1}:
 (0, 0, 0)
 (1, 0, 0)
 (2, 0, 0)
 (3, 0, 0)
 (4, 0, 0)
 (5, 0, 0)
 (6, 0, 0)
 (0, 1, 0)
 (1, 1, 0)
 (2, 1, 0)
 (3, 1, 0)
 (4, 1, 0)
 (5, 1, 0)
 ⋮    

```

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 28, 2018, 3:29pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/14 "2018-11-28T15:29:11Z")

</div>

Indeed, that is very nice! I am going to use it right. Thank you very much 🙂

Just a note, with this procedure we end up building a very large object to only use a tiny part of it (especially as the dimension increases).  
So the approach used by @stevengj should be the most memory efficient (provided one can fix the loop situation mentioned above). So if we take the 3D case it should be like

```julia
function gentuples(N)
    L = NTuple{3, Int}[]
    for i_3 in 1:N
        for i_2 in 1:N-i_3+1
            for i_1 in 1:N-i_2-i_3+2
                push!(L, tuple(i_1,i_2,i_3))
            end
        end
    end
    return L
end

```

I am surprised that the version proposed by @tim.holy (using the adapted form suggested by @jandehaan) is much faster than the loops and is even more memory efficient. What would explain this behaviour ?

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [November 28, 2018, 3:47pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/15 "2018-11-28T15:47:33Z")

</div>

> [@OlivierHnt](#):
>
> What would explain this behaviour ?

Depends on what you are benchmarking, but when you provide a constant `D` and `N`, the compiler will unroll accordingly. Eg look at

```julia
@code_typed makenodes(3, 5)

```

Note that when `N` and `D` are only known at runtime, you will not get this.

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 28, 2018, 4:00pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/16 "2018-11-28T16:00:38Z")

</div>

I am not sure I know what I am looking at, I will do some readings about `@code_typed`. I simply used `@time` before calling the two functions (choosing `N=500` for instance).  
Could you clarify what happens if I do not declare `N` and `D` as constants and I call the function `makenodes` ?

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [November 28, 2018, 4:42pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/17 "2018-11-28T16:42:31Z")

</div>

I am surprised your `gentuples` version isn’t fast in 3d. Are you sure you timed it correctly, not including the compilation time? (Using BenchmarkTools.jl is recommended.)

For one in arbitrary dimensionality, the question is whether the compiler knows the dimensionality. You can use the `Val` construct illustrated above (see the help `?Val`) to pass them as compile-time constants.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [November 29, 2018, 12:14am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/18 "2018-11-29T00:14:04Z")

</div>

> [@OlivierHnt](#):
>
> Would you know how to pass this to `@nloops` ?

Maybe you want:

```julia
using Base.Cartesian
@generated function gentups(::Val{d}, N::Integer) where {d}
    quote
        L = NTuple{d,Int}[]
        $(Symbol("s_",d+1)) = 0
        @nloops $d i (j->1:(j == $d ? N : N-s_{j+1}+1)) (j->s_{j}=i_{j}+s_{j+1}) begin
            push!(L, @ntuple $d j->i_j)
        end
        return L
    end
end
gentups(d::Integer, N::Integer) = gentups(Val(Int(d)), Int(N))

```

If you do `gentups(3, 5)`, all of the tuples sum to `≤ 6` as I think you want.

---

<div class="post-metadata">

**Author:** ![OlivierHnt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/olivierhnt/32/6227_2.png) [@OlivierHnt](https://discourse.julialang.org/u/OlivierHnt)\
**Post date:** [November 29, 2018, 5:58am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/19 "2018-11-29T05:58:34Z")

</div>

@tim.holy using BenchmarkTools.jl gave better results indeed. I will do some reading about `Val` then! In all cases, your approach works really well.

@stevengj it does seem to work but I have to think twice before using it. Given a N I would expect to have, in particular, the tuples `(1,...,1)` all the way down to `(N,1, ...,1)`. While

```julia
julia> gentups(2,3)
6-element Array{Tuple{Int64,Int64},1}:
 (1, 1)
 (2, 1)
 (3, 1)
 (1, 2)
 (2, 2)
 (1, 3)

```

gives the expected result, if I increase the dimension we get

```julia
gentups(3,3)
4-element Array{Tuple{Int64,Int64,Int64},1}:
 (1, 1, 1)
 (2, 1, 1)
 (1, 2, 1)
 (1, 1, 2)

```

so I need to also increase N to have what I want, that is

```julia
julia> gentups(3,4)
10-element Array{Tuple{Int64,Int64,Int64},1}:
 (1, 1, 1)
 (2, 1, 1)
 (3, 1, 1)
 (1, 2, 1)
 (2, 2, 1)
 (1, 3, 1)
 (1, 1, 2)
 (2, 1, 2)
 (1, 2, 2)
 (1, 1, 3)

```

May I ask some technical questions ?  
Why did you choose to make this function into a macro ? I have done some readings about macros but it is still not clear. Is it just because it is convenient to use macro sometimes ?  
Also, I am not sure what the line `gentups(d::Integer, N::Integer) = gentups(Val(Int(d)), Int(N))` at the end really does.

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [November 29, 2018, 9:24am UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/20 "2018-11-29T09:24:06Z")

</div>

You’ll find more answers in the documentation on [`Base.Cartesian`](https://docs.julialang.org/en/latest/devdocs/cartesian/#Base.Cartesian-1) (it seems a little out of date, it is no longer “one of the few ways to write compact and performant multidimensional code.”)

The answer to your question about a macro is that the number of loops is a variable (`d`), and hence if you were typing it out manually you’d need to write different versions for different values of `d`. The macros of `Base.Cartesian` automate this process—macros operate on code itself rather than on values, and so they are good for this kind of application. The `Val` construct, and making the function `@generated`, ensures that `d` is only used at compile time; the compiled function will, in the end, look essentially exactly like the version you’d type out by hand for each dimensionality individually.

Even though I wrote `Base.Cartesian`, I recommend _against_ using it unless you have a compelling reason—the code is harder to understand, and generated functions have their disadvantages. However, this might be such a case: it would depend on whether the comprehension solution is in practice a performance bottleneck for you. If it isn’t, I would personally recommend keeping it simple, but as you point out in principle you can beat it by limiting the looping to just the values you plan to keep, and the `Base.Cartesian` solution can be leveraged to do that.

---

<div class="post-metadata">

**Author:** ![heflwi](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/heflwi/32/6267_2.png) [@heflwi](https://discourse.julialang.org/u/heflwi)\
**Post date:** [November 29, 2018, 2:25pm UTC](https://discourse.julialang.org/t/elegant-way-to-construct-a-set-of-multi-index/18081/21 "2018-11-29T14:25:14Z")

</div>

It seems what you want is for tuple sums to be `≤ N+d-1`, for which the `@nloops` line in `gentups` would need to be tweaked:

```julia
using Base.Cartesian
@generated function gentups(::Val{d}, N::Integer) where {d}
    quote
        L = NTuple{d,Int}[]
        $(Symbol("s_",d+1)) = 0
        @nloops $d i (j->1:(j == $d ? N : N-s_{j+1}+$d-1)) (j->s_{j}=i_{j}+s_{j+1}) begin
            push!(L, @ntuple $d j->i_j)
        end
        return L
    end
end
gentups(d::Integer, N::Integer) = gentups(Val(Int(d)), Int(N))

```

This then gives

```julia
julia> gentups(3,3)
10-element Array{Tuple{Int64,Int64,Int64},1}:
 (1, 1, 1)
 (2, 1, 1)
 (3, 1, 1)
 (1, 2, 1)
 (2, 2, 1)
 (1, 3, 1)
 (1, 1, 2)
 (2, 1, 2)
 (1, 2, 2)
 (1, 1, 3)

julia> gentups(4,2)
5-element Array{NTuple{4,Int64},1}:
 (1, 1, 1, 1)
 (2, 1, 1, 1)
 (1, 2, 1, 1)
 (1, 1, 2, 1)
 (1, 1, 1, 2)

```

as I guess you’re expecting.