2154的加强版
可以发现 由莫比乌斯函数组成的积性函数 在 prime[j]|i 转移时 多出来的u(p*p*a)=0部分不记的 只需考虑F(x) 中的x 乘上了p对之前的影响
#include<cstdio>
#include<cstdlib>
#include<algorithm>
using namespace std;
typedef long long ll;
inline char nc()
{
static char buf[100000],*p1=buf,*p2=buf;
if (p1==p2) { p2=(p1=buf)+fread(buf,1,100000,stdin); if (p1==p2) return EOF; }
return *p1++;
}
inline void read(int &x)
{
char c=nc(),b=1;
for (;!(c>='0' && c<='9');c=nc()) if (c=='-') b=-1;
for (x=0;c>='0' && c<='9';x=x*10+c-'0',c=nc()); x*=b;
}
const int maxn=10000000,P=100000009;
int prime[maxn],vst[maxn],tot;
ll sum[maxn];
inline void Init(){
sum[1]=1;
for (int i=2;i<=maxn;i++)
{
if (!vst[i]){
prime[++tot]=i; sum[i]=(ll)i*(i-1)%P; sum[i]=(P-sum[i])%P;
}
for (int j=1;(ll)prime[j]*i<=maxn && j<=tot;j++)
{
vst[i*prime[j]]=1;
if (i%prime[j]==0){
sum[i*prime[j]]=sum[i]*prime[j]%P;
break;
}
sum[i*prime[j]]=sum[i]*sum[prime[j]]%P;
}
}
for (int i=1;i<=maxn;i++)
(sum[i]+=sum[i-1])%=P;
}
inline ll Sum(ll x,ll y){
return (x*(x+1)/2%P)*(y*(y+1)/2%P)%P;
}
inline void Solve(int n,int m){
if (n>m) swap(n,m);
ll ans=0; int last;
for (int i=1;i<=n;i=last+1)
{
last=min(n/(n/i),m/(m/i));
ans+=Sum(n/i,m/i)*(sum[last]-sum[i-1]);
ans=(ans%P+P)%P;
}
printf("%lld\n",ans);
}
int main()
{
int T,n,m;
Init();
read(T);
while (T--){
read(n); read(m);
Solve(n,m);
}
return 0;
}