题解
花了一上午+一中午终于把这道题A了
首先,我们要求的是bi互不相同的合法方案数
我们可以枚举一个a的集合S,来强制里面的b全部都相同,然后其它的随便放
由于这个题的n的约数非常多,我们可以把它质因数分解一下再来做
那么质因数分解之后怎么来算贡献呢?
我们可以强制每一个质因子都在S中分配相同的幂次
由于题目中的a在相等的情况下是无序的,而在不等的情况下又是有序的
所以为了简化问题,我们干脆把ai都看作是不同的
这样我们枚举集合就可以变成O(2^m)了(好像复杂度更高了。。。)
假设我们所有的质因子的可以随意的分配(什么都不考虑)
那么这样的方案数如果用生成函数来表示就是
(最后的[x^{e_j}]表示取x^{e_j}的系数)
显然这个生成函数是可以利用完全背包来计算的(m个物品,每个物品大小为ai,求花费背包空间为ej的方案数)
如果我们强制了一个ai的集合S必须选同样b
那么计算答案的式子就是(为了简便,我们把生成函数写成闭形式)
(相当于把S强制捆绑成一个大物品,物品大小为\sum ai)
然后我们自然而然地在外面套上了一个容斥系数(-1)^{|S|-1}
以为这样就可以拿到60分了
然而发现这样过不了样例。。。
因为我们容斥的时候只关注了bi互不相同,而这种情况
2 2 2 2 ({}中的元素表示该元素位于S中,里边的数表示集合中的ai,注意,ai相等但是位置不同我们也把它视为不同S)
按照我们的计算方法
+{2} 2 2 2+2 {2} 2 2+2 2 {2} 2+2 2 2 {2}
- {2 2} 2 2 - {2} 2 {2} 2 - {2} 2 2 {2} - 2 {2 2} 2 - 2 {2} 2 {2} - 2 2 {2 2}
+{2 2 2} 2 +{2 2} 2 {2}+{2} 2 {2 2}+2 {2 2 2}
-{2 2 2 2}
我们把S强制捆绑成一个大物品后,计算的答案就是(用物品大小的可重集[ ]来代表背包的结果)
假设tot=1,p1=2,e1=8
+[2 2 2 2]*4=35*4
-[4 2 2]*6=-9*6
+[6 2]*4=2*4
-[8]*1=-1
=93
而原问题的答案显然是0(即怎么分配2都会出现bi相等)
到底是哪里多算呢了?
我们发现,枚举的S集合(这其实是一个可重集),虽说内部的ai编号是不同的,但是S本身是可能被重复枚举到的
也就是说,我们把相同构成的可重集S多算了一遍
这道题其实还需要对可重集S进行容斥
我们枚举一个对a的集合划分P={S1,S2,...S_{|P|}}(这里就相当于在枚举一个集合的集合,只不过这些集合之间无交集,且他们的并集为全集a)
那么我们也可以写出答案的式子
然后又自然而然地乘上一个(-1)^{|P|-1}的容斥系数?
还是过不了样例啊
我们觉得这个容斥系数很烦人,于是先设它为F[P],把它放一边
略加思考
略加思考
略加思考
深入思考
胡乱思考
这个F[P]好像并不好算啊
那我们自己来构造一个F[P]吧
我们想一想,什么时候答案才可以+1呢?
很显然,所有集合的大小都为1的时候
那么如果P中存在一个集合Si大小大于1,则我们要它的贡献为0
(这里我们似乎发现了容斥的本质:通过对所有贡献的加权来获取我们需要的贡献)
于是我们可以定义容斥系数F[P]
这样就完了???
当然没有
一个集合S可能会被计算多次,我们需要保证的是S在累加完所有的贡献之后的系数为[n=1]
即[n=1]=f[n]+(???*f[???]+??*f[??]+???)
我们可以枚举一下S中1号点的所属集合的大小来计算后面的一堆问号
(后面的[n-i=1]表示子问题n-i的系数)
(关于定1号点已经是组合数学的老套路了)
所以最终的式子就是
推导一下(令i=n-1)
所以f[1]=1,f[n]=-(n-1)*(-(n-2))*(-(n-3))*...*(-1)=(-1)^(n-1)*(n-1)!
最终我们完成了容斥系数的构造
所以最终的答案就是
考虑一下,我们可以枚举每一个划分来计算答案
然而划分方案太多了,我们需要简化一下
我们发现对于一个划分P={S1,S2,...,Sk}
我们用到的只有每个集合Si的ai之和,以及每个集合的大小|Si|,把这些ai之和拿来做背包,把Si拿来计算容斥系数
(这里的意思是把sumai与Si的组成相同的P都压缩到一起来求解,最后乘上组成这种P的方案数即可)
(我们可以用结构体+hash来表示一个划分P,利用map来存一个结构体的DP值)
(整个程序就分为四个部分:
1、预处理
2、枚举划分P,并把可以一起计算的划分压缩在一起,构成一个压缩划分集
3、求出每一组压缩划分集的容斥系数(由于Si是一样的,所以容斥系数也相等),并与该压缩划分集的大小相乘
4、对每一组压缩划分集求出背包数组,并算出每个质因子的答案,与它所带的系数相乘,求和就是答案)
我们还可以稍微优化第四步:
把sumai相同的压缩划分集再压缩一次,构成一个二次压缩划分集(可以用multiset来存)
实测压缩划分集最坏大概有57000左右,二次压缩划分集最多只有[30的整数划分=5604]个
所以我们的平均加入的物品个数为O(5604*E(|P|))
每一个物品我们都要花O(maxtj)的复杂度加入背包
所以总复杂度是O(5604*E(|P|)maxtj)
实测E(|P|)≈10.3,maxtj=9995
T到飞起
我们稍微运用一下ctime,发现大部分时间都花在了背包上,而背包的计算次数高达5亿次!!
怎么办呢
开O2:3.4s >> 2.8s
加法取模优化:2.8s >> 2.67s
循环展开:2.67s >> 2.55s(为什么只优化了这么点儿啊。。。)
开小结构体的内部数组:TLE >> WA(不知道哪里出了问题。。。估计成功了也没有什么大优化。。。)
后来中午翘掉午休到机房继续肝这道题,想到了一个奇妙的做法
我把每一个二次压缩划分集的物品大小输出了出来
发现有很多物品是重复的
好像我们有很多重复的计算的呢。。
如果我们把物品大小序列看成一个字符串,把它们插入到一棵Trie树中,dfs一遍,每走一条边就加入一个物品
那么一些前缀相同的字符串的可以只计算一次该前缀
这样就省去了不少重复加入的物品
由于一个二次压缩划分集最多只有30个物品
我们可以对每一个深度开一个DP,每次更新的时候memcpy之前的DP状态就可以了
这样做的物品数大概为12000+,时间复杂度为12000*10000左右,2s+O2是可以通过的
代码1:(TLE,注释满满的代码)
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>
#include<set>
using namespace std;
#define M 10005
#define N 32
#define LL unsigned long long
#define fi first
#define se second
const int mod=1000000007;
bool vis[M];
int prime[M],tot,fac[N],finv[N];
int cnt[M];//mei zhong zhi yin zi de chu xian ci shu
int ksm(int x,int y)
{
int ret=1;
while(y){
if(y&1)ret=1ll*ret*x%mod;
y>>=1;x=1ll*x*x%mod;
}
return ret;
}
void shai(int n)
{
int i,j;vis[1]=1;
for(i=2;i<=n;i++){
if(!vis[i])prime[++tot]=i;
for(j=1;j<=tot;j++){
int tmp=i*prime[j];if(tmp>n)break;
vis[tmp]=1;if(i%prime[j]==0)break;
}
}
for(i=1;i<=tot;i++)
for(j=prime[i];j<=n;j*=prime[i])
cnt[i]+=n/j;
fac[0]=fac[1]=finv[0]=finv[1]=1;
for(i=2;i<=30;i++)fac[i]=1ll*fac[i-1]*i%mod;
finv[30]=ksm(fac[30],mod-2);
for(i=30;i>=1;i--)finv[i-1]=1ll*finv[i]*i%mod;
}
int a[N],tong[N];//tong:mei zhong a de chu xian ci shu(yong lai zui hou /fac[tong[i]])
struct P{
//ya suo hou de hua fen:ba suma he size qu zhi xiang tong de hua fen fang dao yi qi
int len;LL hh1,hh2;
pair<int,int> a[32];//fi:suma,se:size
P(){len=0;hh1=0;memset(a,0,sizeof(a));}
void gethh(){
sort(a+1,a+len+1);hh1=0;
for(int i=1;i<=len;i++)//{
hh1=hh1*3993991+(LL)a[i].fi+(LL)a[i].se*1997;
//hh2=hh2;
//}
}
bool operator < (const P &t)const{
//if(len!=t.len)return len<t.len;
return hh1<t.hh1;
//if(hh2!=t.hh2)return hh2<t.hh2;
//for(int i=1;i<=len;i++){
// if(a[i].fi!=t.a[i].fi)return a[i].fi<t.a[i].fi;
// if(a[i].se!=t.a[i].se)return a[i].se<t.a[i].se;
//}
//return 0;
}
}pa,pb;
map<P,int> mp[2];//tong ji mei zhong ya suo hua fen de chu xian ci shu
map<P,int>::iterator it;
multiset<int> S;//ba suma xiang tong de zai dan du ya suo yi xia
//zhe yang ke yi jia kuai bei bao d ci shu
multiset<int>::iterator its;
map<multiset<int>,int> F;//hua fen de rcxs deng yu ji he de rcxs zhi ji
//F[P]=\prod (-1)^(|Si|-1)*(|Si|-1)!
map<multiset<int>,int>::iterator itf;
int rc(int x){return x&1?mod-1:1;}
int f[M];//zuo wan quan bei bao:O(n*m)
//hao xiang ke yi duo xiang shi qiu ni zai juan ji zuo dao O(m*nlogn)---->sb
#include<ctime>
int main()
{
int n,m,i,j;
scanf("%d%d",&n,&m);
for(i=1;i<=m;i++){scanf("%d",&a[i]);tong[a[i]]++;}
shai(n);sort(a+1,a+m+1);
int now=0;
mp[now][pa]=1;
double c1=clock();
for(i=1;i<=m;i++){
mp[now^1].clear();
for(it=mp[now].begin();it!=mp[now].end();it++){
pa=(*it).fi;
for(j=1;j<=pa.len+1;j++){
pb=pa;
pb.a[j].fi+=a[i];pb.a[j].se++;
if(j==pa.len+1)pb.len++;
pb.gethh();
mp[now^1][pb]=(mp[now^1][pb]+(*it).se)%mod;
}
}
now^=1;
}
printf("%d:%.3f\n",mp[now].size(),(clock()-c1)/1000);
//int ttmp[55]={0},tcnt=0;
for(it=mp[now].begin();it!=mp[now].end();it++){
pa=(*it).fi;int mulf=(*it).se;
//ttmp[++tcnt]=mulf;
S.clear();
for(i=1;i<=pa.len;i++){
S.insert(pa.a[i].fi);
mulf=1ll*mulf*rc(pa.a[i].se-1)%mod*fac[pa.a[i].se-1]%mod;
}
F[S]=(F[S]+mulf)%mod;//lei jia deng xiao ji he de rcxs
//printf("%d\n",F[S]);
}
printf("%d:%.3f\n",F.size(),(clock()-c1)/1000);
//int T = 0;
int ans=0;//ct=0;
for(itf=F.begin();itf!=F.end();itf++){
//ct++;
S=(*itf).fi;
//ttmp[++tcnt]=(*itf).se;
memset(f,0,sizeof(f));f[0]=1;
for(its=S.begin();its!=S.end();its++){
int v=*its;
for(i=v;i<=cnt[1];i++){
f[i]+=f[i-v];
if(f[i]>=mod)f[i]-=mod;
}
//T+=cnt[1]-v;
}
int mulf=1;
for(i=1;i<=tot;i++)
mulf=1ll*mulf*f[cnt[i]]%mod;
ans=(1ll*ans+1ll*(*itf).se*mulf)%mod;
//printf("%d: %d\n",ct,ans);
}
//printf("-- %d \n",T);
//sort(ttmp+1,ttmp+tcnt+1);
//for(i=1;i<=tcnt;i++)
// printf("%d\n",ttmp[i]);
for(i=1;i<=30;i++)if(tong[i])
ans=1ll*ans*finv[tong[i]]%mod;
printf("%d\n",ans);
printf("%.3f\n",(clock()-c1)/1000);
}
代码2:(TLE,加入了无数卡常技巧的代码)
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>
#include<set>
using namespace std;
#define M 10005
#define N 31
#define LL unsigned long long
#define fi first
#define se second
const int mod=1000000007;
bool vis[M];
int prime[M],tot,fac[N],finv[N];
int cnt[M];
int ksm(int x,int y)
{
int ret=1;
while(y){
if(y&1)ret=1ll*ret*x%mod;
y>>=1;x=1ll*x*x%mod;
}
return ret;
}
void shai(int n)
{
int i,j;vis[1]=1;
for(i=2;i<=n;i++){
if(!vis[i])prime[++tot]=i;
for(j=1;j<=tot;j++){
int tmp=i*prime[j];if(tmp>n)break;
vis[tmp]=1;if(i%prime[j]==0)break;
}
}
for(i=1;i<=tot;i++)
for(j=prime[i];j<=n;j*=prime[i])
cnt[i]+=n/j;
fac[0]=fac[1]=finv[0]=finv[1]=1;
for(i=2;i<=30;i++)fac[i]=1ll*fac[i-1]*i%mod;
finv[30]=ksm(fac[30],mod-2);
for(i=30;i>=1;i--)finv[i-1]=1ll*finv[i]*i%mod;
}
int a[N],tong[N];
struct P{
int len;LL hh1,hh2;
pair<int,int> a[N];
P(){len=0;hh1=0;memset(a,0,sizeof(a));}
void gethh(){
sort(a+1,a+len+1);hh1=0;
for(int i=1;i<=len;i++)
hh1=hh1*3993991+a[i].fi+1997ll*a[i].se;
}
bool operator < (const P &t)const{return hh1<t.hh1;}
}pa,pb;
map<P,int> mp[2];
map<P,int>::iterator it;
multiset<int> S;
multiset<int>::iterator its;
map<multiset<int>,int> F;
map<multiset<int>,int>::iterator itf;
//int rc(int x){return x&1?mod-1:1;}
//int f[M];
//#include<ctime>
int ans;
struct trie{
int ch[31],k;
}tr[100005];
int trtot,rt;
int tmp[31],tcnt;
void insert(int &i,int pos,int k)
{
if(!i)i=++trtot;
if(pos>tcnt){
tr[i].k=k;
return;
}
insert(tr[i].ch[tmp[pos]],pos+1,k);
}
int f[32][M];
void BAG(int i,int dep)
{
if(tr[i].k){
int mulf=1;
for(int j=1;j<=tot;j++)
mulf=1ll*mulf*f[dep-1][cnt[j]]%mod;
ans=(1ll*ans+1ll*tr[i].k*mulf)%mod;
return;
}
for(int j=1;j<=30;j++){
if(!tr[i].ch[j])continue;
memcpy(f[dep],f[dep-1],sizeof(f[dep-1]));
for(int k=j;k<=cnt[1];k++){
f[dep][k]+=f[dep][k-j];
if(f[dep][k]>=mod)f[dep][k]-=mod;
}
BAG(tr[i].ch[j],dep+1);
}
}
int main()
{
//freopen("B.out","w",stdout);
int n,m,i,j;
scanf("%d%d",&n,&m);
for(i=1;i<=m;i++){scanf("%d",&a[i]);tong[a[i]]++;}
shai(n);sort(a+1,a+m+1);
int now=0;
mp[now][pa]=1;
//double c1=clock();
for(i=1;i<=m;i++){
mp[now^1].clear();
for(it=mp[now].begin();it!=mp[now].end();it++){
pa=(*it).fi;
for(j=1;j<=pa.len+1;j++){
pb=pa;
pb.a[j].fi+=a[i];pb.a[j].se++;
if(j==pa.len+1)pb.len++;
pb.gethh();
mp[now^1][pb]=(mp[now^1][pb]+(*it).se)%mod;
}
}
now^=1;
}
//printf("%d:%.3f\n",mp[now].size(),(clock()-c1)/1000);
for(it=mp[now].begin();it!=mp[now].end();it++){
pa=(*it).fi;int mulf=(*it).se;
S.clear();
for(i=1;i<=pa.len;i++){
S.insert(pa.a[i].fi);
mulf=1ll*mulf*fac[pa.a[i].se-1]%mod;
if(!(pa.a[i].se&1))mulf=1ll*mulf*(mod-1)%mod;
}
F[S]=(F[S]+mulf)%mod;
}
//printf("%d:%.3f\n",F.size(),(clock()-c1)/1000);
//double T1=clock(),T2=0;
rt=0;trtot=0;
for(itf=F.begin();itf!=F.end();itf++){
S=(*itf).fi;tcnt=0;
for(its=S.begin();its!=S.end();its++)
tmp[++tcnt]=*its;
insert(rt,1,(*itf).se);
}
f[0][0]=1;
BAG(rt,1);
/*int ans=0;int con=0;
for(itf=F.begin();itf!=F.end();itf++){
S=(*itf).fi;
memset(f,0,sizeof(f));f[0]=1;
for(its=S.begin();its!=S.end();its++){
int v=*its;
double cd=clock();
printf("%d ",v);
//con++;
for(i=v;i<=cnt[1];i++){
int t=i-v;
f[i]+=f[t];if(f[i]>=mod)f[i]-=mod;
if(i+8<=cnt[1]){
f[i+1]+=f[t+1];if(f[i+1]>=mod)f[i+1]-=mod;
f[i+2]+=f[t+2];if(f[i+2]>=mod)f[i+2]-=mod;
f[i+3]+=f[t+3];if(f[i+3]>=mod)f[i+3]-=mod;
f[i+4]+=f[t+4];if(f[i+4]>=mod)f[i+4]-=mod;
f[i+5]+=f[t+5];if(f[i+5]>=mod)f[i+5]-=mod;
f[i+6]+=f[t+6];if(f[i+6]>=mod)f[i+6]-=mod;
f[i+7]+=f[t+7];if(f[i+7]>=mod)f[i+7]-=mod;
f[i+8]+=f[t+8];if(f[i+8]>=mod)f[i+8]-=mod;
i+=8;
}
}
T2+=clock()-cd;
}
printf("\n");
int mulf=1;
for(i=1;i<=tot;i++)
mulf=1ll*mulf*f[cnt[i]]%mod;
ans=(1ll*ans+1ll*(*itf).se*mulf)%mod;
}*/
for(i=1;i<=30;i++)if(tong[i])
ans=1ll*ans*finv[tong[i]]%mod;
printf("%d\n",ans);
//printf("%d\n",con);
//printf("%.3f\n",(clock()-c1)/1000);
//printf("other:%.3f\nbag:%.3f\n",(clock()-T1-T2)/1000,T2/2);
}
代码3:(AC 总时间:6000+ms最快代码,并删掉了许多注释)
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>
#include<set>
using namespace std;
#define M 10005
#define N 31
#define LL unsigned long long
#define fi first
#define se second
const int mod=1000000007;
bool vis[M];
int prime[M],tot,fac[N],finv[N];
int cnt[M];
int ksm(int x,int y)
{
int ret=1;
while(y){
if(y&1)ret=1ll*ret*x%mod;
y>>=1;x=1ll*x*x%mod;
}
return ret;
}
void shai(int n)
{
int i,j;vis[1]=1;
for(i=2;i<=n;i++){
if(!vis[i])prime[++tot]=i;
for(j=1;j<=tot;j++){
int tmp=i*prime[j];if(tmp>n)break;
vis[tmp]=1;if(i%prime[j]==0)break;
}
}
for(i=1;i<=tot;i++)
for(j=prime[i];j<=n;j*=prime[i])
cnt[i]+=n/j;
fac[0]=fac[1]=finv[0]=finv[1]=1;
for(i=2;i<=30;i++)fac[i]=1ll*fac[i-1]*i%mod;
finv[30]=ksm(fac[30],mod-2);
for(i=30;i>=1;i--)finv[i-1]=1ll*finv[i]*i%mod;
}
int a[N],tong[N];
struct P{
int len;LL hh1,hh2;
pair<int,int> a[N];
P(){len=0;hh1=0;memset(a,0,sizeof(a));}
void gethh(){
sort(a+1,a+len+1);hh1=0;
for(int i=1;i<=len;i++)
hh1=hh1*3993991+a[i].fi+1997ll*a[i].se;
}
bool operator < (const P &t)const{return hh1<t.hh1;}
}pa,pb;
map<P,int> mp[2];
map<P,int>::iterator it;
multiset<int> S;
multiset<int>::iterator its;
map<multiset<int>,int> F;
map<multiset<int>,int>::iterator itf;
int ans;
struct trie{
int ch[31],k;
}tr[60005];
int trtot,rt;
int tmp[31],tcnt;
void insert(int &i,int pos,int k)
{
if(!i)i=++trtot;
if(pos>tcnt){tr[i].k=k;return;}
insert(tr[i].ch[tmp[pos]],pos+1,k);
}
int f[32][M];
void BAG(int i,int dep)
{
if(tr[i].k){
int mulf=1;
for(int j=1;j<=tot;j++)
mulf=1ll*mulf*f[dep-1][cnt[j]]%mod;
ans=(1ll*ans+1ll*tr[i].k*mulf)%mod;
return;
}
for(int j=1;j<=30;j++){
if(!tr[i].ch[j])continue;
memcpy(f[dep],f[dep-1],sizeof(f[dep-1]));
for(int k=j;k<=cnt[1];k++){
f[dep][k]+=f[dep][k-j];
if(f[dep][k]>=mod)f[dep][k]-=mod;
}
BAG(tr[i].ch[j],dep+1);
}
}
int main()
{
int n,m,i,j;
scanf("%d%d",&n,&m);
for(i=1;i<=m;i++){scanf("%d",&a[i]);tong[a[i]]++;}
shai(n);sort(a+1,a+m+1);
int now=0;
mp[now][pa]=1;
for(i=1;i<=m;i++){
mp[now^1].clear();
for(it=mp[now].begin();it!=mp[now].end();it++){
pa=(*it).fi;
for(j=1;j<=pa.len+1;j++){
pb=pa;
pb.a[j].fi+=a[i];pb.a[j].se++;
if(j==pa.len+1)pb.len++;
pb.gethh();
mp[now^1][pb]=(mp[now^1][pb]+(*it).se)%mod;
}
}
now^=1;
}
for(it=mp[now].begin();it!=mp[now].end();it++){
pa=(*it).fi;int mulf=(*it).se;
S.clear();
for(i=1;i<=pa.len;i++){
S.insert(pa.a[i].fi);
mulf=1ll*mulf*fac[pa.a[i].se-1]%mod;
if(!(pa.a[i].se&1))mulf=1ll*mulf*(mod-1)%mod;
}
F[S]=(F[S]+mulf)%mod;
}
rt=0;trtot=0;
for(itf=F.begin();itf!=F.end();itf++){
S=(*itf).fi;tcnt=0;
for(its=S.begin();its!=S.end();its++)
tmp[++tcnt]=*its;
insert(rt,1,(*itf).se);
}
f[0][0]=1;
BAG(rt,1);
for(i=1;i<=30;i++)if(tong[i])
ans=1ll*ans*finv[tong[i]]%mod;
printf("%d\n",ans);
}
\(^o^)/~\(^o^)/~\(^o^)/~完结撒花\(^o^)/~\(^o^)/~\(^o^)/~