「M-SOLUTIONS プロコンオープン」E問題
2019年6月1日開催のAtCoder M-SOLUTIONS プロコンオープンに参加した.結果はA, B 2完.このnoteでは時間内に解けなかったE問題の私なりの理解,解説を紹介する.
E問題 Product of Arithmetic Progression
問題はこちら.配点は600点.問題の概要をざっくり言うと,指示された等差数列の積を計算し,1000003で割ったあまりを出力する,という計算をたくさん(最大10^5個)する問題.内容は単純だが,実行時間制限: 2 sec以内に最大10^5個の等差数列の積を計算する必要があるため,工夫を要する.
mod p の世界での割り算
まず,mod p の世界での割り算を知る.これに関しては前回のnoteでも紹介したが,Kensuke Otsukiさんの大変詳しい記事 があるのでそちらを参照してほしい.
等差数列の積(mod p)の計算
ここからmod p (p = 1000003)での等差数列の積の計算について話す(内容は公式のpdf解説と同じ).求める等差数列の積は,
x * (x+d) * (x+2d) * (x+3d) * ... * (x+(n-1)d)
である.愚直に,一回掛け算をするごとにpで割ったあまりを計算することでも答えは出るが,制限時間内にはおさまらない.そこで,d ≠0のとき,各項をdで割り,最初にd^nをかけることで
≡ d^n * (x/d) * (x/d + 1) * (x/d + 2) * (x/d +3) * ... * (x/d + n-1) (mod p)
となる.等差数列の積の部分
(x/d) * (x/d + 1) * (x/d + 2) * (x/d +3) * ... * (x/d + n-1) (mod p)
は,公差 1 の等差数列の積になるので,
≡ (x/d + n - 1)! / (x/d - 1)! (mod p)
となる.分母・分子は,整数N ! であるから,あらかじめ必要な大きさのNまでのN! を計算しておくことで,Q個の等差数列の積の計算で共通化するため,計算時間が短縮できる.さらに,mod 1000003 なので,x/d + n - 1 が1000003以上の時は,即座に(x/d + n - 1)! ≡ 0 つまり,求める等差数列の積 ≡ 0 とできる.(つまり 1 ! から1000003! (mod p)まで計算しておけばよい)
d = 0 のときは x^n (mod p) を計算する.
C++での実装は次のようになる.(#include は省略)
#define rep(i,n) for(int i=0;i<n;i++)
using namespace std;
long long factorial(long long i){ // 階乗
if(i > 0)
return i * factorial(i-1);
else
return 1;
}
// a^n mod
long long modpow(long long a, long long n, long long mod) {
long long res = 1;
while (n>0) {
if(n & 1) res = res * a % mod;
a = a * a % mod;
n >>= 1;
}
return res;
}
// a^(-1) mod
long long modinv(long long a, long long mod) {
return modpow(a, mod-2, mod);
}
int main() {
long long Q;
const long long mod = 1000003;
cin >> Q;
// あらかじめ必要な階乗(mod)を計算しておく
vector<long long> factmod;
factmod.push_back(1); // 0! = 1
long long tmp = 1;
for(int i=1; i<mod+10; i++) { // 1000003!(mod 1,000,003) まで計算しておく ちょっと多めに
tmp = (tmp*i)%mod;
factmod.push_back(tmp);
}
long long ans = 1;
rep(i, Q) {
long long x, d, n;
cin >> x >> d >> n;
if(d!=0) {
// 求める積は d^n * (x/d) * (x/d + 1) * ... * (x/d + (n-1))
// つまり d^n * (x/d + (n-1))! / (x/d - 1)!
long long xbyd = x * modinv(d, mod) %mod; // x/d (mod)
if(xbyd+n-1 > mod) {
// x/d+n-1 が modより大きい時 (x/d + (n-1))! がゼロになる.
ans = 0;
} else {
if(xbyd==0) xbyd += mod; // xbyd-1 が負になるのを防ぐ
ans = modpow(d, n, mod) * factmod[xbyd+n-1] %mod * modinv(factmod[(xbyd-1)%mod], mod) %mod;
}
} else { // d = 0
// 求める積は x^n
ans = modpow(x, n, mod);
}
cout << ans << endl;
}
}
公式pdf解説にない注意点としては,N>1000003のとき,N! = 0なので,1000003 ! まで計算しておけば十分.
この記事が気に入ったらサポートをしてみませんか?