結果

問題 No.563 超高速一人かるた large
コンテスト
ユーザー koyumeishi
提出日時 2017-06-18 15:15:30
言語 C++14
(gcc 13.3.0 + boost 1.87.0)
結果
AC  
実行時間 279 ms / 3,000 ms
コード長 11,299 bytes
コンパイル時間 2,354 ms
コンパイル使用メモリ 145,288 KB
実行使用メモリ 33,624 KB
最終ジャッジ日時 2024-10-01 20:00:32
合計ジャッジ時間 5,250 ms
ジャッジサーバーID
（参考情報） 		judge3 / judge1
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
sample AC * 3
other AC * 18
権限があれば一括ダウンロードができます

ソースコード

diff # 

/*
O(n^2 log n) dp + NTT
*/

#include <iostream>
#include <string>
#include <vector>
#include <cassert>
#include <random>
#include <algorithm>
#include <set>
#include <map>
#include <chrono>
using namespace std;

void timer(){
  static chrono::steady_clock::time_point start = chrono::steady_clock::now();
  auto now = chrono::steady_clock::now();
  cerr << "elapsed time : " << chrono::duration_cast<chrono::milliseconds>(now-start).count() << "[ms]" << endl;
}

const long long mod = 1000000007;

long long ans[10001];
long long fact[10001];
long long fact_inv[10001];

long long mod_pow(long long x, long long y, long long mod){  //(x^y) % mod
  if(x==0 && y!=0) return 0;
  long long ret=1LL;
  while(y>0LL){
    if(y&1LL) ret = (ret * x) % mod;
    x = (x*x) % mod;
    y >>= 1LL;
  }
  return ret;
}

long long extgcd(long long a, long long b, long long &x, long long &y){
  long long d=a;
  if(b!=0){
    d = extgcd(b, a%b, y, x);
    y -= (a/b) * x;
  }else{
    x = 1;
    y = 0;
  }
  return d;
}

long long mod_inverse(long long a, long long m){
  long long x,y;
  extgcd(a,m,x,y);
  return (m+x%m)%m;
}
template<typename T = long long>
class Number_Theoretic_Transform {
  // return the vector of F[t] or f[x] where
  // F[t] = sum of { f[x] * exp(-j * theta * t * x) } in x = 0 .. N-1 (FFT)
  // f(x) = 1/N * sum of { F[t] * exp(j * theta * t * x) } in t = 0 .. N-1 (inverse FFT)
  // where theta = 2*PI / N
  // N == 2^k
  // use the rotater as (primitive-root of mod) ^ t in NTT, which is used as exp(j*theta)^t in FFT

  //事前に計算した 単位回転子 rotater (MOD mod 上での位数が2^kとなる数) を 引数に与える。
  //逆変換のときには rotater^-1 (MOD mod) を rotaterに与える

  vector< T > do_NTT(vector< T > A, bool inverse){
    int N = A.size();

    //bit reverse
    for(int bit=0; bit<N; bit++){
      int rbit = 0;
      for(int i=0, tmp_bit = bit; i<k-1; i++, rbit <<= 1, tmp_bit >>=  1){
        rbit |= tmp_bit & 1;
      }
      rbit >>= 1;
      if(bit < rbit){
        swap(A[bit], A[rbit]);
      }
    }


    int dist = 1;
    vector<T>& theta = (inverse?inv_theta_v:theta_v);

    T t = k-1;
    T half = theta[k-1];  //半回転

    for(int level = 1; level<k; level++){
      T W_n = theta[t]; //rotater ^ theta (MOD mod)
      T W = 1;              //rotater
      for(int x=0; x < (1<<(level-1)); x++){
        for(int y=x; y+dist < N; y += (dist<<1)){
          T tmp = A[y+dist]*W;
          if(tmp >= mod) tmp %= mod;

          A[y+dist] = A[y] + (tmp*half) % mod;
          if(A[y+dist] >= mod) A[y+dist] %= mod;

          A[y] = A[y] + tmp;
          if(A[y] >= mod) A[y] %= mod;
        }
        W = W*W_n;
        if(W>=mod) W%=mod;
      }
      dist <<= 1;
      t -= 1;
    }

    if(inverse){
      for(int i=0; i<N; i++){
        A[i] = z * A[i];
        if(A[i] >= mod) A[i] %= mod;
      }
    }
    

    return A;
  }

public:
  const T mod;
  const T rotater;
  const T inv_rotater;
  const T k;
  vector<T> theta_v;
  vector<T> inv_theta_v;
  const T z;

