记录编号	98410	评测结果	AAAAAA
题目名称	[SPOJ 422] 转置更好玩	最终得分	100
用户昵称	cstdio	是否通过	通过
代码语言	C++	运行时间	5.368 s
提交时间	2014-04-22 22:58:55	内存使用	0.32 MiB

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<vector>
#include<cstring>
#define MOD 1000003
using namespace std;
typedef long long ll;
const int SIZEN=1010,SIZEP=60;
bool flag[SIZEN];
int prime[SIZEN]={0},pn=0;
int factr[SIZEP]={0};//底数
int facti[SIZEP]={0};//指数
int factv[SIZEP]={0};//值
int fn;
void getprime(int n){
	for(int i=2;i<=n;i++){
		if(!flag[i]){
			prime[pn++]=i;
			for(int j=i*i;j<=n;j+=i){
				flag[j]=true;
			}
		}
	}
}
int gcd(int a,int b){
	return !b?a:gcd(b,a%b);
}
void extend_gcd(ll a,ll b,ll &x,ll &y){//ax+by=(a,b)
	if(!b){x=1,y=0;}
	else{
		extend_gcd(b,a%b,y,x);
		y-=a/b*x;
	}
}
int pow2[33]={0};
void initpow(void){
	pow2[0]=2;
	for(int i=1;i<32;i++) pow2[i]=1ll*pow2[i-1]*pow2[i-1]%MOD;
}
int quickpow2(int n){
	int ans=1;
	for(int i=0;n>0;n>>=1,i++){
		if(n&1) ans=1ll*ans*pow2[i]%MOD;
	}
	return ans;
}
void getfact(int n){//做分解
	fn=0;
	for(int i=0;i<pn;i++){
		if(n%prime[i]==0){
			factr[fn]=prime[i];
			facti[fn]=0;
			factv[fn]=1;
			while(n%prime[i]==0){
				facti[fn]++;
				factv[fn]*=prime[i];
				n/=prime[i];
			}
			fn++;
		}
	}
	if(n>1){//本身是一个大素数
		factr[fn]=n;
		facti[fn]=1;
		factv[fn]=n;
		fn++;
	}
}
void DFS(int &ans,int g,int phi,int dep){//gcd值是g的情况
	//则情况数量为phi((a+b)/g),这里的phi是Euler函数
	//但变量phi只储存phi((a+b)/g)的表达式中只涉及少于dep的素因子的部分的结果
	//例如,如果最终的Euler函数值是((2-1)*(2^3))*((5-1)*(5^6))*((7-1)*(7^10)),那么phi有可能仅仅是((2-1)*(2^3))
	if(dep==fn){
		//如果这里情况中的某个循环一次走x步,那么(x,a+b)=g
		//可以看到,对每个x,合法的"保持不变"染色个数是2^g,一共phi*2^g个
		ans=(1ll*phi*quickpow2(g)+ans)%MOD;
		return;
	}
	int i=0,p=factv[dep];
	while(i<=facti[dep]){
		DFS(ans,factv[dep]/p*g,phi*(p-p/factr[dep]),dep+1);
		//factv[dep]/p*g是把当前的g值乘上当前枚举的素因子的一个幂
		//具体地说,是乘上factr[dep]^i
		//于是涉及factr[dep]的phi值应当是phi*(factr[dep]-1)*factr[dep]^(facti[dep]-i-1)(根据Euler函数的积性)
		i++,p/=factr[dep];
	}
}
int calc(int a,int b){
	if(!a||!b) return 0;
	int g=gcd(a,a+b);
	getfact((a+b)/g);
	int ans=0;
	DFS(ans,g,1,0);//一"小块"的长度是g,因此合法的置换一次都要走g的倍数
	ll x,y;
	extend_gcd((a+b)/g,MOD,x,y);//此处求逆元,因为要除以总的置换数量(a+b)/g
	x%=MOD;
	ans=((quickpow2(a+b)-x*ans)%MOD+MOD)%MOD;
	return ans;
}
int main(){
	freopen("transp2.in","r",stdin);
	freopen("transp2.out","w",stdout);
	initpow();
	getprime(SIZEN-1);
	int T;
	scanf("%d",&T);
	while(T--){
		int a,b;
		scanf("%d%d",&a,&b);
		printf("%d\n",calc(a,b));
	}
	return 0;
}