【BZOJ4652】【NOI2016】循环之美

传送门： BZOJ4652

推导1

首先有一个神奇的结论，x/y （x，y互质）在k 进制下为纯循环小数当且仅当y与k互质，理由如下

当x/y纯循环，最短循环节长度为L，那么第二个循环节开头数字是x*(k^L)%y的结果,第一个循环节开头是x%y的结果

因此有x*(k^L)=x (mod y)

因为x与y互质，可以消去x，得k^L=1 (mod y)

容易得出(k,y)=1

推导2

（莫比乌斯反演，下面是Word手打图片）

最后一步中把y/d当做一个变量，可以看做是t=y/d，枚举t=1~m/d

按n/d和m/d可以将d分成O(sqrt(n))块

推导3

设

当使n/d=a,m/d=b的最大d为q，最小d为p//就是说分块出来是p到q

要向答案累加a*[S(q)-S(p-1)]*T(b)

S(x)和T(x)要通过递归求得

设k=p1^a1*p2^a2*...*pn^an

设s(i,x)为小于等于x的，与k的前i种质因子的乘积互质的自然数的μ函数的和

t(i,x)为小于等于x的，与k的前i种质因子的乘积互质的自然数的个数

于是S(x)=s(n,x)  T(x)=t(n,x)//这里n是k的质因子种数

递推方法s(i,x)=s(i−1,x)−μ(pi)s(i,⌊xpi⌋)

即减去与k的gcd为pi的幂的自然数的μ函数值之和

其中如果gcd的pi指数大于等于2，那么对应的数的μ函数为0

所以对函数的和有贡献的数只有除以一个pi以后与pi互质的，故由积性函数性质，可以由除pi以后的μ直接乘上pi的μ值（显然是-1）

所以化简一下其实是s(i,x)=s(i-1,x)+s(i,⌊xpi⌋)

类比但是更简单的，t的递推式
t(i,x)=t(i−1,x)−t(i−1,⌊xpi⌋)

其中t（i-1,x）是不考虑i及后面的质因子的互质数个数，t(i−1,⌊xpi⌋)是1到x中与前i种质因子乘积的gcd仅为pi的数的个数

p.s.递推公式参考philipsweng的博客

对于求过的s和t函数，建议map一下，否则可能会超时

递归底层i=0，t(0,x)=x，s需要求μ的前缀和

n<=1e9还挺大的，但是恰好有个东西叫杜教筛，专门求积性函数的前缀和，效率n^(2/3)，完结撒花

杜教筛讲解

附蒟蒻代码

#include<cstdio>
#include<algorithm>
#include<cstring>
#include<map>
using namespace std;
#define long long long
#define to i*prime[j]
#define maxn 1000005 
long n,m,k;
int p;
int miu[maxn];
int prime[maxn];
int sum[maxn];
int k_prf[10];
bool sift[maxn];
map<int,int> mp;
map<int,int> ms[10];
map<int,int> mt[10];
void init()
{
	miu[1]=sum[1]=1;
	for(int i=2;i<maxn;i++)
	{
		if(!sift[i])
		{
			prime[++p]=i;
			miu[i]=-1;
		}
		sum[i]=sum[i-1]+miu[i];
		for(int j=1;j<=p&&to<maxn;j++)
		{
			sift[to]=true;
			if(i%prime[j]==0) 
			{
				miu[to]=0;
				break;
			}
			else miu[to]=-miu[i];
		}
	}
} 
int sum_miu(int x)
{
	if(x<maxn) return sum[x];
	if(mp.count(x)) return mp[x];
	int ans=1;
	for(int i=2,p;p<x;i=p+1)
	{
		p=x/(x/i);
		ans-=sum_miu(x/i)*(p-i+1);
	}
	return mp[x]=ans;
}
long S(int x,int r)//sum_miu of 1~r coprime with the first x-th prmfactors of k
{
	if(!x) return sum_miu(r);
	if(r<=1) return r;
	if(ms[x].count(r)) return ms[x][r];
	return ms[x][r]=S(x-1,r)+S(x,r/k_prf[x]);
}
int T(int x,int r)//numbers of 1~r ...
{
	if(!x) return r;
	if(r<=1) return r;
	if(mt[x].count(r)) return mt[x][r];
	return mt[x][r]=T(x-1,r)-T(x-1,r/k_prf[x]);
}
int main()
{
 	init();
 	scanf("%lld%lld%lld",&n,&m,&k);
 	int kk=k;
 	p=0;
 	for(int i=1;kk>1;i++)
 	{
 		if(kk%prime[i]==0) 
 		{
 			k_prf[++p]=prime[i];
 			while(kk%prime[i]==0) 
 				kk/=prime[i];
		}
	 }
	 long ans=0ll;
	 for(int i=1,t;i<=n&&i<=m;i=t+1)
	 {
	 	t=min(m/(m/i),n/(n/i));
	 	ans+=(n/i)*(S(p,t)-S(p,i-1))*T(p,m/i);
	 }
	 printf("%lld",ans);
}

对于求过的s和t函数，建议map一下，否则可能会超时