  Number_Theoretic_Transform(T mod_, T rotater_, T k_) : 
    mod(mod_), 
    rotater(rotater_), 
    k(k_), 
    inv_rotater(mod_inverse(rotater_, mod)),
    z(mod_inverse(1<<(k-1), mod)) // 1/Nを掛けるところなので N^-1 MOD modを掛けたいところだけど何故か (N/2)^-1 で上手く行く……
  {

    theta_v = vector<T>(k+1,1);
    theta_v[0] = rotater;
    for(int i=1; i<=k; i++){
      theta_v[i] = theta_v[i-1] * theta_v[i-1];
      if(theta_v[i] >= mod) theta_v[i] %= mod;
    }

    inv_theta_v = vector<T>(k+1,1);
    inv_theta_v[0] = inv_rotater;
    for(int i=1; i<=k; i++){
      inv_theta_v[i] = inv_theta_v[i-1] * inv_theta_v[i-1];
      if(inv_theta_v[i] >= mod) inv_theta_v[i] %= mod;
    }
    
    

  };

  vector< T > NTT(vector< T > A){
    return do_NTT(A, false);
  }

  vector< T > INTT(vector< T > A){
    return do_NTT(A, true);
  }

  // vector<double> C | C[i] = ΣA[i]B[size-1-j]
  vector<T> convolution(vector<T> &A, vector<T> &B){
    int n = A.size();

    assert(A.size() == B.size());
    assert( n == (1<<k) );  //Nは2^kである必要がある(事前にresize)

    auto A_NTT = NTT(A);
    auto B_NTT = NTT(B);

    for(int i=0; i<n; i++){
      A_NTT[i] *= B_NTT[i];
      if(A_NTT[i] >= mod) A_NTT[i] %= mod;
    }
    return INTT(A_NTT);
  }
};


long long gcd(long long a, long long b){
  if(b==0) return a;
  return gcd(b, a%b);
}

long long lcm(long long a, long long b){
  if(a<b) swap(a,b);
  if(b==1) return a;
  return a * (b/gcd(a,b));
}

// Z % Yi = Xi であるようなZを求める。Garnerのアルゴリズム O(N^2)
// 参考 http://techtipshoge.blogspot.jp/2015/02/blog-post_15.html
// http://yukicoder.me/problems/448
long long Chinese_Remainder_Theorem_Garner(vector<long long> x, vector<long long> y, long long MOD){
  int N = x.size();

  bool valid = true;
  //前処理
  //gcd(Yi,Yj) == 1 (i!=j) でなくてはならないので、
  //共通の因数 g = gcd(Yi,Yj) を見つけたら片側に寄せてしまう
  for(int i=0; i<N; i++){
    for(int j=i+1; j<N; j++){
      if(i == j) continue;
      long long g = gcd(y[i], y[j]);

      if( x[i]%g != x[j]%g ) valid = false; //解が存在しない

      if(g != 1){
        y[i] /= g; y[j] /= g;
        long long g_ = gcd(y[i], g);
        while(g_ != 1){
          y[i] *= g_;
          g /= g_;
          g_ = gcd(y[i], g);
        }
        y[j] *= g;

        x[i] %= y[i];
        x[j] %= y[j];
      }
    }
  }

  if(!valid){
    cerr << -1 << endl;
    abort();
    return 0;
  }

  //Garner's algorithm
  vector<long long> z(N);
  for(int i=0; i<N; i++){
    z[i] = x[i];
    for(int j=0; j<i; j++){
      z[i] = mod_inverse(y[j], y[i]) % y[i] * (z[i] - z[j]) % y[i];

      z[i] = (z[i]+y[i])%y[i];
    }
  }

  long long ans = 0;
  long long tmp = 1;
  for(int i=0; i<N; i++){
    ans = (ans + z[i] * tmp)%MOD;

    tmp = (tmp * y[i])%MOD;
  }

  return ans;
}


