HDU 5780 gcd

不会...看了题解这个公式就简单了：gcd(x^a-1,x^b-1)=x^gcd(a,b)-1。

那么只要枚举一下d，[1,n]内 gcd(a,b)=d的有几对，对答案做出的贡献为 pair[d]*(x^d-1)。

pair[d]容易计算：即[1,n/d]内 gcd(a,b)=1的有几对，那么就是2*(phi[n/1]+phi[n/2]+...+phi[n/d])-1。

因为n有1000000，T有300，这样做依然超时，还需优化。

因为计算的是[1,n/d]内 gcd(a,b)=1的对数，当n/d相同时，pair[d]相同，可以分块计算。

#pragma comment(linker, "/STACK:1024000000,1024000000")
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<vector>
#include<map>
#include<set>
#include<queue>
#include<stack>
#include<ctime>
#include<iostream>
using namespace std;
typedef long long LL;
const double pi=acos(-1.0),eps=1e-8;
void File()
{
    freopen("D:\\in.txt","r",stdin);
    freopen("D:\\out.txt","w",stdout);
}
inline int read()
{
    char c = getchar();  while(!isdigit(c)) c = getchar();
    int x = 0;
    while(isdigit(c)) { x = x * 10 + c - '0'; c = getchar(); }
    return x;
}

const LL mod=1000000007;
const int maxn=1000010;
LL phi[maxn],prime[maxn],sumPhi[maxn],sz;
bool check[maxn];

LL pow(LL a,LL b,LL p)
{
    LL res=1; a%=p;
    while(b)
    {
        if(b&1) res=(res*a)%p;
        a=(a*a)%p, b/=2;
    }
    return res;
}

void Eorue()
{
    phi[1]=1; sz=0;
    for(int i=2;i<maxn;i++)
    {
        if(!check[i]) { phi[i]=i-1; prime[sz++]=i; }
        for(int j=0;j<sz;j++)
        {
            if(i*prime[j]>=maxn)break;
            check[i*prime[j]]=true;
            if(i%prime[j]==0) { phi[i*prime[j]]=phi[i]*prime[j];break; }
            else phi[i*prime[j]]=phi[i]*(prime[j]-1);
        }
    }
    sumPhi[0]=0;
    for(int i=1;i<maxn;i++) sumPhi[i]=(sumPhi[i-1]+(LL)phi[i])%mod;
}

LL extend_gcd(LL a,LL b,LL &x,LL &y)
{
    if(a==0&&b==0) return -1;
    if(b==0){x=1;y=0;return a;}
    LL d=extend_gcd(b,a%b,y,x);
    y-=a/b*x;
    return d;
}

LL mod_reverse(LL a,LL n)
{
    LL x,y,d=extend_gcd(a,n,x,y);
    if(d==1) return (x%n+n)%n;
    else return -1;
}

int main()
{
    Eorue() ; int T; scanf("%d",&T);
    while(T--)
    {
        LL x; int n; scanf("%lld%d",&x,&n);
        if(x==1) { printf("0\n"); continue; }
        LL ans=0;
        int p=1;
        while(p<=n)
        {
            int pos=n/(n/p);
            LL F=((((pow(x,pos+1,mod)-pow(x,p,mod)+mod)%mod)*mod_reverse(x-1,mod))%mod-(pos-p+1)+mod)%mod;
            LL G=(2*sumPhi[n/p]-1)%mod;
            ans=(ans+F*G%mod)%mod;
            p=pos+1;
        }
        printf("%lld\n",ans);
    }
    return 0;
}