BZOJ3516 国王奇遇记加强版 <扰动法>

Posted on Symbols count in article: 598 Reading time ≈ 2 mins.

Problem

国王奇遇记加强版

Time  Limit:  1  Sec\mathrm{Time\;Limit:\;1\;Sec}
Memory  Limit:  512  MB\mathrm{Memory\;Limit:\;512\;MB}

Description

Input

共一行，包括两个正整数NN和MM。

Output

共一行，为所求表达式的值对109+710^9+7取模的值。

Sample Input

1
5 3

Sample Output

1
36363

Hint

1≤N≤1091\le N\le10^9，1≤m≤10001\le m\le1000

标签：扰动法

Solution

扰动法转和式为递推式。
令

f(k)=∑i=1nikmif(k)=\sum_{i=1}^{n}i^km^i

用扰动法化简

(m−1)f(k)=∑i=1nikmi+1−∑i=1nikmi=∑i=1n+1(i−1)kmi−∑i=1nikmi=nkmn+1+∑i=1n((i−1)k−ik)mi=nkmn+1+∑i=1n((∑j=0k(kj)(−1)k−jij)−ik)mi=nkmn+1+∑i=1n∑j=0k−1(kj)(−1)k−jijmi=nkmn+1+∑j=0k−1(kj)(−1)k−j∑i=1nijmi=nkmn+1+∑j=0k−1(kj)(−1)k−jf(j)\begin{aligned} (m-1)f(k) &=\sum_{i=1}^{n}i^km^{i+1}-\sum_{i=1}^{n}i^km^i\\ &=\sum_{i=1}^{n+1}(i-1)^km^{i}-\sum_{i=1}^{n}i^km^i\\ &=n^km^{n+1}+\sum_{i=1}^{n}((i-1)^k-i^k)m^i\\ &=n^km^{n+1}+\sum_{i=1}^{n}((\sum_{j=0}^{k}\binom{k}{j}(-1)^{k-j}i^j)-i^k)m^i\\ &=n^km^{n+1}+\sum_{i=1}^{n}\sum_{j=0}^{k-1}\binom{k}{j}(-1)^{k-j}i^jm^i\\ &=n^km^{n+1}+\sum_{j=0}^{k-1}\binom{k}{j}(-1)^{k-j}\sum_{i=1}^{n}i^jm^i\\ &=n^km^{n+1}+\sum_{j=0}^{k-1}\binom{k}{j}(-1)^{k-j}f(j)\\ \end{aligned}

于是可以进行O(m2)O(m^2)的递推，注意特判m=1m=1的情况。

简单版见BZOJ3157，加强版见BZOJ4126。

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
#include <bits/stdc++.h>
#define MX 200
#define P 1000000007
using namespace std;
typedef long long lnt;
template <class T> inline void read(T &x) {
    x = 0; int c = getchar(), f = 1;
    for (; !isdigit(c); c = getchar()) if (c == 45) f = -1;
    for (; isdigit(c); c = getchar()) (x *= 10) += f*(c-'0');
}
int n, m; lnt f[60][MX+5], p[MX+5], fac[MX+5], inv[MX+5];
void init() {
    fac[0] = inv[0] = inv[1] = 1;
    for (int i = 1; i <= MX; i++) fac[i] = fac[i-1]*i%P;
    for (int i = 2; i <= MX; i++) inv[i] = (P-P/i*inv[P%i]%P)%P;
    for (int i = 2; i <= MX; i++) inv[i] = inv[i-1]*inv[i]%P;
}
lnt Pow(lnt x, int k) {
    lnt ret = 1;
    for (; k; k >>= 1, x = x*x%P)
        if (k&1) ret = ret*x%P;
    return ret;
}
lnt C(int n, int m) {return fac[n]*inv[m]%P*inv[n-m]%P;}
void calc(int n, int stp) {
    if (!n) return;
    if (n&1) {
        calc(n-1, stp+1);
        for (int i = 0; i <= m; i++) f[stp][i] = m;
        for (int i = 0; i <= m; i++) for (int j = 0; j <= i; j++)
            f[stp][i] = (f[stp][i]+m*C(i, j)%P*f[stp+1][j]%P)%P;
    } else {
        calc(n/2, stp+1); lnt pw = Pow(m, n/2); p[0] = 1;
        for (int i = 1; i <= m; i++) p[i] = p[i-1]*(n/2)%P;
        for (int i = 0; i <= m; i++) f[stp][i] = f[stp+1][i];
        for (int i = 0; i <= m; i++) for (int j = 0; j <= i; j++)
            f[stp][i] = (f[stp][i]+pw*C(i, j)%P*p[i-j]%P*f[stp+1][j]%P)%P;
    }
}
int main() {
    read(n), read(m);
    init(), calc(n, 0);
    printf("%lld\n", f[0][m]);
    return 0;
}