【地址】http://www.spoj.pl/problems/PGCD/
【题目大意】
10组数据,给定(m,n),问gcd(i,j)=素数的对数(1<=i<=m,1<=j<=n) m,n<=10^7
【算法】
不妨设m
先筛出素数和miu函数。
然后暴力容斥。
for (int d=2;d<=m;d++) if (prime[d])
for (int j=1;j<=m/d;j++)
ans+=miu(j)*(m/d/j)*(n/d/j)
这样是 m lg m的
miu是容斥系数。
特别地miu[1]=1。
miu[i]=1 i为Square-Free-Number并且i有偶数个素因子
miu[i]=-1 i为Square-Free-Number并且i有奇数个素因子
miu[i]=0 other
这样会TLE。
那么看我们所求。
sigma{ miu[j]*( m/(d*j) )*n/( n/(d*j) ) } 其中d*j<=m && d为素数
尝试枚举d*j的值。令x=d*j,由于d为素数,令x=p1^k1 * p2^k2 * p3^k3 … 其中pi为素数。
那么现在看看会出什么结果。
1、
假如ki全为1。那么从里面拿出一个素数d,会剩下一个Square-Free-Number。
设x有t个素因子。
j含有的素因子个数必然为t-1。所以所有的j的miu值都一样。而取法有t种于是ans+= ( ((t-1)&1)?-1:1 ) * t * (m/x) * (n/x)
2、
假如ki有一个是2,其它全是1,那么只能拿出那个ki=2的,否则剩下的不是SFN,miu值=0,无须计算。
这样设拿走了那个ki=2的一个素因子以后,剩下t个素因子。那么显然 ans+=( (t&1)?-1:1 )*(m/x)*(n/x)
3、
其它的不可能拿走一个素数以后得到SFN了。于是miu值都=0,不必考虑。
于是。。。枚举x。然后按照上面的情况分别搞就行了。
实际上都是P[x]*(m/x)*(n/x)的形式。由于上面的分类是针对x这个数值的,所以可以通过筛选法把P[x]预处理出来。
然后for (int x=2;x<=m;x++) ans+=P[x]*(m/x)*(n/x)就行了。
但是这样还不够快。
注意到m/x和n/x都只有sqrt级别个结果。于是随着x的递增,(m/x)*(n/x)只有sqrt(m)+sqrt(n)级别个结果。
通过对P[]建立前缀和数组,可以迅速知道区间的P[]之和。那么对于(m/x)*(n/x)相同的部分,直接跳过一整段,
ans+=(S_P[r]-S_P[l-1])*(m/x)*(n/x) 其中(m/l) * (n/l) == (m/r) * (n/r) 。
【时间复杂度】
一开始的筛O(n)的
后面每个询问sqrt(n)的。
【其它】
出现了很诡异的事件- -,我用O(n)的筛,45S。
用O(n lg n)的筛,31S。
不明真相中。
【CODE】
#include
#include
#include
#include
#include
using namespace std;
const int N=(int)1e7;
int F[N+5];
int cnt_1[N+5];
int cnt_2[N+5];
int S[N+5];
void Get_Prime(){
S[2]=0;
for (int i=2;i<=N;i++){
if (cnt_1[i]==0){
F[i]=i;
cnt_1[i]=1;
//这个筛是n lg n的- -,但不知道为什么,居然比线性筛快。求指导。。。
/*
for (int j=i+i,c=2;j<=N;j+=i,c++){
cnt_1[j]++;
if (c==i){
cnt_2[j]++;
if (cnt_2[j]==1 && j/i/i%i==0) cnt_2[j]++;
c=0;
}
}
*/
}
for (int j=2,k=i+i;j<=F[i];j++,k+=i){
if (k>N) break;
F[k]=j;
cnt_2[k]=cnt_2[i]+(j==F[i]); //计算k有多少个i^2形式的因子
cnt_1[k]=cnt_1[i]+1; //计算k有多少个因子。
}
int t;
int &c1=cnt_1[i];
int &c2=cnt_2[i];
if (c2==0)
if (c1&1)
t=c1;
else
t=-c1;
else
if (c2==1) t=( (c1&1)?1:-1 );
else t=0;
S[i+1]=S[i]+t;
}
}
int main(){
Get_Prime();
int Tc,m,n,d1,d2,Next_1,Next_2,Next;
scanf("%d",&Tc);
while (Tc–){
scanf("%d%d",&m,&n);
if (m>n) swap(m,n);
long long ans=0;
for (int i=2;i<=m;){
d1=m/i;
d2=n/i;
Next_1=m/d1+1; //Next_1>i && m/Next_1>m/i
Next_2=n/d2+1; //they are the same
Next=Next_1
ans+=(long long)(S[Next]-S[i])*d1*d2;
i=Next;
}
printf("%lldn",ans);
}
}