With_repetitions_permutations v0.6 VS v1

Hello,

I am trying to update to v1.0 the code with_repetitions_permutations found on
rosettacode that was build for Julia v0.6.

This function simply gives all possible permutations with repetition (function that I did not found in the Combinatorics package).

The two codes works on their respective versions but somehow the v0.6 version is a bit faster (has less allocations) when I compare them with @time.
For example (on the same computer)
v0.6:

@time collect(with_repetitions_permutations([0,1,2],15))
  4.427327 seconds (43.05 M allocations: 3.528 GiB, 56.30% gc time)
14348907-element Array{Array{Int64,1},1}:

 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 [1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 [2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 [0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 ⋮
 [1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
 [2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]

For the v1.0:

@time collect(with_repetitions_permutations([0,1,2],15))
  7.878192 seconds (57.40 M allocations: 6.308 GiB, 63.05% gc time)
14348907-element Array{Array{Int64,1},1}:
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 [1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 [2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 [0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 ⋮
 [1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
 [2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2]

The original code v0.6 is

struct WithRepetitionsPermutations{T}
    a::T
    t::Int
end

with_repetitions_permutations(elements::T, len::Integer) where T =
    WithRepetitionsPermutations{T}(unique(elements), len)

Base.iteratorsize(::WithRepetitionsPermutations) = Base.HasLength()
Base.length(p::WithRepetitionsPermutations) = length(p.a) ^ p.t
Base.iteratoreltype(::WithRepetitionsPermutations) = Base.HasEltype()
Base.eltype(::WithRepetitionsPermutations{T}) where T = T
Base.start(p::WithRepetitionsPermutations) = ones(Int, p.t)
Base.done(p::WithRepetitionsPermutations, s::Vector{Int}) = s[end] > endof(p.a)
function Base.next(p::WithRepetitionsPermutations, s::Vector{Int})
    cur = p.a[s]
    s[1] += 1
    local i = 1
    while i < endof(s) && s[i] > length(p.a)
        s[i] = 1
        s[i+1] += 1
        i += 1
    end
    return cur, s
end

while the code for v1.0 that I wrote while being inspired by the Writing Iterators in Julia 0.7 basically I gave the new name to the iterator and insert the three method start, next, done in the main part.

struct WithRepetitionsPermutations{T}
    a::T
    t::Int
end

with_repetitions_permutations(elements::T, len::Integer) where T =
    WithRepetitionsPermutations{T}(unique(elements), len)

Base.IteratorSize(::WithRepetitionsPermutations) = Base.HasLength()
Base.length(p::WithRepetitionsPermutations) = length(p.a) ^ p.t
Base.IteratorEltype(::WithRepetitionsPermutations) = Base.HasEltype()
Base.eltype(::WithRepetitionsPermutations{T}) where T = T


function Base.iterate(p::WithRepetitionsPermutations, state=ones(Int, p.t))
    if s[end] > lastindex(p.a)
        return nothing
    end

    cur = p.a[s]
    s[1] += 1
    local i = 1
    while i < lastindex(s) && s[i] > length(p.a)
        s[i] = 1
        s[i+1] += 1
        i += 1
    end
    return cur, s
end

My question is then did I did something wrong or not optimal in re writing the code that would explain the speed (allocation) difference? I am very new to Julia and I am not entirely sure to understand this iterator function so it very possible.
Thanks

edit: The problem here was that state should have been replaced by s

Is the question why your two versions differ in performance, or how to solve this in an efficient way? This is how I’d do it:

julia> it = Iterators.product(ntuple(_ -> 0:2, 15)...)
Base.Iterators.ProductIterator{NTuple{15,UnitRange{Int64}}}((0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2, 0:2))

julia> @time p = collect(it);
  0.879027 seconds (7 allocations: 1.604 GiB, 0.77% gc time)

julia> length(p)
14348907

The question was first about the difference between the two.
But I am also interested by your code that is very fast.
BTW with your code how can I simply (and quickly) go from

3×3×3 Array{Tuple{Int64,Int64,Int64},3}:
[:, :, 1] =
 (0, 0, 0)  (0, 1, 0)  (0, 2, 0)
 (1, 0, 0)  (1, 1, 0)  (1, 2, 0)
 (2, 0, 0)  (2, 1, 0)  (2, 2, 0)

[:, :, 2] =
 (0, 0, 1)  (0, 1, 1)  (0, 2, 1)
 (1, 0, 1)  (1, 1, 1)  (1, 2, 1)
 (2, 0, 1)  (2, 1, 1)  (2, 2, 1)

[:, :, 3] =
 (0, 0, 2)  (0, 1, 2)  (0, 2, 2)
 (1, 0, 2)  (1, 1, 2)  (1, 2, 2)
 (2, 0, 2)  (2, 1, 2)  (2, 2, 2)

to

27-element Array{Array{Int64,1},1}:
 [0, 0, 0]
 [1, 0, 0]
 [0, 1, 0]
 [0, 0, 1]
 [1, 1, 0]
 [1, 0, 1]
 [0, 1, 1]
 [2, 0, 0]
 [0, 2, 0]
 [0, 0, 2]
 ⋮        
 [1, 2, 1]
 [1, 1, 2]
 [2, 2, 0]
 [2, 0, 2]
 [0, 2, 2]
 [2, 2, 1]
 [2, 1, 2]
 [1, 2, 2]
 [2, 2, 2]

here I show

Iterators.product(ntuple(_ -> 0:2, 3)...)

The ultimate goal being to sort the vectors as (0,0,0),(1,0,0),(0,1,0),(0,0,1),(1,1,0), etc

An easy way would be:

julia> collect.(Iterators.product(ntuple(_ -> 0:2, 3)...))[:]
27-element Array{Array{Int64,1},1}:
 [0, 0, 0]
 [1, 0, 0]
 [2, 0, 0]
 [0, 1, 0]

But this is not the intended usage of Iterators.product, nor what that rosettacode task is about. If you want to create and store all permutations, you can do better. The point of using an iterator is to lazily iterate over the elements, without needing to store them all.

With regards to your existing code, you have:

I think you want the argument to be named s, not state. This is accessing some global, which isn’t doing what you want (and is probably a type-instability).

@mbauman Thanks for your answer, indeed replacing state by s makes the program v1.0 a tiny bit faster than the v0.6!
@bennedich Thanks, now the version I have without the error in the code is faster than doing collect.(Iterators.product(ntuple(_ -> 0:2, 3)...))[:].
I am not sure to understand the comparison with the rosettacode, both are showing all permutations with repetitions no?

You most likely should not be using collect here. At the top of that rosettacode page, you’ll see this:

Do not store all the intermediate values of the sequence, rather generate them as required, and pass the intermediate result to a deciding routine for combinations selection and/or early generator termination.

When you call collect, you force the immediate generation and storage of all values. If that’s indeed what you want to do, there are other ways of solving this problem that may be more efficient than using an iterator.

So what’s the point of not generating them all? Well, let’s say that you want to find out if there’s any permutation among all of these 14 million where the digits add up to 15. You can do it like this:

julia> it = Iterators.product(ntuple(_ -> 0:2, 15)...);

julia> any(v -> sum(v) == 15, it)
true

Let’s test performance.

julia> using BenchmarkTools

julia> @btime any(v -> sum(v) == 15, $it);
  10.179 μs (0 allocations: 0 bytes)

So instead of taking several seconds, and using many GB of memory, you got the answer in 10 microseconds, and didn’t allocate a single byte! That is possible since the iterator lazily creates one permutation at a time and tests it, and once any permutation passes the test (digits add up to 15), the iteration is stopped.

1 Like

I see, but indeed here I want to generate and sort all the possible vectors and save them in a dictionary.
So even for sorting I have to use collectat some point (because sort needs arrays).
I don’t know enough Julia to think of more efficient methods.

Perhaps you could describe your use case in a bit more detail. How many permutations will you have? (Enough that performance will be an issue?) What will they consist of? How do you want to sort them? What are the keys and values in that dictionary? How will you use the dictionary, and why do you need it sorted?

What you’ve described so far makes little sense to me. First, sorting is a fairly slow operation compared to generating the permutations (it’s not linear in time), so you’re better off creating your permutations in the correct order to begin with, if possible.

Second, dictionaries are not ordered in Julia. You can use OrderedDict or SortedDict from DataStructures.jl if you need it ordered or sorted, but I would first make sure that this is indeed the correct data structure for the problem you’re solving.

1 Like

Ideally I would like to store and sort the permutation with replacement of 0:p in a size N-set.
For example in a dictionary it would look for p=2 and N=2
0 =>(0,0), 1 =>(1,0), 2 =>(0,1), 3 =>(2,0), 4 =>(0,2), 5 =>(1,1), 6 =>(2,1), 7 =>(1,2), 8 =>(2,2)
I know that dictionaries are disorganized but I would like to assign to each vector a label (going from 0 to (p+1)^N) in accordance to the ‘natural’ sorting (norm 1 of the vector).

I don’t know the most efficient and Julia-like way to do that…

Why does this have to be a dictionary, and not just an array, where the keys would correspond to the array index (1-based instead of 0-based)? A dictionary is slower and uses more space than an array.

I.e. the sum of all elements? How would you break ties? In this example (2,0) comes before (1,1), but in your previous example, (1,1,0) came before (2,0,0).

As I wrote in my previous post, if the performance of the current approach is unacceptable, I would try to generate the permutations in the correct order to begin with (in a pre-allocated array), instead of first creating them all, then sorting them.

1 Like

Thanks for all the comments. I think I found the best way to do that, it was already implemented in the multiexponents(m,n) of the Combinaorics package. For example (taken from the package page)

 multiexponents(m, n)
Returns the exponents in the multinomial expansion (x₁ + x₂ + ... + xₘ)ⁿ.
For example, the expansion (x₁ + x₂ + x₃)² = x₁² + x₁x₂ + x₁x₃ + ...
has the exponents:
    julia> collect(multiexponents(3, 2))
    6-element Array{Any,1}:
     [2, 0, 0]
     [1, 1, 0]
     [1, 0, 1]
     [0, 2, 0]
     [0, 1, 1]
     [0, 0, 2]

It displays and sorts exactly what I was looking for in a very fast way.

If this solves your problem – great! But I’m a bit confused. Multiexponents is not the same as permutations with repetition, or am I missing something? How would you generate the example set of 27 permutations of 0:2 that you posted above, using multiexponents? And what about the dictionary?