整数拆分与五边形数定理-优快云博客

本文链接：https://blog.youkuaiyun.com/iamhpp/article/details/114989481

本文介绍了整数拆分的概念，给出了拆分数的生成函数及其与欧拉函数的关系，并通过五边形数定理阐述了递推公式。讨论了HDU 4658问题，要求在限制条件下计算整数的拆分数，提出利用生成函数和递推的方法解决，复杂度为O(n^2)和O(nlogn)。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

定义

对于一个非负整数 $n$ ，它的一个拆分定义为若干个正整数 $a_1,a_2,...,a_k$ ，满足：

$a_1+a_2+...+a_k=n$
$a1≤a2≤...≤aka_1\le a_2\le ...\le a_k$

其中，第二个条件是为了避免重复计算，如 $6 = 1 + 2 + 3$ 和 $6 = 3 + 2 + 1$ 为同一方案。
令 $p_n$ 表示 $n$ 不同的拆分方案数。例如， $4 = 3 + 1 = 2 + 2 = 2 + 1 + 1 = 1 + 1 + 1 + 14 = 4 = 3 + 1 = 2 + 2 = 2 + 1 + 1 = 1 + 1 + 1 + 1$ ，于是 $p_4=5$ 。
对于 $p_n$ 的生成函数 $P (x)$ ，我们可以通过枚举 $1$ 的数量，枚举 $2$ 的数量…，枚举 $n$ 的数量得到，即 $P(x)=(∑i=0∞xi)(∑i=0∞(x2)i)...(∑i=0∞(xn)i)=11−x11−x2...11−xn=∏k=1∞11−xkP(x)=(\sum\limits_{i=0}^{\infty}x^i)(\sum\limits_{i=0}^{\infty}(x^2)^i)...(\sum\limits_{i=0}^{\infty}(x^n)^i)=\frac{1}{1-x}\frac{1}{1-x^2}...\frac{1}{1-x^n}=\prod\limits_{k=1}^{\infty}\frac{1}{1-x^k}$ 。

五边形数定理

