重複円順列の個数をJuliaで解決！

はじめに

2023年5月7日にJuliaのバージョンが1.9.0となりました！こういう時には何かやってみたくなりますね。
https://julialang.org/downloads/

現在，中学3年生の数学を教えています。分野や高校数学Aの「場合の数・確率」です。

Twitterなどでは定期的に円順列の話題が挙がってきます。その中で @ysmemoirsさんのツイートで次のようなものがありました。

https://twitter.com/ysmemoirs/status/1661278744630394880

これは，やらねば！と思い，まず手計算してツイートしました。

https://twitter.com/dannchu/status/1661335029895802883

なかなか検討が進まないので，Juliaで確認することにしました。

同じものを含む円順列

以前，サイトの方で同じものを含む円順列については，まとめました。

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

コード1

同じものの個数がわかっていれば，その円順列の個数を求めることができます。

enkan.jl

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

例：白玉4個，黒玉6個を円形に並べる総数

enkan.jl

enkan([4,6])

22

と求めることができます。

コード2

上記のコード1とは別に全ての順列のリストを求めることもできます。

circperm.jl

# seq（同じものでもOK）のリストからk個選んで円順列を作り，そのリストを返す。
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

例：白玉4個，黒玉6個を円形に並べる順列のリスト

circperm.jl

circperm([1,1,1,1,2,2,2,2,2,2],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]

重複円順列

10個の数字から重複を許して4個選び，円形に並べる

まず，Twitterであった，0〜9の10個の数字から重複を許して4個選び，円形に並べることを考えましす。前セクションのコード2を使えば，とりあえず全てのリストが得られるので，そちらを試しました。

circperm.jl

k = [0,1,2,3,4,5,6,7,8,9,
0,1,2,3,4,5,6,7,8,9,
0,1,2,3,4,5,6,7,8,9,
0,1,2,3,4,5,6,7,8,9]
circperm(k,4)

2530-element Vector{Vector{Int64}}:
[0, 1, 2, 3]
[0, 1, 2, 4]
[0, 1, 2, 5]
[0, 1, 2, 6]
[0, 1, 2, 7]
[0, 1, 2, 8]
[0, 1, 2, 9]
[0, 1, 2, 0]
[0, 1, 2, 1]
[0, 1, 2, 2]
⋮
[7, 9, 9, 7]
[7, 9, 9, 9]
[7, 7, 7, 7]
[8, 9, 8, 9]
[8, 9, 8, 8]
[8, 9, 9, 8]
[8, 9, 9, 9]
[8, 8, 8, 8]
[9, 9, 9, 9]

2530通りのようです。

10個の数字から重複を許して6個選び，円形に並べる

@ysmemoirsさんのツイートの続きを見ると，「授業の流れ上必要になったのは，実はこれの6桁版」とあったので，さらにやってみることにしました。

https://twitter.com/ysmemoirs/status/1661286233828409440

最初は手計算でしたものをツイートしました。

https://twitter.com/dannchu/status/1661335804701196288

さすがに，これは自信がないので，「誰かやってくれないかな。。。」と思っていたのですが，あまり反応もなかったので，Juliaで検算することにしました。しかし，残念ながら，コード2のcircpermでは時間がかかりすぎるので，コード1のenkanを活かすことにしました。分割の様子がわかれば，こちらが利用できそうです。

手計算で行った時のまとめ

まずは，10個の数字から重複を許して4個選び，円形に並べる場合の手計算で分類した時のものを分割を網羅する形でまとめてみました。

4の分割数を考える。

(4) 文字の選び方は10通り，そのそれぞれに対して並べ方は1通り。$10\times 1=10$

(1-3) 文字の選び方は${}_{10}\text{P}_{2}=90$通り，そのそれぞれに対して並べ方は1通り。$90\times 1=90$

(2-2) 文字の選び方は${}_{10}\text{C}_{2}=45$通り，そのそれぞれに対して並べ方は2通り。$45\times 2=90$

