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 ! まで計算しておけば十分．