CCPCネットワーク試合HDU-6706 huntian oy(モビウス反転+デュ教ふるい+sum(i*phi(i))テンプレート)


リンク:http://acm.hdu.edu.cn/showproblem.php?pid=6706
題意:f(n,a,b)=\sum_{i=1}^n \sum_{j=1}^i gcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]\%(10^9+7)を求めます
公式の問題解:

CCPC网络赛 HDU-6706 huntian oy(莫比乌斯反演+杜教筛+sum(i*phi(i))模板)_第1张图片
 
杜教篩強推博客:https://www.cnblogs.com/peng-ym/p/9446555.html
 
#include  
#include 
#define ll long long
#define ull unsigned long long 
using namespace std;
const int N = 3e6+10;
const ll mod = 1e9+7;
ll phi[N],prime[N],cnt,ni2,ni6;
ll pre1[N];
tr1::unordered_map  mp1; 
templateinline void read(T &x)
{
    x=0;
    static int p;p=1;
    static char c;c=getchar();
    while(!isdigit(c)){if(c=='-')p=-1;c=getchar();}
    while(isdigit(c)) {x=(x<<1)+(x<<3)+(c-48);c=getchar();}
    x*=p;
}
void Init()
{
	cnt=0;
	phi[1]=1;
	for(int i=2;i<=N-10;i++)
	{
		if(!phi[i])
		{
			prime[cnt++]=i;
			phi[i]=i-1;
		}
		for(int j=0;j=0&&l<=n;l=r+1)
	{
		r=n/(n/l);
		ans-=(r-l+1)*(l+r)%mod*ni2%mod*S(n/l)%mod; 
		ans=(ans%mod+mod)%mod;
	} 
	return mp1[n]=ans; 
}
ll poww(ll a,ll b)
{
	ll ans=1;
	while(b)
	{
		if(b&1) ans=ans*a%mod;
		a=a*a%mod;	
		b>>=1;
	}
	return ans;
}
int main(void)
{
	Init();
	int t,n,a,b;
	ni2=poww(2,mod-2);
	ni6=poww(6,mod-2);
	scanf("%d",&t);
	while(t--)
	{
		read(n),read(a),read(b); 
		printf("%lld
",(S(n)-1+mod)%mod*ni2%mod); } return 0; }