(1-1-2) 文字の選び方は$\dfrac{10\times9\times 8}{2}=360$通り，そのそれぞれに対して並べ方は$\dfrac{4!}{2!\times 4}=3$通り。$360\times 3=1080$

(1-1-1-1) 文字の選び方は${}_{10}\text{C}_{4}=210$通り，そのそれぞれに対して並べ方は$3!=6通り。$210\times 6=1260$

$$\therefore 10+90+90+1080+1260=2530通り$$

次に5個選ぶ場合も考えてみました。

5の分割数を考える。

(5) 文字の選び方は10通り，そのそれぞれに対して並べ方は1通り。$10\times 1=10$

(1-4) 文字の選び方は${}_{10}\text{P}_{2}=90$通り，そのそれぞれに対して並べ方は$\dfrac{{}_{5}\text{C}_{1}}{5}=1$通り。$90\times 1=90$

(2-3) 文字の選び方は${}_{10}\text{P}_{2}=90$通り，そのそれぞれに対して並べ方は$\dfrac{{}_{5}\text{C}_{2}}{5}=2$通り。$90\times 2=180$

(1-1-3) 文字の選び方は$\dfrac{10\times9\times 8}{2}=360$通り，そのそれぞれに対して並べ方は$\dfrac{5!}{3!\times 5}=4$通り。$360\times 4=1440$

(1-2-2) 文字の選び方は$\dfrac{10\times9\times 8}{2}=360$通り，そのそれぞれに対して並べ方は$\dfrac{5!}{2!2!\times 5}=6$通り。$360\times 6=2160$

(1-1-1-2) 文字の選び方は$\dfrac{10\times9\times 8\times 7}{3!}=840$通り，そのそれぞれに対して並べ方は$\dfrac{5!}{2!\times 5}=12$通り。$840\times 12=10080$

(1-1-1-1-1) 文字の選び方は$\dfrac{10\times9\times 8\times 7\times 6}{5!}=252$通り，そのそれぞれに対して並べ方は$\dfrac{5!}{5}=24$通り。$252\times 24=6048$

$$\therefore 10+90+180+1440+2160+10080+6048=20008通り$$

分割の中の数の最大公約数が5のときは1で割る。最大公約数が１のときは5で割る。

ツイートでの考え方は下記の方法になります。簡単に述べると，分割の中の数の最大公約数が１のときは円形に並べる方法は1列に並べる総数を個数で割れば求まります。最大公約数が1ではないときは色々調べる必要があります。

$$\dfrac{10}{1}+\dfrac{10^5-10}{5}=10+\dfrac{99990}5=10+19998=20008$$

6の分割数を考える。

(6)文字の選び方は$10$通り，そのそれぞれに対して並べ方は1通り。$10\times 1=10$

(1-5) 文字の選び方は${}_{10}\text{P}_{2}=90$通り，そのそれぞれに対して並べ方は1通り。$90\times 1=90$

(2-4) 文字の選び方は${}_{10}\text{P}_{2}=90$通り，そのそれぞれに対して並べ方はenkan([2,4])$=3$通り。$90\times 3=270$

(3-3) 文字の選び方は${}_{10}\text{C}_{2}=45$通り，そのそれぞれに対して並べ方はenkan([3,3])$=4$通り。$45\times 4=180$

(1-1-4) 文字の選び方は$\dfrac{10\times9\times 8}{2!}=360$通り，そのそれぞれに対して並べ方はenkan([1,1,4])$=5$通り。$360\times 5=1800$

(1-2-3) 文字の選び方は$10\times9\times 8=720$通り，そのそれぞれに対して並べ方はenkan([1,2,3])$=10$通り。$720\times 10=7200$

(2-2-2) 文字の選び方は$\dfrac{10\times9\times 8}{3!}=120$通り，そのそれぞれに対して並べ方はenkan([2,2,2])$=16$通り。$120\times 16=1920$

(1-1-1-3) 文字の選び方は$\dfrac{10\times9\times 8\times 7}{3!}=840$通り，そのそれぞれに対して並べ方はenkan([1,1,1,3])$=20$通り。$840\times 20=16800$

