bzoj 4454: C Language Practice 最大公约数

       最大公约数的O(N)预处理O(1)求法。

       首先令m=N^0.5，那么我们可以通过记忆化搜索的形式得到所有的gcd(i,j)，其中i,j<=m，这样状态数和时间都是O(m*m)=O(N)的。

       那么考虑求两个不超过N的数x，y的gcd；我们一定可以把x分解成三个数a,b,c，每个数都满足要么<=m，要么是质数，简单的证明如下：

       假设存在关于某个x的a不满足条件，那么显然有a为合数且a>m。那么假设a的一种分解方式为a=uv（不妨令u>v），由于a>m，而abc=x<N=m^2，那么bc<m，因此x也存在另一种分解方式u,v,bc，如果此时u<=m或u为质数，分解完毕；否则可以对u不断进行类似的操作指导满足条件。

       那么考虑对a,b,c和y分别求gcd然后再并起来；假设对a操作。通过预处理我们可以知道a是否是质数，如果是质数，那么判断是否有y mod a=0，满足就让y/=a,ans*=a，否则一定有a<=m，那么由gcd(m,n)=gcd(n,m mod n)，可以得到gcd(y,a)=gcd(a,y mod a)，显然a和y mod a都<=m，那么这部分预处理已经得到了。

       像这样就可以做到O(N)预处理O(1)。注意这道题目卡内存，只能精确取m=N^0.5，实践中可以把m取得大一点。

AC代码如下：

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;

int n,m,cnt,a[2005],b[2005],f[1000005][3],g[1005][1005],p[1000005],c[100005];
int read(){
	int x=0; char ch=getchar();
	while (ch<'0' || ch>'9') ch=getchar();
	while (ch>='0' && ch<='9'){ x=x*10+ch-'0'; ch=getchar(); }
	return x;
}
int get_g(int x,int y){
	return (x && y)?g[y][x%y]:x|y;
}
int gcd(int x,int y){
	if (x<=1000 && y<=1000) return g[x][y];
	if (!x || !y) return x|y;
	int d=1,i;
	for (i=0; i<3; i++) if (f[x][i]>1){
		int t=f[x][i];
		if (p[t]==t){
			if (!(y%t)){
				y/=t; d*=t;
			}
		} else{
			int z=g[t][y%t];
			y/=z; d*=z;
		}
	}
	return d;
}
int main(){
	int i,j,cas=read();
	for (i=0; i<=1000; i++)
		for (j=0; j<=i; j++) g[i][j]=g[j][i]=get_g(i,j);
	p[1]=1;
	for (i=2; i<=1000000; i++){
		if (!p[i]){
			p[i]=i; c[++cnt]=i;
		}
		for (j=1; j<=cnt && i*c[j]<=1000000; j++){
			p[i*c[j]]=c[j];
			if (!(i%c[j])) break;
		}
	}
	f[1][0]=f[1][1]=f[1][2]=1;
	for (i=2; i<=1000000; i++){
		memcpy(f[i],f[i/p[i]],sizeof(f[i]));
		if (f[i][0]*p[i]<=1000) f[i][0]*=p[i];
		else if (f[i][1]*p[i]<=1000) f[i][1]*=p[i];
		else f[i][2]*=p[i];
	}
	while (cas--){
		m=read(); n=read();
		for (i=0; i<m; i++) a[i]=read();
		for (i=0; i<n; i++) b[i]=read();
		unsigned int ans=0;
		for (i=0; i<m; i++) for (j=0; j<n; j++){
			ans+=gcd(a[i],b[j])^i^j;
		}
		printf("%u\n",ans);
	}
	return 0;
}

by lych

2016.3.30