Juliaで順列・組み合わせ・円順列・数珠順列を求めてみる！

はじめに

TwitterでZennのJulia関連のものが流れてきました。

https://zenn.dev/ohno/articles/03e65bfa028baa

順列と組み合わせで全てのリストを作成するプログラムでした。（Shuhei Ohnoさん）

円順列と数珠順列で『誰かコメントください！』とあったので，『これはやらねば！』と思い作ってみました。

リストが同じものを含む場合

Shuhei Ohnoさんのコードを見てみて思ったのですが，[1,1,2,3,4]などの同じものを含む順列・組み合わせなどには対応していなかったので，そちらも確認したいと思います。

順列（1列に並べる）

まずは，パッケージCombinatorics.jlを使います。

順列（同じものを含まない）

using Combinatorics
seq = [1, 2, 3, 4, 5]
p = union(permutations(seq, 3)) 
# 順列（同じものでもOK）この例は集合から3個とって1列に並べる。

60-element Vector{Vector{Int64}}:
 [1, 2, 3]
 [1, 2, 4]
 [1, 2, 5]
 [1, 3, 2]
 [1, 3, 4]
 [1, 3, 5]
 [1, 4, 2]
 [1, 4, 3]
 [1, 4, 5]
 [1, 5, 2]
 ⋮
 [5, 2, 1]
 [5, 2, 3]
 [5, 2, 4]
 [5, 3, 1]
 [5, 3, 2]
 [5, 3, 4]
 [5, 4, 1]
 [5, 4, 2]
 [5, 4, 3]

個数も表示されて，60個であることがわかります！

順列（同じものを含む）

同じものを含む順列でもやってみましょう。

seq2 = [1, 1, 1, 2, 2, 3]
p2 = union(permutations(seq2, 3)) 

19-element Vector{Vector{Int64}}:
 [1, 1, 1]
 [1, 1, 2]
 [1, 1, 3]
 [1, 2, 1]
 [1, 2, 2]
 [1, 2, 3]
 [1, 3, 1]
 [1, 3, 2]
 [2, 1, 1]
 [2, 1, 2]
 [2, 1, 3]
 [2, 2, 1]
 [2, 2, 3]
 [2, 3, 1]
 [2, 3, 2]
 [3, 1, 1]
 [3, 1, 2]
 [3, 2, 1]
 [3, 2, 2]

19個です。

組み合わせ

続いて，組み合わせです。

組み合わせ（同じものを含まない）

seq = [1, 2, 3, 4, 5]
c = union(combinations(seq, 3)) 
# 組み合わせ（同じものでもOK）この例は集合から3個とってくる組み合わせ。

10-element Vector{Any}:
 [1, 2, 3]
 [1, 2, 4]
 [1, 2, 5]
 [1, 3, 4]
 [1, 3, 5]
 [1, 4, 5]
 [2, 3, 4]
 [2, 3, 5]
 [2, 4, 5]
 [3, 4, 5]

${}_5\text{C}_3=10$個ですね。

組み合わせ（同じものを含む）

それでは，同じもの含む場合です。

seq2 = [1, 1, 1, 2, 2, 3]
c = union(combinations(seq2, 3)) 
# 組み合わせ（同じものでもOK）この例は集合から3個とってくる組み合わせ。

6-element Vector{Any}:
 [1, 1, 1]
 [1, 1, 2]
 [1, 1, 3]
 [1, 2, 2]
 [1, 2, 3]
 [2, 2, 3]

6個です。高校数学の教科書などでは例題で『樹形図で考えてみよう！』となってます。

円順列

さて，円順列です。私自身かつて，『同じもの含む円順列の総数』は公式を駆使して計算したことがあります。

https://zenn.dev/dannchu/articles/2a27581585a338

今回は実際にそのリストを求めます。前回の結果とも一致するか確認したいです。

今回の方針は

ターゲットの列に対して，それ以下の列をシフトしてマッチしたものを取り除いていく。

です。適宜，unionで重複を除いています。

円順列（同じものを含まない）

# 円順列（同じものでもOK）この例は集合から3個とってくる円順列。
function circperm(seq, k)
    p = union(permutations(seq, k))
    n = length(p)
    d = []
    for i = 1:n-1, j = i+1:n, t = 1:k-1
        if p[i] == circshift(p[j], t)
            push!(d, j)
        end
    end
    deleteat!(p, sort!(union!(d)))
end

circperm(seq, 3)

20-element Vector{Vector{Int64}}:
 [1, 2, 3]
 [1, 2, 4]
 [1, 2, 5]
 [1, 3, 2]
 [1, 3, 4]
 [1, 3, 5]
 [1, 4, 2]
 [1, 4, 3]
 [1, 4, 5]
 [1, 5, 2]
 [1, 5, 3]
 [1, 5, 4]
 [2, 3, 4]
 [2, 3, 5]
 [2, 4, 3]
 [2, 4, 5]
 [2, 5, 3]
 [2, 5, 4]
 [3, 4, 5]
 [3, 5, 4]

${}_5\text{C}_3\times \dfrac{3!}{3}=20$個ですね。

円順列（同じものを含む）

seq2 = [1, 1, 1, 2, 2,3]
circperm(seq2, 3)

7-element Vector{Vector{Int64}}:
 [1, 1, 1]
 [1, 1, 2]
 [1, 1, 3]
 [1, 2, 2]
 [1, 2, 3]
 [1, 3, 2]
 [2, 2, 3]

7個です。

念の為，以前作ったものでも確かめてみます。

seq3 = [1, 1, 1, 1, 2,2,2,2,2,2]
circperm(seq3, 10)