(1-1-2-2) 文字の選び方は$\dfrac{10\times9\times 8\times 7}{2!2!}=1260$通り，そのそれぞれに対して並べ方はenkan([1,1,2,2])$=30$通り。$1260\times 30=37800$

(1-1-1-1-2) 文字の選び方は$\dfrac{10\times9\times 8\times 7\times 6}{4!}=1260$通り，そのそれぞれに対して並べ方はenkan([1,1,1,1,2])$=60$通り。$1260\times 60=75600$

(1-1-1-1-1-1) 文字の選び方は$\dfrac{10\times9\times 8\times 7\times 6\times 5}{6!}=210$通り，そのそれぞれに対して並べ方はenkan([1,1,1,1,1,1])$=120$通り。$210\times 120=25200$

$$\therefore 10+90+270+180+1800+7200+1920+16800+37800+75600+25200=166870通り$$

分割の中の数の最大公約数に着目する。

私のツイートでは分割の中の数の最大公約数が１のときは円形に並べる方法は1列に並べる総数を個数6で割れば求まります。最大公約数が1ではないものは(6)，(2-4)，(3-3)，(2-2-2)の4種類あります。その中で「6個で1つ」でないようなものを除きます。数が多くなると大変です。

コード化まとめ

分割数について

分割数はCombinatorics.jlのpartitionsを使いました。分割数だけでなく，どのような分割なのかも与えられます。便利です。

https://juliamath.github.io/Combinatorics.jl/dev/api/#Combinatorics.partitions-Tuple{Integer}

partitions

using Combinatorics
collect(partitions(6))

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

分割の様子から数字の選び方を定める

手計算で行うときに前半この総数が必要になるのですけど，初歩的に「配列の中から数字をそれぞれカウントする」ということが最初「どうしよう？？？」と思ったのですけれど，DataStructures.jlのcounterで実現できることがわかりました。

https://discourse.julialang.org/t/counting-number-of-occurences-in-an-array/31613

using DataStructures
c = counter([1,2,3,3,4,4,4,6])

Accumulator{Int64, Int64} with 5 entries:
4 => 3
6 => 1
2 => 1
3 => 2
1 => 1

keys(c)

KeySet for a Dict{Int64, Int64} with 5 entries. Keys:
4
6
2
3
1

values(c)

ValueIterator for a Dict{Int64, Int64} with 5 entries. Values:
3
1
1
2
1

関数の作成（iro,circperm3）

関数を作成してチェックです。

iro.jl

function iro(k,n)
   p = factorial(n)÷factorial(n-length(k))
   c = counter(k)
   for i in values(c)
       p = p÷factorial(i)
   end
   p
end

circperm3.j

function circperm3(m,n)
    k = collect(partitions(m))
    p = 0
    q = n^m
    for i in k
        p += enkan(i)*iro(i,n)
    end
    p
end

circperm3(6,10)

166870

バッチリです。

そのほかも確認です。

for i=1:10
   println("0〜9の10個の数字から重複を許して$i 個選び円形に並べる方法は",circperm3(i,10),"通りです。")
end

0〜9の10個の数字から重複を許して1 個選び円形に並べる方法は10通りです。
0〜9の10個の数字から重複を許して2 個選び円形に並べる方法は55通りです。
0〜9の10個の数字から重複を許して3 個選び円形に並べる方法は340通りです。
0〜9の10個の数字から重複を許して4 個選び円形に並べる方法は2530通りです。
0〜9の10個の数字から重複を許して5 個選び円形に並べる方法は20008通りです。
0〜9の10個の数字から重複を許して6 個選び円形に並べる方法は166870通りです。
0〜9の10個の数字から重複を許して7 個選び円形に並べる方法は1428580通りです。
0〜9の10個の数字から重複を許して8 個選び円形に並べる方法は12501280通りです。
0〜9の10個の数字から重複を許して9 個選び円形に並べる方法は111111340通りです。
0〜9の10個の数字から重複を許して10 個選び円形に並べる方法は1000010044通りです。