%%%tourist
问有多少个符合条件的集合,使
an/g−n
和
an/g
,
an−n
奇偶性不一样
要使
an/g−n
和
an/g
奇偶性不一样,只要n是奇数就行了
要使
an/g−n
和
an−n
奇偶性不一样,就是要使
an/g
和
an
奇偶性不一样,就是说,集合里最大的数在去掉集合的gcd后改变了奇偶性,因为奇数去掉因子不可能变成偶数,所以只能是偶数变成了奇数
所以我们可以枚举集合里最大的
an
恰好被
2k
整除,这时因为
g
要包含所有2,所以集合里的所有数都要被
令
f(b,n,p)
表示在
[1,b]
中选n个数,最大的数奇偶性是p的方案数
f(b,0,0)=1,f(b,0,1)=0
考虑将他转移到
f(2b,n,p)
可以得到
直接转移每个n的 f(b,n,p) ,一次的复杂度是 O(n2) 的,我们可以这样
其中上面的 F∗F 是卷积,用FFT优化(模数比较小,用FFT没有问题),可以做到 O(nlogn) 转移
对于 f(b,n,p) 到 f(b+1,n,p) 的转移,可以做到 O(n)
那么对于每个k,我们需要
O(nlogn logmaxa)
计算他的答案
于是总复杂度是
O(nlognlog2maxa)
的
不能通过这道题
但是我们可以发现,我们在计算一个k的同时也是在计算其他k的答案
于是有个
logmaxa
就省掉了
于是复杂度就变成
O(nlogn logmaxa)
了
code:
#include<set>
#include<map>
#include<deque>
#include<queue>
#include<stack>
#include<cmath>
#include<ctime>
#include<bitset>
#include<string>
#include<vector>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<climits>
#include<complex>
#include<iostream>
#include<algorithm>
#define ll long long
using namespace std;
const int maxn = 130000;
const double pi = acos(-1);
struct E
{
double x,y;
E(){}
E(const double _x,const double _y){x=_x;y=_y;}
}w[maxn],t0[maxn],t1[maxn]; int id[maxn];
inline E operator +(const E &x,const E &y){return E(x.x+y.x,x.y+y.y);}
inline E operator -(const E &x,const E &y){return E(x.x-y.x,x.y-y.y);}
inline E operator *(const E &x,const E &y){return E(x.x*y.x-x.y*y.y,x.y*y.x+x.x*y.y);}
ll Mod;
int b;
int N,ln,n;
void DFT(E s[],const int sig)
{
for(int i=0;i<N;i++) if(i<id[i]) swap(s[i],s[id[i]]);
for(int m=2;m<=N;m<<=1)
{
int t=m>>1,tt=N/m;
for(int j=0;j<N;j+=m)
{
for(int i=0;i<t;i++)
{
E wn=sig==1?w[i*tt]:w[N-i*tt];
E tx=s[j+i],ty=s[j+i+t]*wn;
s[j+i]=tx+ty;
s[j+i+t]=tx-ty;
}
}
}
if(sig==-1) for(int i=0;i<N;i++) s[i].x/=(double)N;
}
ll re=0;
ll f0[maxn],f1[maxn],temp[maxn];
ll nf0[maxn],nf1[maxn];
void solve()
{
f0[0]=1;
int now=0;
for(int i=30;i>=1;i--)
{
if(now&1)
{
for(int j=0;j<=n;j++) t0[j]=E((f0[j]+f1[j])%Mod,0);
for(int j=0;j<=n;j++) t1[j]=E(f1[j],0);
for(int j=n+1;j<N;j++) t0[j]=t1[j]=E(0,0);
DFT(t0,1); DFT(t1,1);
for(int j=0;j<N;j++) t0[j]=t0[j]*t1[j];
DFT(t0,-1);
for(int j=0;j<=n;j++) nf0[j]=(f0[j]+((ll)(t0[j].x+0.5)))%Mod;
for(int j=0;j<=n;j++) t0[j]=E((f0[j]+f1[j])%Mod,0);
for(int j=0;j<=n;j++) t1[j]=E(f0[j],0);
for(int j=n+1;j<N;j++) t0[j]=t1[j]=E(0,0);
DFT(t0,1); DFT(t1,1);
for(int j=0;j<N;j++) t0[j]=t0[j]*t1[j];
DFT(t0,-1);
for(int j=0;j<=n;j++) nf1[j]=(/*f1[j]+*/((ll)(t0[j].x+0.5))%Mod-f0[j]+Mod)%Mod;
}
else
{
for(int j=0;j<=n;j++) t0[j]=E((f0[j]+f1[j])%Mod,0);
for(int j=0;j<=n;j++) t1[j]=E(f0[j],0);
for(int j=n+1;j<N;j++) t0[j]=t1[j]=E(0,0);
DFT(t0,1); DFT(t1,1);
for(int j=0;j<N;j++) t0[j]=t0[j]*t1[j];
DFT(t0,-1);
for(int j=0;j<=n;j++) nf0[j]=(((ll)(t0[j].x+0.5))%Mod-f1[j]+Mod)%Mod;
for(int j=0;j<=n;j++) t0[j]=E((f0[j]+f1[j])%Mod,0);
for(int j=0;j<=n;j++) t1[j]=E(f1[j],0);
for(int j=n+1;j<N;j++) t0[j]=t1[j]=E(0,0);
DFT(t0,1); DFT(t1,1);
for(int j=0;j<N;j++) t0[j]=t0[j]*t1[j];
DFT(t0,-1);
for(int j=0;j<=n;j++) nf1[j]=(f1[j]+((ll)(t0[j].x+0.5))%Mod)%Mod;
}
for(int j=0;j<=n;j++) f0[j]=nf0[j],f1[j]=nf1[j];
now<<=1;
if((b&(1<<i))>0)
{
for(int j=n;j>=1;j--) f1[j]=(f1[j]+f0[j-1]+f1[j-1])%Mod;
now|=1;
}
for(int j=1;j<=n;j+=2) re+=f1[j];
re%=Mod;
}
}
int main()
{
scanf("%d%d%I64d",&n,&b,&Mod);
N=1,ln=0; while(N<=2*n) N<<=1,ln++;
for(int i=1;i<N;i++) id[i]=(id[i>>1]>>1)|((i&1)<<ln-1);
for(int m=2;m<=N;m<<=1)
{
int t=m>>1,tt=N/m;
for(int i=0;i<t;i++)
{
w[i*tt]=E(cos(2*pi*i/m),sin(2*pi*i/m));
w[N-i*tt]=E(cos(2*pi*i/m),sin(-2*pi*i/m));
}
}
solve();
printf("%I64d\n",re);
return 0;
}