1到11的任意一个排列都可以分解为一些不相交的循环子排列（以下称为状态），11共有56种不同的分割方式，因此这$11!$种状态可以分割为56组，每组分解为循环子排列后的结构相同。
特别地，目标状态由11个长度为1的循环子排列组成，该组中只有这一个状态，为叙述方便，规定其为第一组。
选择的3个数有种可能：分布在同一个循环子排列中，分布在三个不同的循环子排列中，以及有两个数在同一个循环子排列，一个数在另一个子排列中。
交换两个位于同一循环子排列的数所得的新状态与原状态在同一组中，而交换两个不同循环子排列的数将会合并这两个子排列。
因此先计算出每组的大小，然后枚举所有三个数的选法，对每一种选法检验6个排列是如何在组间转化的，如此可以计算出组间的状态转移矩阵A。
计算$I+A+A^2+\ldots=(I-A)^{-1}$第一列的总和（第一列可以理解为从第一组有序状态出发，平均分别需要多少次能转化到其他组）。
最终结果是$48271207$。
注：以下代码改编在官方论坛，且在实现中进一步略去了A的第一行和第一列。此外，因用到Numpy库求逆，故为Python 3代码。
注2：本题涉及的随机过程知识，一般在数学系本科高年级及以上的相应课程中开设。
from numpy import matrix
from math import factorial,comb
from itertools import combinations
target=11
partitions=((1,1,1,1,1,1,1,1,1,2),(1,1,1,1,1,1,1,1,3),(1,1,1,1,1,1,1,2,2),(1,1,1,1,1,1,1,4),(1,1,1,1,1,1,2,3),(1,1,1,1,1,1,5),(1,1,1,1,1,2,2,2),(1,1,1,1,1,2,4),(1,1,1,1,1,3,3),(1,1,1,1,1,6),(1,1,1,1,2,2,3),(1,1,1,1,2,5),(1,1,1,1,3,4),(1,1,1,1,7),(1,1,1,2,2,2,2),(1,1,1,2,2,4),(1,1,1,2,3,3),(1,1,1,2,6),(1,1,1,3,5),(1,1,1,4,4),(1,1,1,8),(1,1,2,2,2,3),(1,1,2,2,5),(1,1,2,3,4),(1,1,2,7),(1,1,3,3,3),(1,1,3,6),(1,1,4,5),(1,1,9),(1,2,2,2,2,2),(1,2,2,2,4),(1,2,2,3,3),(1,2,2,6),(1,2,3,5),(1,2,4,4),(1,2,8),(1,3,3,4),(1,3,7),(1,4,6),(1,5,5),(1,10),(2,2,2,2,3),(2,2,2,5),(2,2,3,4),(2,2,7),(2,3,3,3),(2,3,6),(2,4,5),(2,9),(3,3,5),(3,4,4),(3,8),(4,7),(5,6),(11,))
m=len(partitions)
a=[]
v=[1]*m
for k,p in enumerate(partitions):
p=partitions[k]
x=target
c=[0]*(target+1)
for i in p:
c[i]+=1
v[k]*=comb(x,i)
x-=i
for i in range(1,target+1):
v[k]*=factorial(i-1)**c[i]
v[k]//=factorial(c[i])
row=[0]*m
row[k]=target*(target-1)*(target-2)
x=0
perm=[]
for i in p:
for j in range(x+1,x+i):
perm.append(j)
perm.append(x)
x+=i
for i0,i1,i2 in combinations(range(target),3):
for j0,j1,j2 in ((i0,i1,i2),(i0,i2,i1),(i1,i0,i2),(i1,i2,i0),(i2,i0,i1),(i2,i1,i0)):
swapped=perm[:]
swapped[j0],swapped[j1],swapped[j2]=perm[i0],perm[i1],perm[i2]
q=[]
untouched=[True]*target
while True:
i=0
while i < target and not untouched[i]:
i+=1
if i==target:
break
x=0
while untouched[i]:
untouched[i]=False
i=swapped[i]
x+=1
q.append(x)
try:
row[partitions.index(tuple(sorted(q)))]-=1
except:
pass
a.append(row)
print(round((matrix(v)*matrix(a).I*matrix([[1]]*m))[0,0]/factorial(target-3)))
|