22-element Vector{Vector{Int64}}:
 [1, 1, 1, 1, 2, 2, 2, 2, 2, 2]
 [1, 1, 1, 2, 1, 2, 2, 2, 2, 2]
 [1, 1, 1, 2, 2, 1, 2, 2, 2, 2]
 [1, 1, 1, 2, 2, 2, 1, 2, 2, 2]
 [1, 1, 1, 2, 2, 2, 2, 1, 2, 2]
 [1, 1, 1, 2, 2, 2, 2, 2, 1, 2]
 [1, 1, 2, 1, 1, 2, 2, 2, 2, 2]
 [1, 1, 2, 1, 2, 1, 2, 2, 2, 2]
 [1, 1, 2, 1, 2, 2, 1, 2, 2, 2]
 [1, 1, 2, 1, 2, 2, 2, 1, 2, 2]
 ⋮
 [1, 1, 2, 2, 1, 2, 2, 1, 2, 2]
 [1, 1, 2, 2, 1, 2, 2, 2, 1, 2]
 [1, 1, 2, 2, 2, 1, 1, 2, 2, 2]
 [1, 1, 2, 2, 2, 1, 2, 1, 2, 2]
 [1, 1, 2, 2, 2, 1, 2, 2, 1, 2]
 [1, 1, 2, 2, 2, 2, 1, 2, 1, 2]
 [1, 2, 1, 2, 1, 2, 1, 2, 2, 2]
 [1, 2, 1, 2, 1, 2, 2, 1, 2, 2]
 [1, 2, 1, 2, 2, 1, 2, 1, 2, 2]

22個です。

以前作ったのコードで確かめてみると，

using Combinatorics #multinomial(a)多項係数を求める
using Primes    #totient(i)　オイラーのトーシェント関数

function divisors(n)    #約数のリストを求める関数
   X=[]
   for i=1:n
       if n % i==0
           X =push!(X,i)
       end
   end
   X
end

function enkan(a)
   l=gcd(a)    #リスト（配列）aの最大公約数を求める。
   N=sum(a)    #aの総和
   A=divisors(l)
   p=0
   for k in A
       q=map(x -> x÷k,a) 
       p +=totient(k)*multinomial(q...)
   end
   p÷N
end

enkan([4 6])

22

バッチリです！

数珠順列

最近の若い人は数珠（じゅず）って知らないよね。仏教でお葬式のときに手に持っているものです。

今回の方針は

円順列のリストを作っておいて，ターゲットの列に対して，それ以下の列をリバースしてマッチしたものを取り除いていく。

です。こちらも適宜，unionで重複を除いています。

数珠順列（同じものを含まない）

# 数珠順列（同じものでもOK）この例は集合から3個とってくる数珠順列。
function ringperm(seq, k)
   p = circperm(seq, k)
   n = length(p)
   d = []
   for i = 1:n-1, j = i+1:n, t = 1:k-1
       if p[i] == circshift(reverse(p[j]),t)
           push!(d, j)
       end
   end
   deleteat!(p, sort!(union!(d)))
end

ringperm(seq, 3)

10-element Vector{Vector{Int64}}:
 [1, 2, 3]
 [1, 2, 4]
 [1, 2, 5]
 [1, 3, 4]
 [1, 3, 5]
 [1, 4, 5]
 [2, 3, 4]
 [2, 3, 5]
 [2, 4, 5]
 [3, 4, 5]

${}_5\text{C}_3\times \dfrac{3!}{3}\times2=10$個ですね。

数珠順列（同じものを含む）

あまりやったことはありません。数が多いと考えるだけで大変です。

普通は裏返して，同じになる並びにペアを考えて2で割ります。同じものを含まない場合は必ずペアが存在するので，単純に2で割ればいいですね。

seq2 = [1, 1, 1, 2, 2,3]
ringperm(seq2, 3)

6-element Vector{Vector{Int64}}:
 [1, 1, 1]
 [1, 1, 2]
 [1, 1, 3]
 [1, 2, 2]
 [1, 2, 3]
 [2, 2, 3]

6個です。円順列は7個あったのに，全然減らないですね。7個中5個は裏返しても自分自身で，
[1, 2, 3]と[1, 3, 2]だけがペアとなります。

もう一つ例をやってみます。

seq3 = [1, 1, 1, 1, 2,2,2,2,2,2]
ringperm(seq3, 10)

16-element Vector{Vector{Int64}}:
 [1, 1, 1, 1, 2, 2, 2, 2, 2, 2]
 [1, 1, 1, 2, 1, 2, 2, 2, 2, 2]
 [1, 1, 1, 2, 2, 1, 2, 2, 2, 2]
 [1, 1, 1, 2, 2, 2, 1, 2, 2, 2]
 [1, 1, 2, 1, 1, 2, 2, 2, 2, 2]
 [1, 1, 2, 1, 2, 1, 2, 2, 2, 2]
 [1, 1, 2, 1, 2, 2, 1, 2, 2, 2]
 [1, 1, 2, 1, 2, 2, 2, 1, 2, 2]
 [1, 1, 2, 1, 2, 2, 2, 2, 1, 2]
 [1, 1, 2, 2, 1, 1, 2, 2, 2, 2]
 [1, 1, 2, 2, 1, 2, 1, 2, 2, 2]
 [1, 1, 2, 2, 1, 2, 2, 1, 2, 2]
 [1, 1, 2, 2, 2, 1, 1, 2, 2, 2]
 [1, 2, 1, 2, 1, 2, 1, 2, 2, 2]
 [1, 2, 1, 2, 1, 2, 2, 1, 2, 2]
 [1, 2, 1, 2, 2, 1, 2, 1, 2, 2]

16個です。円順列は22個でした。