struct trie_node{
  map<char, int> child;
  int count;
};
vector<trie_node> trie;

vector<long long> s[200001];
vector<long long> w[200001];

vector<long long> convolution(vector<long long> a, vector<long long> b){
  int x = a.size();
  int y = b.size();
  if(x*y <= (x+y)*30 ){
    vector<long long> res(x+y-1);
    for(int i=0; i<x; i++){
      for(int j=0; j<y; j++){
        (res[i+j] += a[i] * b[j] % mod) %= mod;
      }
    }
    return res;
  }

  long long k = 5;
  while( (1<<k) <= x+y ) k++;
  k++;
  assert( k < 16 );

  a.resize(1<<k, 0);
  b.resize(1<<k, 0);

  const static vector<vector<vector<int>>> param = {
    {
    },
    {
    },
    {
      {2, 17011, 4989209},
      {2, 7027, 4937873},
      {2, 12941, 4925573},
    },
    {
      {3, 13003, 4827737},
      {3, 9521, 4813577},
      {3, 10079, 4767457},
    },
    {
      {4, 15859, 4899329},
      {4, 4679, 4895057},
      {4, 15749, 4866209},
    },
    {
      {5, 15817, 4994401},
      {5, 15227, 4988449},
      {5, 14947, 4917217},
    },
    {
      {6, 8117, 4989121},
      {6, 17123, 4974913},
      {6, 761, 4952257},
    },
    {
      {7, 14551, 4928513},
      {7, 9781, 4832257},
      {7, 4663, 4797953},
    },
    {
      {8, 4973, 4933889},
      {8, 10159, 4856321},
      {8, 16229, 4817921},
    },
    {
      {9, 10567, 4996609},
      {9, 12143, 4996097},
      {9, 13249, 4984321},
    },
    {
      {10, 631, 4995073},
      {10, 15727, 4979713},
      {10, 14029, 4953089},
    },
    {
      {11, 3041, 4933633},
      {11, 4201, 4909057},
      {11, 569, 93493249},
    },
    {
      {12, 16633, 4902913},
      {12, 5527, 4845569},
      {12, 1009, 97497089},
    },
    {
      {13, 13687, 48701441},
      {13, 709, 4866049},
      {13, 163, 4816897},
    },
    {
      {14, 24763, 99500033},
      {14, 26249, 99106817},
      {14, 14867, 98992129},
    },
    {
      {15, 9859, 98861057},
      {15, 6427, 97615873},
      {15, 11393, 97320961},
    },
  };

  vector<long long> P;
  vector<long long> R;

  vector<vector<long long>> tmp(3);
  for(int i=0; i<3; i++){
    R.push_back( param[k][i][1] );
    P.push_back( param[k][i][2] );

    assert( mod_pow(R[i], 1<<k, P[i]) == 1  );
    assert( mod_pow(R[i], 1<<k-1, P[i]) != 1  );

    Number_Theoretic_Transform<long long> NTT(P[i], R[i], k);
    tmp[i] = NTT.convolution(a,b);
  }

  vector<long long> res(x+y-1);
  for(int i=0; i<res.size(); i++){
    vector<long long> w = {tmp[0][i], tmp[1][i], tmp[2][i]};
    res[i] = Chinese_Remainder_Theorem_Garner(w,P, mod);
  }

  return res;
}

