csp模拟轻飘飘的时间【杜教筛】【卢卡斯】

最新推荐文章于 2022-06-23 22:46:43 发布

gigo_64

最新推荐文章于 2022-06-23 22:46:43 发布

阅读量132

点赞数

分类专栏：莫比乌斯反演卢卡斯定理筛法

本文链接：https://blog.csdn.net/qq_42835815/article/details/102941895

版权

莫比乌斯反演同时被 3 个专栏收录

6 篇文章 0 订阅

订阅专栏

筛法

3 篇文章 0 订阅

订阅专栏

卢卡斯定理

2 篇文章 0 订阅

订阅专栏

传送门

杜教筛板子，只要学了并且还会推导就能随便做（除了卡常）。

$\sum_{i=1}^n\sum_{j=1}^m\binom{gcd(i,j)}{B}$

n,m范围1e10，B范围mod以内，mod9990017。

其实这个mod如此与众不同是在提示你可以卢卡斯。

$\sum_{d=1}\binom{d}{B}\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor}\sum_{k|gcd(i,j)}\mu(k)$

$\sum_{d=1}\binom{d}{B}\sum_{k=1}\mu(k)\left \lfloor \frac{n}{kd} \right \rfloor\left \lfloor \frac{m}{kd} \right \rfloor$

枚举T=kd

$\sum_{T=1}\sum_{d|T}\binom{d}{B}\mu(\frac{T}{d})\left \lfloor \frac{n}{T} \right \rfloor\left \lfloor \frac{m}{T} \right \rfloor$

发现 $\sum_{d|T}\binom{d}{B}\mu(\frac{T}{d})$ 实际上是组合数和莫比乌斯函数的狄利克雷卷积。

令 $F(T)=\sum_{d|T}\binom{d}{B}\mu(\frac{T}{d})$

$\sum_{T=1}F(T)\left \lfloor \frac{n}{T} \right \rfloor\left \lfloor \frac{m}{T} \right \rfloor$

所以目标是F(T)的前缀和。这东西都说了有莫比乌斯函数的卷积，根据套路 $\mu*I=\epsilon$ 我们直接卷上一个I。

提示：I函数任何时候取值都是1。

根据杜教筛公式 $S(n)=\sum_{i=1}^nG(i)-\sum_{i=2}^nS(\left \lfloor \frac{n}{i} \right \rfloor)$

其中 $G(n)=\binom{n}{B}$

所以我们要求一个组合数的一列和。这个一列和有个很好推导的公式就是 $\sum_{i=1}^nG(i)=\binom{n+1}{B+1}$

那公式就弄完了。根据杜教筛预处理最优复杂度证明，预处理 $n^{\frac{2}{3}}$ 的F（直接埃氏筛），然后记忆化搜索。

少取模，多卡常，最烦卡常了。

#include<bits/stdc++.h>
#include<tr1/unordered_map>
using namespace std;
using namespace tr1; 
#define in read()
#define LL long long
LL in{
	LL cnt=0,f=1;char ch=0;
	while(!isdigit(ch)){
		ch=getchar();
	}
	while(isdigit(ch)){
		cnt=cnt*10ll+ch-48ll;
		ch=getchar();
	}return cnt;
}
const int mod=9990017;
const int N=(int)pow(1e10,2.0/3);
int F0[N];
unordered_map<LL,int> F;
LL n,m,B; int fac[mod+1],ifac[mod+1];
inline int c(int x,int y){
	if(x<y)return 0;
	return (LL)fac[x]*ifac[y]%mod*ifac[x-y]%mod;
}
inline int C(int x,int y){ 
	if(x<y)return 0;
	return (LL)c(x%mod,y%mod)*((x>=mod)?C(x/mod,y/mod):1)%mod;
}
inline int ksm(int a,int b){
	int gu=1;
	while(b){
		if(b&1)
		gu=(LL)gu*a%mod;a=(LL)a*a%mod;b>>=1; 
	}return gu; 
}
inline void prepare(){
	fac[0]=1;
	for(register int i=1;i<mod;++i)fac[i]=(LL)fac[i-1]*i%mod;
	ifac[mod-1]=ksm(fac[mod-1],mod-2);
	for(register int i=mod-2;i>=1;i--)ifac[i]=(LL)ifac[i+1]*(i+1)%mod;ifac[0]=1;
	for(register int i=1;i<=N;++i)F0[i]=C(i,B);
	for(register int i=1;i<=N;++i){
		for(register int j=i*2;j<=N;j+=i){
			F0[j]=(F0[j]-F0[i]);if(F0[j]<0)F0[j]+=mod;

		}
	}
	for(register int i=1;i<=N;++i)F0[i]=(F0[i-1]+F0[i])%mod;
}
inline int S(LL x){
	if(x<=N)return F0[x];
	if(F.count(x))return F[x];
	int sum=C(x+1,B+1);
	for(register LL l=2,r;l<=x;l=r+1){
		r=x/(x/l);
		sum-=(LL)(r-l+1)*S(x/l)%mod;
		if(sum<0)sum+=mod;
	}return F[x]=sum;
}
signed main(){
	
	n=in;m=in;B=in;if(n>m)swap(n,m);
	prepare();
	int ans=0;
	for(register LL l=1,r;l<=n;l=r+1){
		r=min(n/(n/l),m/(m/l));
		ans=(ans+(LL)(S(r)-S(l-1))%mod*(n/l)%mod*(m/l)%mod)%mod;
	}printf("%d", (ans+mod)%mod);
		return 0;
}