Description
Input
一个正整数T表示数据组数
接下来T行 每行两个正整数 表示N、M
Output
T行 每行一个整数 表示第i组数据的结果
Sample Input
1
4 5
Sample Output
122
HINT
T <= 10000
N, M<=10000000
Solution
这题就是上一题（BZOJ2154）的升级版
我们接着推
[ans=sum_{d=1}^Ndcdot sum_{i=1}^{N}mu(i)cdot i^2cdot s(lfloor frac{N}{id}
floor,lfloor frac{M}{id}
floor)
]
把(id)提出来
[ans=sum_{T=1}^Ns(lfloor frac{N}{T}
floor,lfloor frac{M}{T}
floor)sum_{d|T}frac{T}{d}cdot d^2cdot mu(d)
]
第一个sigma整除分段，我们主要来考虑第二个sigma
设(H(n)=sum_{d|n}frac{n}{d}cdot d^2cdot mu(d))，(F(i)=i)，(G(i)=i^2mu(i))
易知：(F)为积性函数，而(G(ab)=(ab)^2mu(ab)=a^2mu(a)b^2mu(b))（(mu)也是积性函数）(=G(a)G(b))，所以(G)也是积性函数
那么它们卷起来的(H)就也会是积性函数
- (H(1)=1)
- (H(p,p is a prime)=p-p^2)（直接带进去）
- (H(p^a,p is a prime)=sum_{i=0}^ap^{a-i}cdot p^{2i}cdot mu(p^i))
当(i>1)时，根据mu的定义，(mu(p^i))为0，所以(H(p^a,p is a prime)=p^a-p^{a+1})
这也就是说在素数筛时，对于一个(i*prime[j])，如果(i\%prime[j]!=0)，说明(prime[j])是一个新的(i*prime[j])的质数约数，那么(H[i*prime[j]])就加上一个(H[prime[j]]) - (H(prod_{i=1} P_i^{a_i} )=prod_{i=1} H(P_i^{a_i}))
当(i\%pime[j]=0)，也就是(i*prime[j])中有一质数约数的指数大于1,这个时候因为第三点，它的(mu)是0，不会有额外的贡献，只是(H(i*prime[j]))会乘一个(prime[j])
(H(P_1^{a_1+1}prod_{i=2}P_i^{a_i})=H(P_1^{a_1+1})H(prod_{i=2}P_i^{a_i}))（积性函数）
( =prod_{i=2}(P_i^{a_i}-P_i^{a_i+1})(P_1^{a_1+1}-P_1^{a_1+2}))
( =prod_{i=2}(P_i^{a_i}-P_i^{a_i+1})(P_1^{a_1}-P_1^{a_1+1})P_1)
( =prod_{i=1}(P_i^{a_i}-P_i^{a_i+1})P_1)
( =H(prod_{i=1}P_i^{a_i})P_1)
后面素数筛搞完后，整除分块，做完
哦，还有个天大的坑。。。
mod的数是1e8+9，它也不是质数。。。
#include<bits/stdc++.h>
#define ll long long
const int Mod=1e8+9,MAXN=10000000+10;
ll s[MAXN],H[MAXN];
int prime[MAXN],cnt,vis[MAXN];
template<typename T> inline void read(T &x)
{
T data=0,w=1;
char ch=0;
while(ch!='-'&&(ch<'0'||ch>'9'))ch=getchar();
if(ch=='-')w=-1,ch=getchar();
while(ch>='0'&&ch<='9')data=((T)data<<3)+((T)data<<1)+(ch^'0'),ch=getchar();
x=data*w;
}
template<typename T> inline void write(T x,char c=' ')
{
if(x<0)putchar('-'),x=-x;
if(x>9)write(x/10);
putchar(x%10+'0');
if(c!=' ')putchar(c);
}
template<typename T> inline void chkmin(T &x,T y){x=(y<x?y:x);}
template<typename T> inline void chkmax(T &x,T y){x=(y>x?y:x);}
template<typename T> inline T min(T x,T y){return x<y?x:y;}
template<typename T> inline T max(T x,T y){return x>y?x:y;}
inline void init()
{
memset(vis,1,sizeof(vis));
vis[0]=vis[1]=0;
H[1]=1;
for(register ll i=2;i<MAXN;++i)
{
if(vis[i])
{
prime[++cnt]=i;
H[i]=(i-(ll)i*i)%Mod;
}
for(register ll j=1;j<=cnt&&i*prime[j]<MAXN;++j)
{
vis[i*prime[j]]=0;
if(i%prime[j])H[i*prime[j]]=H[i]*H[prime[j]]%Mod;
else
{
H[i*prime[j]]=H[i]*(ll)prime[j]%Mod;
break;
}
}
}
for(register ll i=1;i<MAXN;++i)s[i]=(s[i-1]+H[i])%Mod;
}
inline ll S(ll x,ll y)
{
return ((x+1)*x/2)%Mod*(((y+1)*y/2)%Mod)%Mod;
}
inline ll solve(ll N,ll M)
{
ll res=0;
if(N>M)std::swap(N,M);
for(register ll i=1;;)
{
if(i>N)break;
ll j=min(N/(N/i),M/(M/i));
(res+=S(N/i,M/i)*(s[j]-s[i-1])%Mod)%=Mod;
i=j+1;
}
return (res+Mod)%Mod;
}
int main()
{
init();
int T;
read(T);
while(T--)
{
ll N,M;
read(N);read(M);
write(solve(N,M),'
');
}
return 0;
}