// dp[i][x] := i 頂点で合計長 x の組み合わせ
void dfs(int pos, long long d){
  if( trie[pos].count == 1 ){
    w[pos].resize(2);
    s[pos].resize(2);
    w[pos][0] = w[pos][1] = 1;
    s[pos][0] = 0; s[pos][1] = 1;
    return;
  }
  if( trie[pos].child.size() == 1 ){
    int nx = trie[pos].child.begin()->second;
    dfs(nx, d+1);
    swap( s[pos], s[nx] );
    swap( w[pos], w[nx] );
    return;
  }

  w[pos].resize(1);
  s[pos].resize(1);
  w[pos][0] = 1;
  s[pos][0] = 0;

  for(auto p : trie[pos].child){
    if( p.first == '}' ) continue;

    int nx = p.second;
    dfs(nx, 1);

    for(int i=1; i<w[pos].size(); i++){
      w[pos][i] = w[pos][i] * fact_inv[i] % mod;
      s[pos][i] = s[pos][i] * fact_inv[i] % mod;
    }
    for(int i=1; i<w[nx].size(); i++){
      w[nx][i] = w[nx][i] * fact_inv[i] % mod;
      s[nx][i] = s[nx][i] * fact_inv[i] % mod;
    }

    auto w_ = convolution(w[pos], w[nx]);
    auto s1 = convolution(s[pos], w[nx]);
    auto s2 = convolution(w[pos], s[nx]);
    vector<long long> s_(w_.size());
    for(int i=0; i<s1.size(); i++){
      s_[i] = s1[i] + s2[i];
      s_[i] %= mod;
    }

    // for(int i=0; i<w[pos].size(); i++){
    //   for(int j=0; j<w[nx].size(); j++){
    //     if(i+j >= w_.size()){
    //       w_.resize(i+j+1);
    //       s_.resize(i+j+1);
    //     }
    //     (w_[i+j] += w[pos][i] * w[nx][j] % mod) %= mod;
    //     (s_[i+j] += s[pos][i] * w[nx][j] % mod) %= mod;
    //     (s_[i+j] += w[pos][i] * s[nx][j] % mod) %= mod;
    //   }
    // }

    swap(w[pos], w_);
    swap(s[pos], s_);

    for(int i=1; i<w[pos].size(); i++){
      w[pos][i] = w[pos][i] * fact[i] % mod;
      s[pos][i] = s[pos][i] * fact[i] % mod;
    }
  }

  w[pos].resize( trie[pos].count + 1 );
  s[pos].resize( trie[pos].count + 1 );
  for(long long i=1; i<w[pos].size(); i++){
    (s[pos][i] += d*min<long long>(i, trie[pos].count-1)%mod * w[pos][i]) %= mod;
  }
}

int main(){
  timer();

  fact[0] = fact_inv[0] = 1;
  for(int i=1; i<10001; i++){
    fact[i] = fact[i-1] * i % mod;
    fact_inv[i] = mod_pow(fact[i], mod-2, mod);
  }

  int n;
  cin >> n;
  vector<string> s_(n);
  for(int i=0; i<n; i++){
    cin >> s_[i];
    s_[i] += '{';
  }
  s_.push_back("}");

  trie.push_back( trie_node{{}, 0} );
  for(int i=0; i<s_.size(); i++){
    int pos = 0;
    for(int j=0; j<s_[i].size(); j++){
      trie[pos].count++;
      if( trie[pos].child.count( s_[i][j] ) == 0 ){
        trie.push_back( trie_node{ {}, 0} );
        trie[pos].child[s_[i][j]] = trie.size()-1;
      }
      pos = trie[pos].child[s_[i][j]];
    }
    trie[pos].count++;
  }

  dfs(0, 0);

  for(int i=1; i<=n; i++){
    ans[i-1] = s[0][i];
  }

  for(int i=0; i<n; i++) cout << ans[i]%mod << endl;

  timer();
  return 0;
}
0