原创:HAM
10人分4组,有几种分法,当然不考虑组编号了,比如10人分1组,就只有1种分法,10人分10组也只有一种分法。
假如4人分2组,分法为:
12 , 34
13 , 24
14 , 23
123 , 4
124 , 3
134 , 2
234 , 1
共7种,其余都是重复的。
现在我们考虑10人分4组,像上面那样穷举是肯定不行的,必须找出一种通用的算法来解决。
分组意味着每组至少有一个人,所以我们这样考虑:
第一组选7个人,其余组各1人,则分法有:
Binomial[10,7] * Binomial[3,1] * Binomial[2,1] * 1 / 3!
这里Binomial[n,m] = n!/(m!*(n-m)!)
第一组选6个人,第二组1人,第三组1人,第四组2人:
Binomial[10,6] * Binomial[4,1] * Binomial[3,1] * 1 / 2!
后面除以3!,2!是因为这几组人数相同,会出现相同的分组,不同的排列。
第一、二组选2人,第三、四组选3人:
Binomial[10,2] * Binomial[8,2] * Binomial[6,3] * 1 / 2! / 2!
由此规律,可以看出里面蕴含着一种递归,也就是说,n人分m组(n>=m),如果第一组分a个人(1<=a<=n-m+1),那么第一组分a个人的次数为(暂时不考虑相同分组的组组之间的排列):
T[n,1]::第一组a人 = 1 , n >= 1
T[n,m]::第一组a人 = Binomial[n,a]*T[n-a,m-1] , n >= m >= 2
那么n人分m组的分法就为:
T[n,1] = 1 , n >= 1
T[n,m] = Sum[Binomial[n,i]*T[n-i,m-1], {i,1,n-m+1}] , n >= m >= 2
enum[n,m] = T[n,m]/m! ,这就是最终的分组,/m!是因为除去m组相同分组的m!种排列
下面是我以此模型写出的递归Mathematica程序:
T[n_,m_]:=
Block[{res = 0, i},
If[Mod[n, 1] != 0 || Mod[m, 1] != 0, Return[0]]; (*判断不合法输入*)
If[n < m || n < 1 || m < 1, Return[0]]; (*判断不合法输入*)
If[m != 1,
For [i = n-m+1, i >= 1, i--,
res += Binomial[n,i]*T[n-i,m-1];
];
Return[res]
,
Return[1];
]
]
enum[n_,m_]:=T[n,m]/m!
输入:
enum[10, 4]
输出:
34105
不难看出该递归模型内含有大量的重复子问题,时间效率也是指数次O((n-m+1)^(m-1))的,那么现在的任务是避免重复计算子问题,也就是使用迭代计算,下面就是这样的Mathematica程序,至于为什么这么写,你可以去研究上面的递归程序,从而发现规律:
enum[n_, m_]:=
Block[{i, j, k, l, t, sta},
If[Mod[n, 1] != 0 || Mod[m, 1] != 0, Return[0]]; (*判断不合法输入*)
If[n < m || n < 1 || m < 1, Return[0]]; (*判断不合法输入*)
If[m == 1, Return[1]];
j = 2;
i = n - m + 2;
sta = Table[1,{i-1}];
While[i != n,
For[k = i, k > j-1, k--,
t = 0;
For[l = j-1, l < k, l++,
t += sta[[l-j+2]]*Binomial[k,l];
];
sta[[k-j+1]] = t;
];
i++;
j++;
];
t = 0;
For[l = m-1, l < n, l++,
t += sta[[l-m+2]]*Binomial[n,l];
];
Return[t / m!];
]
输入:
enum[10, 4]
输出:
34105
以上两个程序都是以递归模型来考虑的,只不过上面程序是从顶向下,下面的程序是从底向上,没有多少本质区别。
现在我们以枚举表面不同分组的方式来解决。
为了不使枚举不同分组时出现类似于(2,3,1,4)和(1,3,4,2)这样两种一样的分组情形,我们在枚举时可以采用从小到大排列方式来枚举,例如10人分4组:
1,1,1,7
1,1,2,6
1,1,3,5
1,1,4,4
1,2,2,5
1,2,3,4
1,3,3,3
2,2,2,4
2,2,3,3
下面的这个程序技巧性非常强,你也没有必要去研究这段代码,总之就是通过枚举各种不同的分组来计算,当然因为是枚举的原因,速度也不会快,只会慢。
enum[n_, m_]:=
Block[{i, j, q, p, r, s, t, sta, res},
If[Mod[n, 1] != 0 || Mod[m, 1] != 0, Return[0]]; (*判断不合法输入*)
If[n < m || n < 1 || m < 1, Return[0]]; (*判断不合法输入*)
If[m == 1, Return[1]];
j = 1;
res = 0;
sta = Table[0,{m}];
While[j > 0,
sta[[j]]++;
For[i = j+1, i < m, i++,
sta[[i]] = sta[[j]];
];
t = sta[[1]];
r = 1;
p = 1;
For[i = 2, i < j, i++,
If[sta[[i]] == t,
r++
,
p *= r!;
t = sta[[i]];
r = 1;
];
];
p *= r!*(m-j)!;
s = 0;
For[i = 1, i < m, i++,
s += sta[[i]];
];
If[sta[[ m-1]] <= n-s,
sta[[ m]] = n-s;
t = Binomial[n,sta[[1]]];
q = sta[[1]];
For[i = 2, i < m, i++,
t *= Binomial[n-q,sta[[i]]];
q += sta[[i]];
];
t /= p;
If[sta[[ m-1]] == n-s,
t /= m-j+1;
];
res += t
,
j--;
Continue[];
];
p = p/(m-j);
j = m-1;
t = Binomial[n,sta[[1]]];
q = sta[[1]];
For[i = 2, i < j, i++,
t *= Binomial[n-q,sta[[i]]];
q += sta[[i]];
];
s++;
sta[[j]]++;
While[sta[[j]] <= n-s,
sta[[ m]] = n-s;
r = t * Binomial[n-q,sta[[j]]] / p;
If[sta[[j]] == n-s,
r /= 2;
];
res += r;
s++;
sta[[j]]++;
];
j--;
];
Return[res]
]
输入:
enum[10, 4]
输出:
34105
现在我们继续讨论,有上个递归式中,不难发现,就算是优化了的迭代程序,其中也大量的使用Binomial[]函数,这样其中隐含的进行了大量重复计算,换一种角度考虑,比如现在n人分m组,假如先选出 n-1 人进行分组,再把最后一个人塞进去,如果先选出的 n-1 个人分 m 组的话,分发就有:
m * enum[n-1,m] , n > m > 1
enum[n-1,m]是n-1人分m组的分发,*m是因为最后的的那个人可以有m种方法选择到哪组;
因为每组至少一个人,所以还可以 n-1 人先分到 m-1 组中,而最后一个人就没得选了,只能塞到那个空组里,分发就有:
enum[n-1,m-1] , n > m > 1
总结:
enum[n,1] = 1 , n > 1 ,第一类初始状态
enum[n,m] = 1 , n = m >= 1 ,第二类初始状态
enum[n,m] = enum[n-1,m-1] + m*enum[n-1,m] , n > m > 1
以此递归模型写出Mathematica程序:
enum[n_,m_]:=
Block[{},
If[Mod[n, 1] != 0 || Mod[m, 1] != 0, Return[0]]; (*判断不合法输入*)
If[n < m || n < 1 || m < 1, Return[0]]; (*判断不合法输入*)
If[m != 1 && n > m, (*判断是否是初始状态,如果是就返回1*)
Return[enum[n-1,m-1]+m*enum[n-1,m]] (*根据子问题的解得出本次解*)
,
Return[1]; (*返回所有的初始值,也就是递归的叶子*)
];
]
输入:
enum[10, 4]
输出:
34105
优化使用迭代方法,避免重复计算子问题,至于为什么这么写,你可以去研究上面的递归程序,从而发现规律:
enum[n_,m_]:=
Block[{sta, i, j, c},
If[Mod[n, 1] != 0 || Mod[m, 1] != 0, Return[0]]; (*判断不合法输入*)
If[n < m || n < 1 || m < 1, Return[0]]; (*判断不合法输入*)
c = n-m;
If[c == 0, Return[1]]; (*当n = m时,下列代码无法处理,但结果必定为1*)
sta = Table[1,{c}];
For[i = 2, i <= m, i++,
sta[[1]] += i;
For[j = 2, j <= c, j++,
sta[[j]] += i*sta[[j-1]];
];
];
Return[sta[[c]]];
]
输入:
enum[10, 4]
输出:
34105