kmjp's blog

競技プログラミング参加記です

Codeforces #235 Div2. E. Olympic Games

本番には解けなかった問題。
http://codeforces.com/contest/401/problem/E

問題

(0,0)-(N,M)の格子点(N+1)*(M+1)個がある。
このうち、以下の条件を満たす格子点の対の数を答えよ。

  • 2点を結ぶ線分上に他の格子点がない。
  • 2点のユークリッド距離はL以上R以下である。

解法

まず、上下または左右同士の対が条件を満たすのは、L==1の時だけである。
この際上下の対が(N+1)M組、左右の対がN(M+1)組できる。

残りは斜めに位置する対で条件を満たすものを考える。
ある点(x,y)に対し、その右上の点(x+dx,y+dy)が条件を満たすのは以下の条件を満たすときである。

  •  L \le \sqrt{dx^2+dy^2} \le R
  • gcd(dx,dy)=1

また、(x,y)-(x+dx,y+dy)が(0,0)-(N,M)内に入るような(x,y)の組み合わせは(N+1-dx)*(M+1-dy)通りある。
同様に(x,y)-(x-dx,y+dy)の組み合わせも同数ある。

よって上記条件を満たす(dx,dy)について2*(N+1-dx)*(M+1-dy)の総和を取ればそれが答えとなる。
dyは総当たりで1~min(M,R)までチェックし、条件を満たすdxに対し2*(N+1-dx)の総和を取って(M+1-dy)を足しこんでいけばよい。

dxの条件は2つ。1つは L \le \sqrt{dx^2+dy^2} \le Rなので \sqrt{L^2-dy^2} \le dx \le \sqrt{R^2-dy^2}

問題はgcd(dx,dy)=1の条件である。
これはdyを素因数分解して包除原理を用いればよい。
dxのうち、dyの素因数を偶数個掛けて得られるものの倍数は足し、奇数個かけて得られるもの倍数は引く。
例えばdy=30=2*3*5なら:

  • dxが1の倍数の時の2*(N+1-dx)の総和を足す
  • dxが2・3・5の倍数の時の2*(N+1-dx)の総和を引く
  • dxが2*3、3*5、2*5の倍数の時の2*(N+1-dx)の総和を足す
  • dxが2*3*5の倍数の時の2*(N+1-dx)の総和を引く

なお、dxがある数Pの倍数の時の2*(N+1-dx)の総和は前処理しておくとよい。
これ計算量どの位だろう。dyの総当たりがO(M)、素因数分解がO(√M)、包除原理がO(2^素因数の数)なんだけど素因数分解と包除原理どっちが上かな?

ll N,M,L,R,P;
ll XL[160001],XR[160001];
vector<ll> S[100001];

vector<ll> enumdiv(ll n) {
	vector<ll> V;
	for(ll i=2;i*i<=n;i++) {
		if(n%i==0) V.push_back(i);
		while(n%i==0) n/=i;
	}
	if(n>1) V.push_back(n);
	return V;
}

void solve() {
	int f,i,j,k,l,x,y;
	cin>>N>>M>>L>>R>>P;
	
	ll ret=0;
	// right&&down
	if(L==1) ret=(N*(M+1) + M*(N+1))%P;
	S[1].push_back(0);
	for(i=1;i<=N;i++) S[1].push_back((N+1-i)%P);
	for(i=2;i<=M;i++) {
		S[i].push_back(0);
		for(j=1;j*i<=N;j++) S[i].push_back((S[i].back()+S[1][j*i])%P);
	}
	for(i=1;i<=N;i++) S[1][i]=(S[1][i]+S[1][i-1])%P;
	for(y=R-1;y>0;y--) {
		XL[y]=XL[y+1];
		XR[y]=XR[y+1];
		while(XL[y]*XL[y]+y*(ll)y<L*L) XL[y]++;
		while(XR[y]*XR[y]+y*(ll)y<=R*R) XR[y]++;
		XR[y]--;
	}
	for(y=1;y<=min(M,R);y++) {
		ll line=0;
		if(XL[y]==0) XL[y]++;
		if(XL[y]>N) continue;
		if(XR[y]>N) XR[y]=N;
		
		vector<ll> D = enumdiv(y);
		FOR(i,1<<D.size()) {
			k=1;
			l=(__builtin_popcount(i)%2)?-1:1;
			FOR(j,D.size()) if(i&(1<<j)) k*=D[j];
			line += l*(S[k][XR[y]/k]-S[k][(XL[y]-1)/k]);
			line %= P;
			if(line<0) line+=P;
		}
		
		ret+=(line*2)%P*(M+1-y);
		ret%=P;
		if(ret<0) ret+=P;
		
	}
	cout<<ret<<endl;
}

まとめ

本番、gcd(dx,dy)==1の問題が解消できず解答できなかった。
前処理と素因数を使った包除原理か…覚えておこう。