BZOJ 2693: jzptab [莫比乌斯反演 线性筛]

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;
const int N=1e7+5,MOD=100000009;
inline int read() {
    char c=getchar();
    int x=0,f=1;
    while(c<‘0‘||c>‘9‘){if(c==‘-‘)f=-1;c=getchar();}
    while(c>=‘0‘&&c<=‘9‘){x=x*10+c-‘0‘;c=getchar();}
    return x*f;
}
int n,m;
bool notp[N];int p[N];
ll s[N],mu[N],g[N];
void sieve(){
    mu[1]=1;
    g[1]=1;
    for(int i=2;i<N;i++){
        if(!notp[i]) p[++p[0]]=i,mu[i]=-1,g[i]=i-(ll)i*i;
        for(int j=1;j<=p[0]&&i*p[j]<N;j++){
            int t=i*p[j];
            notp[t]=1;
            if(i%p[j]==0){
                mu[t]=0;
                g[t]=(g[i]*p[j])%MOD;
                break;
            }
            mu[t]=-mu[i];
            g[t]=(g[i]*g[p[j]])%MOD;
        }
    }
    
    for(int i=1;i<N;i++) g[i]=(g[i]+g[i-1])%MOD;    

}
inline ll S(ll x,ll y){
    return ((x*(x+1)/2)%MOD)*((y*(y+1)/2)%MOD)%MOD;
}
int main(){
    sieve();
    int T=read();
    while(T--){
        n=read();
        m=read();
        if(n>m) swap(n,m);
        ll ans=0,r=0;
        for(ll D=1;D<=n;D=r+1){
            r=min(n/(n/D),m/(m/D));
            ans=(ans+S(n/D,m/D)*(g[r]-g[D-1]))%MOD;
        }
        printf("%lld\n",(ans+MOD)%MOD);
    }
}