【NOI2019模拟2019.6.27】B (生成函数+整数划分dp|多项式exp)

时间:2023-03-09 05:43:01
【NOI2019模拟2019.6.27】B (生成函数+整数划分dp|多项式exp)

Description:


【NOI2019模拟2019.6.27】B (生成函数+整数划分dp|多项式exp)

\(1<=n,k<=1e5,mod~1e9+7\)

题解:


考虑最经典的排列dp,每次插入第\(i\)大的数,那么可以增加的逆序对个数是\(0-i-1\)。

不难得到生成函数:

\(Ans=\prod_{i=0}^{n-1}(\sum_{j=0}^ix^j)[x^k]\)

\(=\prod_{i=1}^{n}{1-x^i\over 1-x}[x^k]\)

分母是一个经典的生成函数:

\({1\over 1-x}^n=(\sum_{i>=0}x^i)^n=\sum_{i>=0}C_{i+n-1}^{n-1}\)

那么问题就变为了求:

\(\prod_{i=1}^{n}{1-x^i}\)的前k项。

考虑利用整数划分dp,相当于把k划分成若干不同且<=n的数和,系数是\((-1)^{数的个数}\)。

不难得出dp:

设\(f[i][j]\)表示已经划分了i个数,和为j的所有方案系数和。

转移\(f[i][j]=f[i][j-i]-f[i-1][j-i]+f[i-1][j-(n+1)]\)

由于\(i<=\sqrt {2k}\),所以复杂度是\(O(k\sqrt k)\)。

另一种多项式exp的做法:

不妨对这个式子进行ln最后再exp回去。

我们知道有:

\(ln(1+x)\)

$=\int ~ln(1+x)' $

\(=\int~{1\over 1+x}\)

\(=\int ~ \sum_{i>=0}(-1)^ix^i\)

\(=\sum_{i>=1}{(-1)^{i-1}x^i\over i}\)

\(Ans=exp(\sum_{i=1}^n(ln(1-x^i)-ln(1-x)))[x^k]\)

\(ln(1-x^i)=-\sum_{j>=1}{x^{ij} \over j}\)

所以暴力展开只有调和级数个有用项。

\(ln(1-x)\)同理暴力展开后乘上n即可。

复杂度\(O(n~log~n)\),但是要写MTT,所以跑得巨慢,又难写。

Code:

#include<bits/stdc++.h>
#define fo(i, x, y) for(int i = x, B = y; i <= B; i ++)
#define ff(i, x, y) for(int i = x, B = y; i < B; i ++)
#define fd(i, x, y) for(int i = x, B = y; i >= B; i --)
#define ll long long
#define pp printf
#define hh pp("\n")
using namespace std; const int mo = 1e9 + 7; ll ksm(ll x, ll y) {
ll s = 1;
for(; y; y /= 2, x = x * x % mo)
if(y & 1) s = s * x % mo;
return s;
} const int N = 1e5 + 5; int n, k, m;
ll fac[N * 2], nf[N * 2];
ll f[450][N];
ll g[N];
void calc(int n) {
fac[0] = 1; fo(i, 1, n) fac[i] = fac[i - 1] * i % mo;
nf[n] = ksm(fac[n], mo - 2); fd(i, n, 1) nf[i - 1] = nf[i] * i % mo;
}
ll C(int n, int m) {
return fac[n] * nf[n - m] % mo * nf[m] % mo;
} int main() {
freopen("b.in", "r", stdin);
freopen("b.out", "w", stdout);
calc(200000);
scanf("%d %d", &n, &k);
m = sqrt(2 * k);
f[0][0] = 1;
fo(i, 1, m) {
fo(j, i, k) {
f[i][j] = f[i][j - i] - f[i - 1][j - i];
if(j >= n + 1) f[i][j] += f[i - 1][j - (n + 1)];
f[i][j] %= mo;
}
}
ll ans = 0;
fo(i, 0, k) {
fo(j, 0, m) g[i] += f[j][i];
g[i] %= mo;
ans += g[i] * fac[n - 1 + (k - i)] % mo * nf[k - i] % mo;
}
ans = (ans % mo * nf[n - 1] % mo + mo) % mo;
pp("%lld\n", ans);
}
#include<bits/stdc++.h>
#define fo(i, x, y) for(int i = x, B = y; i <= B; i ++)
#define ff(i, x, y) for(int i = x, B = y; i < B; i ++)
#define fd(i, x, y) for(int i = x, B = y; i >= B; i --)
#define ll long long
#define pp printf
#define hh pp("\n")
#define db double
using namespace std; const int mo = 1e9 + 7; typedef vector<ll> V;
#define si size()
#define pb push_back namespace ntt {
const db pi = acos(-1);
struct P {
db x, y;
P(db _x = 0, db _y = 0) { x = _x, y = _y;}
};
P operator + (P a, P b) { return P(a.x + b.x, a.y + b.y);}
P operator - (P a, P b) { return P(a.x - b.x, a.y - b.y);}
P operator * (P a, P b) { return P(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x);}
const int nm = 1 << 18;
int r[nm];
P w[nm], c0[nm], c1[nm], c2[nm], c3[nm]; void dft(P *a, int n) {
ff(i, 0, n) {
r[i] = r[i / 2] / 2 + (i & 1) * (n / 2);
if(i < r[i]) swap(a[i], a[r[i]]);
} P b;
for(int i = 1; i < n; i *= 2) for(int j = 0; j < n; j += 2 * i)
ff(k, 0, i) b = a[i + j + k] * w[i + k], a[i + j + k] = a[j + k] - b, a[j + k] = a[j + k] + b;
}
void rev(P *a, int n) {
reverse(a + 1, a + n);
ff(i, 0, n) a[i].x /= n, a[i].y /= n;
}
P conj(P a) { return P(a.x, -a.y);}
void fft(V &a, V b) {
#define qz(x) ((ll) round(x))
int n = a.si;
// ff(i, 0, n) c0[i] = P(a[i], 0), c1[i] = P(b[i], 0);
// dft(c0, n); dft(c1, n);
// ff(i, 0, n) c0[i] = c0[i] * c1[i];
// dft(c0, n); rev(c0, n);
// ff(i, 0, n) a[i] = qz(c0[i].x) % mo;
// return;
ff(i, 0, n) {
c0[i] = P(a[i] & 32767, a[i] >> 15);
c1[i] = P(b[i] & 32767, b[i] >> 15);
}
dft(c0, n); dft(c1, n);
ff(i, 0, n) {
int j = (n - i) & (n - 1);
P k, d0, d1, d2, d3;
k = conj(c0[j]);
d0 = (k + c0[i]) * P(0.5, 0);
d1 = (k - c0[i]) * P(0, 0.5);
k = conj(c1[j]);
d2 = (k + c1[i]) * P(0.5, 0);
d3 = (k - c1[i]) * P(0, 0.5);
c2[i] = d0 * d2 + d1 * d3 * P(0, 1);
c3[i] = d1 * d2 + d0 * d3;
}
dft(c2, n); dft(c3, n);
rev(c2, n); rev(c3, n);
ff(i, 0, n) {
a[i] = qz(c2[i].x) + (qz(c2[i].y) % mo << 30) + (qz(c3[i].x) % mo << 15);
a[i] %= mo;
}
}
void build() {
for(int i = 1; i < nm; i *= 2) ff(j, 0, i)
ntt :: w[i + j] = P(cos(pi * j / i), sin(pi * j / i));
}
} ll ksm(ll x, ll y) {
ll s = 1;
for(; y; y /= 2, x = x * x % mo)
if(y & 1) s = s * x % mo;
return s;
} V operator + (V a, V b) {
a.resize(max(a.si, b.si));
ff(i, 0, a.si) a[i] = (a[i] + b[i]) % mo;
return a;
}
V operator - (V a, V b) {
a.resize(max(a.si, b.si));
ff(i, 0, a.si) a[i] = (a[i] - b[i] + mo) % mo;
return a;
}
V operator * (V a, int x) {
ff(i, 0, a.si) a[i] = a[i] * x % mo;
return a;
}
V operator * (V a, V b) {
int sa = a.si + b.si - 1, n = 1;
for(; n <= sa; n *= 2);
a.resize(n); b.resize(n);
ntt :: fft(a, b);
a.resize(sa);
return a;
} V qni(V a) {
int n0 = 1; for(; n0 < a.si; n0 *= 2);
V b; b.resize(1); b[0] = ksm(a[0], mo - 2);
for(int n = 2; n <= n0; n *= 2) {
V c = a; c.resize(n);
c = c * b; c.resize(n); c = c * b; c.resize(n);
b = b * 2 - c;
}
b.resize(a.si);
return b;
} V qd(V a) {
ff(i, 1, a.si) a[i - 1] = a[i] * i % mo;
a.resize(a.si - 1);
return a;
}
V jf(V a) {
a.pb(0);
fd(i, a.si - 1, 1) a[i] = a[i - 1] * ksm(i, mo - 2) % mo;
a[0] = 0;
return a;
}
V ln(V a) {
int sa = a.si;
a = jf(qni(a) * qd(a));
a.resize(sa);
return a;
}
V exp(V a) {
int n0 = 1; for(; n0 < a.si; n0 *= 2);
V b; b.resize(1); b[0] = 1;
for(int n = 2; n <= n0; n *= 2) {
V c = b; c.resize(n);
V d = a; d.resize(n);
c = d - ln(c); c[0] ++;
b = b * c; b.resize(n);
}
b.resize(a.si);
return b;
} V a; const int N = 1e5 + 5; int n, k; ll ni[N]; int main() {
freopen("b.in", "r", stdin);
freopen("b.out", "w", stdout);
ntt :: build();
n = 1e5;
fo(i, 1, n) ni[i] = ksm(i, mo - 2);
scanf("%d %d", &n, &k);
a.clear(); a.resize(k + 1);
fo(j, 1, k) a[j] = ni[j] % mo * (n - 1) % mo;
fo(i, 2, n) {
fo(j, 1, k / i) a[i * j] -= ni[j];
}
fo(j, 1, k) a[j] %= mo;
a = exp(a);
pp("%lld\n", (a[k] + mo) % mo);
}