定义欧拉函数 $ϕ(x)=∏k=1∞(1−xk)\phi(x)=\prod\limits_{k=1}^{\infty}(1-x^k)$ ，可以发现， $ϕ(x)P(x)=1\phi(x)P(x)=1$ 。
聪明的欧拉发现， $ϕ(x)=∑k=−∞∞(−1)kx3k2−k2=∑k=0∞(−1)kx3k2±k2\phi(x)=\sum\limits_{k=-\infty}^{\infty}(-1)^kx^{\frac{3k^2-k}{2}}=\sum\limits_{k=0}^{\infty}(-1)^kx^{\frac{3k^2\pm k}{2}}$ 。
容易证明，每个 $x$ 的系数只能为 $- 1, 0, 1$ ，且除了 $k = 0$ 外，不存在两个 $k_1,k_2$ ，使得 $3k12±k12=3k22±k22{\frac{3k_1^2\pm k_1}{2}}={\frac{3k_2^2\pm k_2}{2}}$ 。
于是， $ϕ(x)=1−x−x2+x5+x7+...\phi(x)=1-x-x^2+x^5+x^7+...$ 。那么考虑递推，根据 $ϕ(x)P(x)=1\phi(x)P(x)=1$ ，我们知道当 $n > 0$ 时， $p (n) - p (n - 1) - p (n - 2) + p (n - 5) + p (n - 7) . . . = 0$ 。
于是我们可以递推。注意到 $ϕ(x)\phi(x)$ 的项数是 $O(nO(\sqrt{n}$ ) 级别的，因此递推的复杂度是 $O(nn)O(n\sqrt{n})$ 的。
递推的代码：

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;

void debug_out(){
    cerr << endl;
}
template<typename Head, typename... Tail>
void debug_out(Head H, Tail... T){
    cerr << " " << to_string(H);
    debug_out(T...);
}
#ifdef local
#define debug(...) cerr<<"["<<#__VA_ARGS__<<"]:",debug_out(__VA_ARGS__)
#else
#define debug(...) 55
#endif

const int N = 1e5 + 5;
int f[N], g[N];
int F(int x){
	return (3 * x * x - x) / 2;
}
int main(){
	ios::sync_with_stdio(false);
	cin.tie(0), cout.tie(0);
	int n = 100000;
	const int mod = 1e9 + 7;
	f[0] = 1 % mod;
	for(int i = 1; i <= n; i++){//求 f[i]
		for(int j = 1; ; j++){//枚举 k
			int cof, x = F(-j), y = F(j);//次数分别为 x 和 y
			if(j % 2) cof = -1;//系数
			else cof = 1;
			if(x <= i) f[i] = (f[i] - cof * f[i - x]) % mod;
			if(y <= i) f[i] = (f[i] - cof * f[i - y]) % mod; 
			if(x > i || y > i) break;
		}
		f[i] = (f[i] + mod) % mod;
	}
	int T;
	cin >> T;
	while(T--){
		cin >> n;
		cout << f[n] << '\n';
	}
	return 0;
}

此外，还可以使用求任意模数多项式的逆来 $O (n l o g n)$ 求得 $P (x)$ 。
代码如下。

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;

void debug_out(){
    cerr << endl;
}
template<typename Head, typename... Tail>
void debug_out(Head H, Tail... T){
    cerr << " " << to_string(H);
    debug_out(T...);
}
#ifdef local
#define debug(...) cerr<<"["<<#__VA_ARGS__<<"]:",debug_out(__VA_ARGS__)
#else
#define debug(...) 55
#endif

namespace fft {
  const double pi = acos(-1.0);
  struct Complex {
    double r, i;
    Complex(double x = 0, double y = 0) : r(x), i(y) {}
    Complex operator+ (const Complex& b) const {
      return Complex(r + b.r, i + b.i);
    }
    Complex operator- (const Complex& b) const {
      return Complex(r - b.r, i - b.i);
    }
    Complex operator* (const Complex& b) const {
      return Complex(r * b.r - i * b.i, r * b.i + i * b.r);
    }
  };
  Complex conj(Complex a) {
    return Complex(a.r, -a.i);
  }

  int base = 1;
  vector<int> rev = { 0, 1 };
  vector<Complex> roots = { { 0, 0 }, { 1, 0 } };

  void ensure_base(int nbase) {
    if (nbase <= base) {
      return;
    }
    rev.resize(1 << nbase);
    for (int i = 0; i < (1 << nbase); i++) {
      rev[i] = (rev[i >> 1] >> 1) + ((i & 1) << (nbase - 1));
    }
    roots.resize(1 << nbase);
    while (base < nbase) {
      double angle = 2 * pi / (1 << (base + 1));
      for (int i = 1 << (base - 1); i < (1 << base); i++) {
        roots[i << 1] = roots[i];
        double angle_i = angle * (2 * i + 1 - (1 << base));
        roots[(i << 1) + 1] = Complex(cos(angle_i), sin(angle_i));
      }
      base++;
    }
  }
  void fft(vector<Complex> &a, int n = -1) {
    if (n == -1) {
      n = a.size();
    }
    assert((n & (n - 1)) == 0);
    int zeros = __builtin_ctz(n);
    ensure_base(zeros);
    int shift = base - zeros;
    for (int i = 0; i < n; i++) {
      if (i < (rev[i] >> shift)) {
        swap(a[i], a[rev[i] >> shift]);
      }
    }
    for (int k = 1; k < n; k <<= 1) {
      for (int i = 0; i < n; i += 2 * k) {
        for (int j = 0; j < k; j++) {
          Complex z = a[i + j + k] * roots[j + k];
          a[i + j + k] = a[i + j] - z;
          a[i + j] = a[i + j] + z;
        }
      }
    }
  }

  vector<Complex> fa, fb;
  vector<int> multiply(const vector<int> &a, const vector<int> &b) {
    int need = a.size() + b.size() - 1;
    int nbase = need > 1 ? 32 - __builtin_clz(need - 1) : 0;
    ensure_base(nbase);
    int sz = 1 << nbase;
    if (sz > (int) fa.size()) {
      fa.resize(sz);
    }
    for (int i = 0; i < sz; i++) {
      int x = (i < (int) a.size() ? a[i] : 0);
      int y = (i < (int) b.size() ? b[i] : 0);
      fa[i] = Complex(x, y);
    }
    fft(fa, sz);
    Complex r(0, -0.25 / sz);
    for (int i = 0; i <= (sz >> 1); i++) {
      int j = (sz - i) & (sz - 1);
      Complex z = (fa[j] * fa[j] - conj(fa[i] * fa[i])) * r;
      if (i != j) {
        fa[j] = (fa[i] * fa[i] - conj(fa[j] * fa[j])) * r;
      }
      fa[i] = z;
    }
    fft(fa, sz);
    vector<int> res(need);
    for (int i = 0; i < need; i++) {
      res[i] = fa[i].r + 0.5; // watch out that fa[i].r < 0
    }
    return res;
  }
  vector<int> multiply_mod(const vector<int> &a, const vector<int> &b, int m, int eq = 0) {
    int need = a.size() + b.size() - 1;
    int nbase = need > 1 ? 32 - __builtin_clz(need - 1) : 0;
    ensure_base(nbase);
    int sz = 1 << nbase;
    if (sz > (int) fa.size()) {
      fa.resize(sz);
    }
    for (int i = 0; i < (int) a.size(); i++) {
      int x = (a[i] % m + m) % m;
      fa[i] = Complex(x & ((1 << 15) - 1), x >> 15);
    }
    fill(fa.begin() + a.size(), fa.begin() + sz, Complex {0, 0});
    fft(fa, sz);
    if (sz > (int) fb.size()) {
      fb.resize(sz);
    }
    if (eq) {
      copy(fa.begin(), fa.begin() + sz, fb.begin());
    } else {
      for (int i = 0; i < (int) b.size(); i++) {
        int x = (b[i] % m + m) % m;
        fb[i] = Complex(x & ((1 << 15) - 1), x >> 15);
      }
      fill(fb.begin() + b.size(), fb.begin() + sz, Complex {0, 0});
      fft(fb, sz);
    }
    double ratio = 0.25 / sz;
    Complex r2(0, -1), r3(ratio, 0), r4(0, -ratio), r5(0, 1);
    for (int i = 0; i <= (sz >> 1); i++) {
      int j = (sz - i) & (sz - 1);
      Complex a1 = (fa[i] + conj(fa[j]));
      Complex a2 = (fa[i] - conj(fa[j])) * r2;
      Complex b1 = (fb[i] + conj(fb[j])) * r3;
      Complex b2 = (fb[i] - conj(fb[j])) * r4;
      if (i != j) {
        Complex c1 = (fa[j] + conj(fa[i]));
        Complex c2 = (fa[j] - conj(fa[i])) * r2;
        Complex d1 = (fb[j] + conj(fb[i])) * r3;
        Complex d2 = (fb[j] - conj(fb[i])) * r4;
        fa[i] = c1 * d1 + c2 * d2 * r5;
        fb[i] = c1 * d2 + c2 * d1;
      }
      fa[j] = a1 * b1 + a2 * b2 * r5;
      fb[j] = a1 * b2 + a2 * b1;
    }
    fft(fa, sz);
    fft(fb, sz);
    vector<int> res(need);
    for (int i = 0; i < need; i++) {
      long long aa = fa[i].r + 0.5;
      long long bb = fb[i].r + 0.5;
      long long cc = fa[i].i + 0.5;
      res[i] = (aa + ((bb % m) << 15) + ((cc % m) << 30)) % m;
    }
    return res;
  }
  vector<int> square_mod(const vector<int> &a, int m) {
    return multiply_mod(a, a, m, 1);
  }
}using namespace fft;

namespace poly {
  // reference exgcd and fft::multiply_mod
  #define Poly vector<int>
  const int MOD = 1e9 + 7;
  void exgcd(int a, int b, int &x, int &y){
  	if(!b){
  		x = 1, y = 0;
  		return;
	  }
	  exgcd(b, a % b, y, x);
	  y -= a / b * x;
  }
  vector<int> inv(const vector<int> &a) {
    if(a.size() == 1) {
      int x, y;
      exgcd(a[0], MOD, x, y);
      const int inv = (x % MOD + MOD) % MOD;
      return vector<int>(1, inv < 0 ? inv + MOD : inv);
    }
    const int na = a.size(), nb = (na + 1) >> 1;
    vector<int> b(a.begin(), a.begin() + nb);
    b = inv(b);
    vector<int> c = fft::multiply_mod(b, b, MOD, 1);
    c.resize(na);
    c = fft::multiply_mod(a, c, MOD);
    b.resize(na), c.resize(na);
    for(int i = 0; i < na; i++) {
      c[i] = (((2ll * b[i] - c[i]) % MOD) + MOD) % MOD;
    }
    return c;
  }

  // make mod function by yourself
  vector<int> divide(const vector<int> &a, const vector<int> &b) {
    const int n = a.size(), m = b.size();
    vector<int> A(a), B(b);
    reverse(A.begin(), A.end()), reverse(B.begin(), B.end());
    A.resize(n - m + 1), B.resize(n - m + 1);
    B = inv(B);
    vector<int> C = fft::multiply_mod(A, B, MOD);
    C.resize(n - m + 1), reverse(C.begin(), C.end());
    return C;
  }

  vector<int> Mod(const vector<int> &a, const vector<int> &b) {
    vector<int> c = divide(a, b);
    c = fft::multiply_mod(b, c, MOD);
    c.resize(b.size() - 1);
    for(int i = 0; i < int(c.size()); i++) c[i] = (a[i] + MOD - c[i]) % MOD;
    return c;
  }
}using namespace poly;
int F(int x){
	return (3 * x * x - x) / 2;
}
int main(){
	ios::sync_with_stdio(false);
	cin.tie(0), cout.tie(0);
	int n = 100000;
	Poly f(n + 1);
	f[0] = 1;
	for(int i = 1; ; i++){
		int x = F(i), y = F(-i), c = i % 2? MOD - 1: 1;
		if(x <= n) f[x] = c;
		if(y <= n) f[y] = c;
		if(x >= n || y >= n) break;
	}
	Poly g = inv(f);
	int T;
	cin >> T;
	while(T--){
		cin >> n;
		cout << g[n] << '\n';
	}
	return 0;
}

[HDU 4658] Integer Partition

有 $T$ 组数据，每组数据给出 $n, k$ ，求 $n$ 的划分方案数，满足 $a_i< k$ 。 $n,k,T≤105n,k,T\le 10^5$ 。

容易得出，生成函数为 $F(x)=∏i=0∞1−xik1−xi=ϕ(xk)P(x)F(x)=\prod\limits_{i=0}^{\infty}\frac{1-x^{ik}}{1-x^i}=\phi(x^k)P(x)$ 。
其中， $ϕ(xk)=(1−xk)(1−x2k)...(1−xnk)=1−xk−x2k+x5k+x7k+...=∑i=0∞(−1)ix(3i2±i)k2\phi(x^k)=(1-x^k)(1-x^{2k})...(1-x^{nk})=1-x^k-x^{2k}+x^{5k}+x^{7k}+...=\sum\limits_{i=0}^{\infty}(-1)^ix^{\frac{(3i^2\pm i)k}{2}}$ 。
然后一开始求出 $P (x)$ ，每次询问求 $ϕ(xk)\phi(x^k)$ 即可。复杂度为 $O(nn)O(n\sqrt{n})$ 。
代码如下。

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;

void debug_out(){
    cerr << endl;
}
template<typename Head, typename... Tail>
void debug_out(Head H, Tail... T){
    cerr << " " << to_string(H);
    debug_out(T...);
}
#ifdef local
#define debug(...) cerr<<"["<<#__VA_ARGS__<<"]:",debug_out(__VA_ARGS__)
#else
#define debug(...) 55
#endif

const int N = 1e5 + 5, mod = 1e9 + 7;
int F(int x){
	return (3 * x * x - x) / 2;
}
int f[N], g[N];
int main(){
	ios::sync_with_stdio(false);
	cin.tie(0), cout.tie(0);
	int n = 100000;
	f[0] = 1;
	for(int i = 1; i <= n; i++){
		for(int j = 1; ; j++){
			int cof, x = F(-j), y = F(j);
			if(j % 2) cof = -1;
			else cof = 1;
			if(x <= i) f[i] = (f[i] - cof * f[i - x]) % mod;
			if(y <= i) f[i] = (f[i] - cof * f[i - y]) % mod; 
			if(x > i || y > i) break;
		}
		f[i] = (f[i] + mod) % mod;
	}
	int T;
	cin >> T;
	while(T--){
		int k, ans = 0;
		cin >> n >> k;
		ans = f[n];
		for(int j = 1; ; j++){
			int cof, x = F(-j) * k, y = F(j) * k;
			if(j % 2) cof = -1;
			else cof = 1;
			if(x <= n) ans = (ans + cof * f[n - x]) % mod;
			if(y <= n) ans = (ans + cof * f[n - y]) % mod; 
			if(x > n || y > n) break;
		}
		cout << (ans + mod) % mod << '\n';
	}
	return 0;
}