博客园 首页 私信博主 显示目录 隐藏目录 管理 动画

ICPC2019 南京. E. Observation(思路)

题目链接


\(Description\)
\(f(d)\)表示空间中到原点距离为\(d\)的整点个数,给定\(L,R,k,p\),求

\[\sum_{d=L}^Rf(d)\ \mathbb{xor}\ k\mod p \]

\(L,R\leq 10^{13},R-L+1\leq 10^6\)

\(Solution\)
必然是找规律。OEIS可以直接找到规律,或者打表。

\(f\)都是\(6\)的倍数,令\(g(x)=\frac{f(x)}{6}\),可以发现\(g(x)\)为积性函数,即\(g(\prod p_i^{k^i})=\prod g(p_i^{k_i})\)
列出\(g(p)\ (p为质数)\)可以知道,若\(p\%4=1\)\(g(p)=p\);若\(p\%4=3\)\(g(p)=p+2\)
再看\(g(p^2)\ (p为质数)\)可以知道,若\(p\%4=1\)\(g(p^2)=p^2\);若\(p\%4=3\)\(g(p^2)=p^2+h(p)\),其中\(h(p)=8k\)\(k\)\(p\)为第几个\(\%4=3\)的质数/whl
考虑如何更好地表示\(h(p)\)
由OEIS结论可发现\((p-1)h(p)=2(p^2-1)\),即\(g(p^2)=p^2+\frac{2(p^2-1)}{p-1}\),且该结论对\(p^k\)同样成立。
所以

\[f(p^k)=\begin{cases}1,&p=2\\p^k,&p\%4=1\\p^k+\frac{2(p^k-1)}{p-1},&p\%4=3\end{cases} \]

求出\(\sqrt n\)以内的质数,对每个数就可以\(O(\log n)\)分解质因数求了。
复杂度\(O(Tn\log n)\)

PS:注意&的运算优先级比==低。。&运算要加括号!
\(ans\)求和过程中最多是\(10^6\)\(2\times 10^{13}\)的数,所以还是要取模!


//1261ms	20MB
#include <bits/stdc++.h>
#define pc putchar
#define gc() getchar()
#define pb emplace_back
typedef long long LL;
const int N=3162300,M=1e6+5;
const LL LIM=2e18;

int cnt,P[230000];
bool not_P[N];

inline LL read()
{
	LL now=0,f=1; char c=gc();
	for(;!isdigit(c);c=='-'&&(f=-1),c=gc());
	for(;isdigit(c);now=now*10+c-48,c=gc());
	return now*f;
}
void Init(const int n)
{
	for(int i=2; i<=n; ++i)
	{
		if(!not_P[i]) P[++cnt]=i;
		for(int j=1; j<=cnt && i*P[j]<=n; ++j)
		{
			not_P[i*P[j]]=1;
			if(!(i%P[j])) break;
		}
	}
}
void Solve()
{
	static LL A[M],f[M];

	LL L=read(),R=read(),K=read(),mod=read(),ans=0;
	if(!L) ans+=K, ++L;
	int n=(int)(R-L+1);
	for(int i=1; i<=n; ++i) A[i]=i+L-1,f[i]=1;

	for(int i=1,lim=::cnt; i<=lim&&P[i]<=R; ++i)
	{
		int p=P[i];
		LL now=L/p*p; now<L&&(now+=p);
		for(now-=L-1; now<=n; now+=p)
		{
			LL pk=1;
			while(!(A[now]%p)) A[now]/=p, pk*=p;
			if(p==2) pk=1;
			else if((p&3)==3) pk+=2*(pk-1)/(p-1);
			f[now]*=pk;
		}
	}
	for(int i=1; i<=n; ++i)
	{
		if((A[i]&3)==3) f[i]*=(A[i]+2);//本身/剩下为大质数 
		else f[i]*=A[i];
		ans+=((f[i]*6)^K)%mod;
		ans>=LIM&&(ans%=mod);//取模! 
	}
	printf("%lld\n",ans%mod);
}

int main()
{
	Init(N-2);
	for(int T=read(); T--; Solve());

	return 0;
}
posted @ 2021-03-13 15:57  SovietPower  阅读(181)  评论(0编辑  收藏  举报