引言
在组合数学中，如何将n个不同的元素划分成k个非空子集？又如何计算n个元素的k阶排列（即恰好使用k个循环的排列）数目？这些问题的答案由斯特林数（Stirling Numbers）给出。作为连接幂函数与下降阶乘、排列与组合的桥梁，斯特林数在计数组合学、概率论和算法分析中扮演着核心角色。本文将系统讲解第一类与第二类斯特林数的定义、递推关系、显式公式，并用Java实现高效计算框架。
数学定义
第二类斯特林数 S(n, k)
将n个不同元素划分为k个非空、无序子集的方法数。记作 S(n, k) 或 {n \brace k}。
例如，S(4, 2) = 7：
– {1,2,3} {4} | {1,2,4} {3} | {1,3,4} {2} | {2,3,4} {1}
– {1,2} {3,4} | {1,3} {2,4} | {1,4} {2,3}
递推公式：S(n, k) = S(n-1, k-1) + k * S(n-1, k)
– S(n-1, k-1)：第n个元素单独成一个子集
– k * S(n-1, k)：第n个元素放入已有的k个子集之一
边界条件：S(n, 0) = 0 (n>0)，S(0, 0) = 1，S(n, n) = 1
第一类斯特林数 s(n, k)
n个不同元素恰好形成k个循环排列（轮换）的排列数。记作 s(n, k) 或 [n \brace k]。
递推公式：s(n, k) = s(n-1, k-1) + (n-1) * s(n-1, k)
– s(n-1, k-1)：第n个元素单独成一个循环
– (n-1) * s(n-1, k)：第n个元素插入已有循环中的某个位置（n-1个位置）
边界条件：s(n, 0) = 0 (n>0)，s(0, 0) = 1，s(n, n) = 1
核心性质与显式公式
第二类斯特林数的显式公式
$$S(n, k) = \frac{1}{k!} \sum_{i=0}^{k} (-1)^i \binom{k}{i} (k-i)^n$$
这可以通过容斥原理推导得出。
斯特林数与幂函数的关系
$$x^n = \sum_{k=0}^{n} S(n, k) \cdot x^{\underline{k}}$$
其中 x^{\underline{k}} = x(x-1)(x-2)...(x-k+1) 为下降阶乘。
这个公式将幂函数表示为下降阶乘的线性组合，在求和变换和有限差分中有重要应用。
Java实现
import java.util.*;
/**
* 斯特林数计算：第一类与第二类斯特林数
* 支持动态规划递推和容斥原理显式公式两种方法
*/
public class StirlingNumbers {
private static final long MOD = 1_000_000_007L;
/**
* 第二类斯特林数：动态规划递推
* S(n, k) = S(n-1, k-1) + k * S(n-1, k)
* 时间复杂度: O(n*k), 空间复杂度: O(k)
*/
public static long stirling2nd(int n, int k) {
if (n < 0 || k < 0 || k > n) return 0;
if (n == 0 && k == 0) return 1;
if (k == 0) return 0;
// 滚动数组优化空间
long[] dp = new long[k + 1];
dp[0] = 1; // S(0, 0) = 1
for (int i = 1; i <= n; i++) {
// 倒序更新避免覆盖
long[] newDp = new long[k + 1];
for (int j = 1; j <= Math.min(i, k); j++) {
newDp[j] = (dp[j - 1] + j * dp[j]) % MOD;
}
dp = newDp;
}
return dp[k];
}
/**
* 第二类斯特林数：显式公式（容斥原理）
* S(n, k) = (1/k!) * sum_{i=0}^{k} (-1)^i * C(k,i) * (k-i)^n
* 适用于k较小的情况
*/
public static long stirling2ndExplicit(int n, int k) {
if (n < 0 || k < 0 || k > n) return 0;
if (n == 0 && k == 0) return 1;
if (k == 0) return 0;
long result = 0;
for (int i = 0; i <= k; i++) {
long sign = (i % 2 == 0) ? 1 : -1;
long comb = combination(k, i);
long power = modPow(k - i, n, MOD);
long term = sign * comb * power % MOD;
result = (result + term + MOD) % MOD;
}
// 除以k!，即乘以k!的模逆元
long invFactK = modInverse(factorial(k), MOD);
return result * invFactK % MOD;
}
/**
* 第一类斯特林数（无符号）：动态规划递推
* c(n, k) = c(n-1, k-1) + (n-1) * c(n-1, k)
* 时间复杂度: O(n*k), 空间复杂度: O(k)
*/
public static long stirling1stUnsigned(int n, int k) {
if (n < 0 || k < 0 || k > n) return 0;
if (n == 0 && k == 0) return 1;
if (k == 0) return 0;
long[] dp = new long[k + 1];
dp[0] = 1;
for (int i = 1; i <= n; i++) {
long[] newDp = new long[k + 1];
for (int j = 1; j <= Math.min(i, k); j++) {
newDp[j] = (dp[j - 1] + (i - 1) * dp[j]) % MOD;
}
dp = newDp;
}
return dp[k];
}
/**
* 计算完整斯特林数三角形（返回前n行的所有值）
*/
public static long[][] stirling2ndTriangle(int n) {
long[][] S = new long[n + 1][n + 1];
S[0][0] = 1;
for (int i = 1; i <= n; i++) {
for (int j = 1; j <= i; j++) {
S[i][j] = (S[i - 1][j - 1] + j * S[i - 1][j]) % MOD;
}
}
return S;
}
/**
* 贝尔数：集合划分的总数，B(n) = sum_{k=0}^{n} S(n, k)
*/
public static long bellNumber(int n) {
long[][] S = stirling2ndTriangle(n);
long sum = 0;
for (int k = 0; k <= n; k++) {
sum = (sum + S[n][k]) % MOD;
}
return sum;
}
/**
* 将幂函数x^n表示为下降阶乘的线性组合（使用斯特林数）
* x^n = sum_{k=0}^{n} S(n,k) * x(x-1)...(x-k+1)
*/
public static long powerToFallingFactorial(long x, int n) {
long[][] S = stirling2ndTriangle(n);
long result = 0;
long falling = 1; // x^0 = 1
for (int k = 0; k <= n; k++) {
if (k > 0) {
falling = falling * (x - k + 1) % MOD;
}
result = (result + S[n][k] * falling) % MOD;
}
return result;
}
// ============ 辅助函数 ============
private static long combination(int n, int k) {
if (k < 0 || k > n) return 0;
long num = 1, den = 1;
for (int i = 0; i < k; i++) {
num = num * (n - i) % MOD;
den = den * (i + 1) % MOD;
}
return num * modInverse(den, MOD) % MOD;
}
private static long factorial(int n) {
long result = 1;
for (int i = 2; i <= n; i++) {
result = result * i % MOD;
}
return result;
}
private static long modPow(long base, long exp, long mod) {
long result = 1 % mod;
base %= mod;
while (exp > 0) {
if ((exp & 1) == 1) result = result * base % mod;
base = base * base % mod;
exp >>= 1;
}
return result;
}
private static long modInverse(long a, long mod) {
return modPow(a, mod - 2, mod); // 费马小定理，要求mod为质数
}
public static void main(String[] args) {
// 测试1：第二类斯特林数三角形
System.out.println("=== 第二类斯特林数 S(n,k) ===");
long[][] S = stirling2ndTriangle(8);
for (int i = 0; i <= 8; i++) {
for (int j = 0; j <= i; j++) {
System.out.printf("%6d ", S[i][j]);
}
System.out.println();
}
// 测试2：验证显式公式与递推的一致性
System.out.println("\n=== 验证显式公式 ===");
for (int n = 5; n <= 10; n++) {
for (int k = 1; k <= n; k++) {
long dp = stirling2nd(n, k);
long exp = stirling2ndExplicit(n, k);
if (dp != exp) {
System.out.println("不一致! n=" + n + " k=" + k);
}
}
}
System.out.println("n=5..10全部一致");
// 测试3：第一类斯特林数
System.out.println("\n=== 第一类斯特林数 c(n,k) ===");
for (int n = 1; n <= 6; n++) {
for (int k = 1; k <= n; k++) {
System.out.printf("%6d ", stirling1stUnsigned(n, k));
}
System.out.println();
}
// 测试4：贝尔数
System.out.println("\n=== 贝尔数 ===");
for (int n = 0; n <= 8; n++) {
System.out.println("B(" + n + ") = " + bellNumber(n));
}
// 测试5：幂转下降阶乘验证
System.out.println("\n=== 幂转下降阶乘验证 ===");
long x = 5;
for (int n = 1; n <= 6; n++) {
long direct = modPow(x, n, MOD);
long viaStirling = powerToFallingFactorial(x, n);
System.out.println(x + "^" + n + " = " + direct + " (直接) vs " + viaStirling + " (斯特林)");
}
}
}
复杂度分析
- 动态规划递推：时间
O(n*k)，空间O(k)（滚动数组优化后）。 - 显式公式（容斥）：时间
O(k log n)（快速幂），空间O(1)，适合k较小（如k ≤ 10^5）而n极大的场景。 - 完整三角形：时间
O(n²)，空间O(n²)，适用于需要同时查询多个值的情况。
应用场景
- 集合划分问题：会议室分配、任务分组的方案计数
- 排列循环结构：分析置换群中的共轭类大小
- 多项式插值：利用斯特林数在幂基与下降阶乘基之间转换
- 算法复杂度分析：分析递归算法调用树的结构数
总结
斯特林数作为组合数学的基石，优雅地刻画了集合划分与循环排列的本质结构。本文从递推定义出发，推导了显式公式，实现了动态规划与容斥原理两种计算方法，并展示了贝尔数、幂-下降阶乘转换等扩展应用。掌握斯特林数，将为解决复杂的组合计数问题提供强大的理论工具。
