diff --git a/base/combinatorics.jl b/base/combinatorics.jl index aa113500dc2be..5041b8044b0ff 100644 --- a/base/combinatorics.jl +++ b/base/combinatorics.jl @@ -252,110 +252,202 @@ end done(p::Permutations, s) = !isempty(s) && s[1] > length(p.a) -# Algorithm H from TAoCP 7.2.1.4 -# Partition n into m parts -function integer_partitions{T<:Integer}(n::T, m::T) - if n < m || m < 2 - throw("Require n >= m >= 2!") - end - # H1 - a = [n - m + 1, ones(T, m), -1] - # H2 - while true - produce(a[1:m]) - if a[2] < a[1] - 1 - # H3 - a[1] = a[1] - 1 - a[2] = a[2] + 1 - continue # to H2 +# Integer Partitions + +immutable IntegerPartitions + n::Int +end + +length(p::IntegerPartitions) = npartitions(p.n) + +partitions(n::Integer) = IntegerPartitions(n) + +start(p::IntegerPartitions) = Int[] +done(p::IntegerPartitions, xs) = length(xs) == p.n +next(p::IntegerPartitions, xs) = (xs = nextpartition(p.n,xs); (xs,xs)) + +function nextpartition(n, as) + if isempty(as); return Int[n]; end + + xs = similar(as,0) + sizehint(xs,length(as)+1) + + for i = 1:length(as)-1 + if as[i+1] == 1 + x = as[i]-1 + push!(xs, x) + n -= x + while n > x + push!(xs, x) + n -= x + end + push!(xs, n) + + return xs + end + push!(xs, as[i]) + n -= as[i] end - # H4 - j = 3 - s = a[1] + a[2] - 1 - if a[j] >= a[1] - 1 - while true - s = s + a[j] - j = j + 1 - if a[j] < a[1] - 1 - break # end loop + push!(xs, as[end]-1) + push!(xs, 1) + + xs +end + +const _npartitions = (Int=>Int)[] +function npartitions(n::Int) + if n < 0 + 0 + elseif n < 2 + 1 + elseif (np = get(_npartitions, n, 0)) > 0 + np + else + np = 0 + sgn = 1 + for k = 1:n + np += sgn * (npartitions(n-k*(3k-1)>>1) + npartitions(n-k*(3k+1)>>1)) + sgn = -sgn end - end + _npartitions[n] = np end - # H5 - if j > m - break # terminate +end + + +# Algorithm H from TAoCP 7.2.1.4 +# Partition n into m parts +# in colex order (lexicographic by reflected sequence) + +immutable FixedPartitions + n::Int + m::Int +end + +length(f::FixedPartitions) = npartitions(f.n,f.m) + +partitions(n::Integer, m::Integer) = (@assert 2 <= m <= n; FixedPartitions(n,m)) + +start(f::FixedPartitions) = Int[] +done(f::FixedPartitions, s::Vector{Int}) = !isempty(s) && s[1]-1 <= s[end] +next(f::FixedPartitions, s::Vector{Int}) = (xs = nextfixedpartition(f.n,f.m,s); (xs,xs)) + +function nextfixedpartition(n, m, bs) + as = copy(bs) + if isempty(as) + # First iteration + as = [n-m+1, ones(Int, m-1)] + elseif as[2] < as[1]-1 + # Most common iteration + as[1] -= 1 + as[2] += 1 + else + # Iterate + local j + s = as[1]+as[2]-1 + for j = 3:m + if as[j] < as[1]-1; break; end + s += as[j] + end + x = as[j] += 1 + for k = j-1:-1:2 + as[k] = x + s -= x + end + as[1] = s end - x = a[j] + 1 - a[j] = x - j = j - 1 - # H6 - while j > 1 - a[j] = x - s = s - x - j = j - 1 + + return as +end + +const _nipartitions = ((Int,Int)=>Int)[] +function npartitions(n::Int,m::Int) + if n < m || m == 0 + 0 + elseif n == m + 1 + elseif (np = get(_nipartitions, (n,m), 0)) > 0 + np + else + _nipartitions[(n,m)] = npartitions(n-1,m-1) + npartitions(n-m,m) end - a[1] = s - end end + # Algorithm H from TAoCP 7.2.1.5 # Set partitions -function partitions{T}(s::AbstractVector{T}) - n = length(s) - n == 0 && return - if n == 1 - produce(Array{T,1}[s]) - return - end - - # H1 - a = zeros(Int,n) - b = ones(Int,n-1) - m = 1 - - while true - # H2 - # convert from restricted growth string a[1:n] to set of sets - temp = [ Array(T,0) for k = 1:n ] - for k = 1:n - push!(temp[a[k]+1], s[k]) - end - result = Array(Array{T,1},0) - for arr in temp - if !isempty(arr) - push!(result, arr) - end - end - #produce(a[1:n]) # this is the string representing set assignment - produce(result) - if a[n] != m - # H3 - a[n] = a[n] + 1 - continue # to H2 - end - # H4 - j = n - 1 - while a[j] == b[j] - j = j - 1 +immutable SetPartitions{T} + s::AbstractVector{T} +end + +length(s::SetPartitions) = npartitions(s) + +partitions(s::AbstractVector) = SetPartitions(s) + +start(p::SetPartitions) = (n = length(p.s); (zeros(Int32, n), ones(Int32, n-1), n, 1)) +done(p::SetPartitions, s) = !isempty(s) && s[1][1] > 0 +next(p::SetPartitions, s) = nextsetpartition(p.s, s...) + +function nextsetpartition(s, a, b, n, m) + function makeparts(s, a, m) + part = {} + for i = 0:m + ss = s[findin(a,i)] + if !isempty(ss) + push!(part, ss) + end + end + part end - # H5 - if j == 1 - break # terminate + + if isempty(s); return ({s}, ([1], Int[], n, 1)); end + + part = makeparts(s,a,m) + + if a[end] != m + a[end] += 1 + else + local j + for j = n-1:-1:1 + if a[j] != b[j] + break + end + end + a[j] += 1 + m = b[j] + (a[j] == b[j]) + for k = j+1:n-1 + a[k] = 0 + b[k] = m + end + a[end] = 0 end - a[j] = a[j] + 1 - # H6 - m = b[j] + (a[j] == b[j]) - j = j + 1 - while j < n - a[j] = 0 - b[j] = m - j = j + 1 + + return (part, (a,b,n,m)) + +end + + +npartitions(p::SetPartitions) = nsetpartitions(length(p.s)) +npartitions(v::AbstractVector) = nsetpartitions(length(v)) + +const _nsetpartitions = (Int=>Int)[] +function nsetpartitions(n::Int) + if n < 0 + 0 + elseif n < 2 + 1 + elseif (wn = get(_nsetpartitions, n, 0)) > 0 + wn + else + wn = 0 + for k = 0:n-1 + wn += binomial(n-1,k)*nsetpartitions(n-1-k) + end + _nsetpartitions[n] = wn end - a[n] = 0 - end end + # For a list of integers i1, i2, i3, find the smallest # i1^n1 * i2^n2 * i3^n3 >= x # for integer n1, n2, n3 diff --git a/base/deprecated.jl b/base/deprecated.jl index 867e553c7cbcf..69881bcddbbc3 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -376,3 +376,7 @@ function addprocs_sge(np::Integer) error("Base.addprocs_sge is discontinued - add package ClusterManagers and then use ClusterManagers.addprocs_sge instead.") end export addprocs_sge + +function integer_partitions(n,m) + error("integer_partitions(n,m) has been renamed to partitions(n,m), and is now an iterator. Please update your code.") +end diff --git a/base/exports.jl b/base/exports.jl index 27db8148b8533..8bcacf1aaedbc 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -517,6 +517,8 @@ export ndims, nnz, nonzeros, + npartitions, + nsetpartitions, nthperm!, nthperm, ones, diff --git a/doc/stdlib/base.rst b/doc/stdlib/base.rst index 3cbd2cbeda496..4f1cd36c82ddb 100644 --- a/doc/stdlib/base.rst +++ b/doc/stdlib/base.rst @@ -3156,7 +3156,7 @@ Combinatorics .. function:: reverse(v [, start=1 [, stop=length(v) ]] ) - Reverse vector ``v``, optionally from start to stop. + Return a copy of ``v`` reversed from start to stop. .. function:: reverse!(v [, start=1 [, stop=length(v) ]]) -> v @@ -3176,20 +3176,45 @@ Combinatorics iterator object. Use ``collect(permutations(a,n))`` to get an array of all permutations. -.. function:: integer_partitions(n, m) +.. function:: partitions(n) + + Generate all integer arrays that sum to ``n``. Because the number of + partitions can be very large, this function returns an iterator + object. Use ``collect(partitions(n))`` to get an array of all + partitions. + +.. function:: partitions(n, m) Generate all arrays of ``m`` integers that sum to ``n``. Because - the number of partitions can be very large, this function runs inside - a Task to produce values on demand. Write - ``c = @task integer_partitions(n,m)``, then iterate ``c`` or call - ``consume`` on it. + the number of partitions can be very large, this function returns an + iterator object. Use ``collect(partitions(n,m))`` to get an array of + all partitions. .. function:: partitions(array) - Generate all set partitions of the elements of an array, represented as - arrays of arrays. Because the number of partitions can be very large, this - function runs inside a Task to produce values on demand. Write - ``c = @task partitions(a)``, then iterate ``c`` or call ``consume`` on it. + Generate all set partitions of the elements of an array, + represented as arrays of arrays. Because the number of partitions + can be very large, this function returns an iterator object. Use + ``collect(partitions(array))`` to get an array of all partitions. + +.. function:: npartitions(n) + + Calculate the number of ways of writing the integer ``n`` as a sum + of positive integers. + +.. function:: npartitions(n, m) + + Calculate the number of ways of writing the integer ``n`` as a sum + of ``m`` positive integers. + +.. function:: npartitions(array) + + Calculate the number of ways to partition ``array``, assuming + unique elements. Equivalent to ``nsetpartitions(length(array))`` + +.. function:: nsetpartitions(n) + + Calculate the number of ways to partition a set of size ``n``. Statistics diff --git a/test/Makefile b/test/Makefile index 32bd7f26728df..d322696104e6c 100644 --- a/test/Makefile +++ b/test/Makefile @@ -6,7 +6,7 @@ TESTS = all core keywordargs numbers strings unicode collections hashing \ math functional bigint sorting statistics spawn parallel arpack file \ git pkg pkg2 resolve resolve2 suitesparse complex version pollfd mpfr \ broadcast socket floatapprox priorityqueue readdlm regex float16 \ - combinations + combinatorics $(TESTS) :: @$(PRINT_JULIA) $(call spawn,$(JULIA_EXECUTABLE)) ./runtests.jl $@ diff --git a/test/combinatorics.jl b/test/combinatorics.jl index 2c342730ad08d..e12acb32f3e4c 100644 --- a/test/combinatorics.jl +++ b/test/combinatorics.jl @@ -1,9 +1,18 @@ -@test factorial(7) = 5040 +@test factorial(7) == 5040 @test factorial(7,3) == 7*6*5*4 @test binomial(5,3) == 10 @test isperm(shuffle([1:1000])) a = randcycle(10) @test ipermute!(permute!([1:10], a),a) == [1:10] @test collect(combinations("abc",2)) == ["ab","ac","bc"] -@test collect(collect(permutations("abc"))) == ["abc","acb","bac","bca","cab","cba"] +@test collect(permutations("abc")) == ["abc","acb","bac","bca","cab","cba"] +@test collect(partitions(4)) == { [4], [3,1], [2,2], [2,1,1], [1,1,1,1] } +@test collect(partitions(8,3)) == { [6,1,1], [5,2,1], [4,3,1], [4,2,2], [3,3,2] } +@test collect(partitions([1,2,3])) == { {[1,2,3]}, {[1,2],[3]}, {[1,3],[2]}, {[1],[2,3]}, {[1],[2],[3]} } + +@test length(collect(partitions(30))) == npartitions(30) == length(partitions(30)) +@test length(collect(partitions(90,4))) == npartitions(90,4) == length(partitions(90,4)) +@test length(collect(partitions('a':'h'))) == npartitions('a':'h') == nsetpartitions(8) == length(partitions('a':'h')) + +