神仙题
题目大意:
定义
f
m
(
n
)
f_m(n)
fm(n)为可重集合S的个数,其中S满足S的元素都是m的非负整数次幂,并且S的元素和是n。
定义
g
m
1
(
n
)
=
f
m
(
n
)
,
g
m
k
(
n
)
=
∑
i
=
0
n
g
m
k
−
1
(
i
)
f
m
(
n
−
i
)
,
k
≥
2
g_m^1(n)=f_m(n),g_m^k(n)=\sum_{i=0}^ng_m^{k-1}(i)f_m(n-i),k\geq2
gm1(n)=fm(n),gmk(n)=∑i=0ngmk−1(i)fm(n−i),k≥2,求:
∑
i
=
0
n
g
m
(
i
)
\sum_{i=0}^ng_m(i)
∑i=0ngm(i)。
n
,
m
≤
1
0
18
,
k
≤
20
n,m\le10^{18},k\le20
n,m≤1018,k≤20
题解:
考虑
f
m
(
n
)
f_m(n)
fm(n)其实是,你有
1
,
m
,
m
2
,
m
3
,
…
1,m,m^2,m^3,\dots
1,m,m2,m3,…,用这些数字拼出
n
n
n的方案数,并且这些数字是无序的(即
m
m
m要出现在
1
1
1之后,
m
2
m^2
m2在
m
m
m之后)。
如果要求 g m k ( n ) g_m^k(n) gmk(n),那么其实是,你先将 n n n表示成 n = ∑ i = 1 k x i n=\sum_{i=1}^k x_i n=∑i=1kxi,然后每个 x i x_i xi用 1 , m , m 2 , m 3 , … 1,m,m^2,m^3,\dots 1,m,m2,m3,…表示出来,加起来就是,要用: 1 , m , m 2 , m 3 , … , 1 , m , m 2 , m 3 , … , 1 , m , m 2 , m 3 , … , … , … 1,m,m^2,m^3,\dots,1,m,m^2,m^3,\dots,1,m,m^2,m^3,\dots,\dots,\dots 1,m,m2,m3,…,1,m,m2,m3,…,1,m,m2,m3,…,…,…(总共k段)拼出n,然后拼的时候元素要是无序的(即尽管有些数字的值相同,但是标号是不同的,视为不同的数字,然后要满足标号不减)。
假定现在要求
g
m
k
(
n
)
g_m^k(n)
gmk(n)(一项)。
显然数字顺序不重要,可以重新写成
1
,
1
,
1
,
…
,
1
,
m
,
m
,
m
,
…
,
m
,
m
2
,
m
2
,
m
2
,
…
,
m
2
,
…
1,1,1,\dots,1,m,m,m,\dots,m,m^2,m^2,m^2,\dots,m^2,\dots
1,1,1,…,1,m,m,m,…,m,m2,m2,m2,…,m2,…,每段有k个。
那么考虑一个普及组级别的dp,设
F
(
i
,
j
)
F(i,j)
F(i,j)表示考虑了前i个数字,和为j的方案数。显然有(记第
i
i
i个数字是
a
i
a_i
ai)枚举当前数字用了几个:
F
(
i
,
j
)
=
∑
k
≥
0
F
(
i
−
1
,
j
−
k
a
i
)
F(i,j)=\sum_{k\geq0}F(i-1,j-ka_i)
F(i,j)=∑k≥0F(i−1,j−kai)。
然后可以注意到,假设我已经算出了
F
(
i
−
1
)
F(i-1)
F(i−1),那么因为我是重新排列过数字的顺序了,所以之后加入的数字都是
a
i
a_i
ai的倍数,那么第二维模
a
i
a_i
ai的余数就不会变了。最终只关心n,所以
F
(
i
,
j
)
F(i,j)
F(i,j)只关心
j
%
a
i
=
n
%
a
i
j\%a_i=n\%a_i
j%ai=n%ai的那些
j
j
j。
也就是说,若
j
j
j不满足上述条件可以从现在开始不管了。
那么既然余数只关心
n
%
a
i
n\%a_i
n%ai,相当于是个定值,因此现在只需要表示前面的倍数即可,即重新令:
H
(
i
,
j
)
=
F
(
i
,
j
a
i
+
n
%
a
i
)
H(i,j)=F(i,ja_i+n\%a_i)
H(i,j)=F(i,jai+n%ai)。
其转移可以直接写出:
H
(
i
,
j
)
=
F
(
i
,
j
a
i
+
n
%
a
i
)
=
∑
k
=
0
j
F
(
i
−
1
,
(
j
−
k
)
a
i
+
n
%
a
i
)
=
∑
k
=
0
j
F
(
i
−
1
,
k
a
i
+
n
%
a
i
)
=
∑
k
=
0
j
H
(
i
−
1
,
k
a
i
a
i
−
1
+
⌊
n
%
a
i
a
i
−
1
⌋
)
H(i,j)=F(i,ja_i+n\%a_i)=\sum_{k=0}^jF(i-1,(j-k)a_i+n\%a_i)=\sum_{k=0}^jF(i-1,ka_i+n\%a_i)\\=\sum_{k=0}^jH\left(i-1,k\frac{a_i}{a_{i-1}}+\left\lfloor\frac{n\%a_i}{a_{i-1}}\right\rfloor\right)
H(i,j)=F(i,jai+n%ai)=k=0∑jF(i−1,(j−k)ai+n%ai)=k=0∑jF(i−1,kai+n%ai)=k=0∑jH(i−1,kai−1ai+⌊ai−1n%ai⌋)
注意由于
{
a
t
}
\{a_{t}\}
{at}是重新排列过的,因此后面的a是前面的倍数,那么
a
i
a
i
−
1
\frac{a_i}{a_{i-1}}
ai−1ai要么是1,要么是m,总之是个只和i有关的常数,后面的那个也是。
也就是
H
(
i
,
j
)
=
∑
k
=
0
j
H
(
i
−
1
,
d
i
k
+
r
i
)
H(i,j)=\sum_{k=0}^jH(i-1,d_ik+r_i)
H(i,j)=∑k=0jH(i−1,dik+ri)。
有个结论是,上面那个东西,如果
H
(
i
−
1
)
H(i-1)
H(i−1)是个关于第二维的
t
t
t次多项式,那么
H
(
i
)
H(i)
H(i)是关于第二维的
t
+
1
t+1
t+1次多项式。
(其实你取
d
i
=
1
,
c
i
=
0
d_i=1,c_i=0
di=1,ci=0,就是经典的“多项式的前缀和还是多项式并且次数加一”)
然后来填一开始的
g
g
g的前缀和的坑。
如果从一开始就理解为,求前缀和是用那些数字拼出一个大小不超过n 的数字s,然后新增
x
k
+
1
=
n
−
s
x_{k+1}=n-s
xk+1=n−s这个数字,其方案数是1,或者说
x
k
+
1
x_{k+1}
xk+1只能使用
1
1
1(而不是
1
,
m
,
m
2
,
m
3
,
…
1,m,m^2,m^3,\dots
1,m,m2,m3,…)来拼成,那么就是用
k
+
1
k+1
k+1个
1
1
1和
k
k
k个
m
i
,
i
>
0
m^i,i>0
mi,i>0拼出n的方案数。
这样每次转移以及最后求答案的时候暴力拉格朗日插值即可。
注意最后要求的是
F
(
c
n
t
,
n
)
=
H
(
c
n
t
,
⌊
n
a
c
n
t
⌋
)
F(cnt,n)=H\left(cnt,\left\lfloor\frac{n}{a_{cnt}}\right\rfloor\right)
F(cnt,n)=H(cnt,⌊acntn⌋)
实现的时候要注意小心爆int/longlong之类的。
#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=a;i<=b;i++)
#define Rep(i,v) rep(i,0,(int)v.size()-1)
#define lint long long
#define mod 1000000007
#define ull unsigned lint
#define db long double
#define pb push_back
#define mp make_pair
#define fir first
#define sec second
#define gc getchar()
#define debug(x) cerr<<#x<<"="<<x
#define sp <<" "
#define ln <<endl
using namespace std;
typedef pair<int,int> pii;
const int N=100010,T=3000;
inline int P(int &x) { return x>=mod?x-=mod:x; }
int fac[T],facinv[T];
inline int fast_pow(int x,int k,int ans=1) { for(;k;k>>=1,x=(lint)x*x%mod) (k&1)?ans=(lint)ans*x%mod:0;return ans; }
inline int sol(int x,int s) { return (s&1)?(x?mod-x:0):x; }
inline int prelude(int n)
{
fac[0]=1;rep(i,1,n) fac[i]=(lint)fac[i-1]*i%mod;
facinv[n]=fast_pow(fac[n],mod-2);
for(int i=n-1;i>=0;i--) facinv[i]=(i+1ll)*facinv[i+1]%mod;
return 0;
}
struct poly{
vector<int> f,pre,suf,xs;int n;
inline int resize(int _n)
{
return n=_n+1,f.resize(n),pre.resize(n),suf.resize(n),xs.resize(n),0;
}
inline int PRE(int x) { return x<0?1:pre[x]; }
inline int SUF(int x) { return x>=n?1:suf[x]; }
inline int getxs()
{
rep(i,0,n-1) xs[i]=(lint)f[i]*facinv[i]%mod*facinv[n-1-i]%mod;
return 0;
}
inline int F(lint x)
{
int ans=0;if(x<n) return f[(int)x];x%=mod;
rep(i,0,n-1) pre[i]=PRE(i-1)*(x-i)%mod;
for(int i=n-1;i>=0;i--) suf[i]=SUF(i+1)*(x-i)%mod;
rep(i,0,n-1) ans+=sol((lint)xs[i]*PRE(i-1)%mod*SUF(i+1)%mod,n-i-1),P(ans);
return ans;
}
}h[T];lint mi[100],a[T];
int main()
{
lint n,m;int k,cnt=0,t=0;mi[0]=1;cin>>n>>m>>k;
while(mi[cnt]<=n/m) mi[cnt+1]=mi[cnt]*m,cnt++;
rep(i,0,cnt) rep(j,1,k) a[++t]=mi[i];
a[0]=1,prelude(t);
h[0].resize(0),h[0].f[0]=1,h[0].getxs();
rep(i,1,t)
{
h[i].resize(i);lint r=n%a[i]/a[i-1];
rep(j,0,i)
h[i].f[j]=(j?h[i].f[j-1]:0)+h[i-1].F(a[i]/a[i-1]*j+r),P(h[i].f[j]);
h[i].getxs();
}
return !printf("%d\n",h[t].F(n/a[t]));
}