bzoj3434 [Wc2014]时空穿梭

式子的推导网上有很多dalao写的都很好，就不再赘述了。

有两点是我没有想到的，再次写一下。

首先就是初始的式子。不知道直角三角形斜边上整数点个数=gcd(a,b)+1，更没有拓展到多维。

第二就是划出来60分的式子之后的优化，把后边的关于T的多项式的系数通过O(n2)时间算了出来$O(n^2*m)->O(n^3*\sqrt{m})$

十分巧妙。

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cmath>
#define N 100500
#define mod 10007
using namespace std;
int C
[25],mu
,h
,nn,prime
,tot,vis
,g[25]
,pw
[25],sum[25][25]
;
void init(){
C[0][0]=1;
for(int i=1;i<=100000;i++){
C[i][0]=1;
for(int j=1;j<=20&&j<=i;j++)
C[i][j]=(C[i-1][j-1]+C[i-1][j])%mod;
}
mu[1]=1;h[++nn]=1;
for(int i=2;i<=100000;i++){
if(!vis[i]){
prime[++tot]=i;
mu[i]=-1;
}
for(int j=1;j<=tot&&i*prime[j]<=100000;j++){
vis[i*prime[j]]=1;
if(i%prime[j]==0){
mu[i*prime[j]]=0;
break;
}
mu[i*prime[j]]=-mu[i];
}
if(mu[i])h[++nn]=i;
}
for(int k=2;k<=20;k++){
for(int i=1;i<=100000;i++)
for(int j=1;j<=nn&&h[j]*i<=100000;j++)
(g[k][i*h[j]]+=mu[h[j]]*C[i-1][k-2]+mod)%=mod;
}
for(int i=1;i<=100000;i++){
pw[i][0]=1;
for(int j=1;j<=11;j++)pw[i][j]=pw[i][j-1]*i%mod;
}
for(int k=1;k<=20;k++)
for(int j=0;j<=11;j++)
for(int i=1;i<=100000;i++)
sum[k][j][i]=(sum[k][j][i-1]+g[k][i]*pw[i][j]%mod)%mod;
}
int T,n,c,m[25],ans,f[25];
void work(int x){
memset(f,0,sizeof f);f[0]=1;
for(int i=1;i<=n;i++){
for(int j=i;j>0;j--)
f[j]=(f[j]*2%mod*m[i]%mod-f[j-1]*(m[i]/x+1)%mod+mod)%mod;
f[0]=f[0]*2*m[i]%mod;
}
}
int main(){
init();
scanf("%d",&T);
while(T--){
scanf("%d%d",&n,&c);ans=0;
for(int i=1;i<=n;i++)scanf("%d",&m[i]);
sort(m+1,m+n+1);
for(int i=1,pos;i<=m[1];i=pos+1){
pos=N;
for(int j=1;j<=n;j++)pos=min(pos,m[j]/(m[j]/i));
int now=1;
for(int j=1;j<=n;j++)
now=(now*(m[j]/i)%mod*5004%mod);
work(i);
for(int j=0;j<=n;j++)
(ans+=now*f[j]%mod*(sum[c][j][pos]-sum[c][j][i-1]+mod)%mod)%=mod;
}
printf("%d\n",ans);
}
return 0;
}

View Code