BZOJ 2693 jzptab 莫比乌斯反演

BZOJ 2693 jzptab 莫比乌斯反演
题目大意：给定n,m，求i从1到n，j从1到m，的i与j的最小公倍数之和。
这题真的是有问题，难想的一批，公式恐惧症无药可救患者。。。。。。
以下我们继续给出PoPoQQQ的PPT中的内容，当中做出一些解释。

注解：这里的F(x,y)要求gcd（i,j）==1，所以我们想要求解F(x,y)需要继续莫比乌斯反演。

注解：这里是一个标准的莫比乌斯反演，前面乘以一个i^2只是为了将除掉的i乘回来。

注解：我们注意到将公式全部打开之后出现了两个整除，那么我们就通过换元来将di变成一个与第二维没有关系的东西从而提出来变成正常的前缀和形式，从而求解。

注解：积性函数即为满足f(i*j)==f(i)*f(j)(gcd(i,j)==1)的函数，而且后面的式子为标准的狄利克雷卷积形式，积性函数乘以积性函数等于积性函数，两个积性函数的狄利克雷卷积也是积性函数，我们可以利用这种性质来线性筛出我们想要求出的那一坨，具体办法是质数手算，一个数中不包含另一个质数时利用两个函数值相乘，一个数中包含另一个数时找规律推出应该是什么（能不能推出来就看脸吧QAQ）

#include <cstdio>  
#include <cstring>  
#include <iostream>  
#include <algorithm>  
#define M 10001000  
#define MOD 100000009  
using namespace std;  
typedef long long ll;  
int prime[1001001],tot;  
bool not_prime[M];  
ll h[M],sum[M];  
void Linear_Shaker()  
{  
    int i,j;  
    h[1]=1;  
    for(i=2;i<M;i++)  
    {  
        if(!not_prime[i])  
        {  
            prime[++tot]=i;  
            h[i]=(i-(ll)i*i)%MOD;  
        }  
        for(j=1;prime[j]*i<M;j++)  
        {  
            not_prime[prime[j]*i]=1;  
            if(i%prime[j]==0)  
            {  
                h[prime[j]*i]=(prime[j]*h[i])%MOD;  
                break;  
            }  
            h[prime[j]*i]=(h[prime[j]]*h[i])%MOD;  
        }  
    }  
    for(i=1;i<M;i++)  
        sum[i]=(sum[i-1]+h[i])%MOD;  
}  
inline ll Sum(ll x,ll y)  
{  
    x%=MOD;y%=MOD;  
    ll re1=(x*(x+1)>>1)%MOD;  
    ll re2=(y*(y+1)>>1)%MOD;  
    return re1*re2%MOD;  
}  
int Query(int n,int m)  
{  
    int i,last,re=0;  
    if(n>m) swap(n,m);  
    for(i=1;i<=n;i=last+1)  
    {  
        last=min(n/(n/i),m/(m/i));  
        re+=Sum(n/i,m/i)*(sum[last]-sum[i-1])%MOD;  
        re%=MOD;  
    }  
    return (re+MOD)%MOD;  
}  
int main()  
{  
    int T,n,m;  
    Linear_Shaker();  
    for(cin>>T;T;T--)  
    {  
        scanf("%d%d",&n,&m);  
        printf("%d\n", Query(n,m) );  
    }  
}