結果

問題 No.2120 場合の数の下８桁
コンテスト
ユーザー emthrm
提出日時 2022-11-04 21:37:17
言語 C++17
(gcc 13.3.0 + boost 1.87.0)
結果
AC  
実行時間 1,139 ms / 2,000 ms
コード長 3,421 bytes
コンパイル時間 2,209 ms
コンパイル使用メモリ 201,616 KB
最終ジャッジ日時 2025-02-08 17:23:57
ジャッジサーバーID
（参考情報） 		judge5 / judge3
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
other AC * 20
権限があれば一括ダウンロードができます

ソースコード

diff # 

// https://ferin-tech.hatenablog.com/entry/2018/01/17/010829
#define __USE_MINGW_ANSI_STDIO 0
#include <bits/stdc++.h>

using namespace std;
typedef long long ll;
#define int ll
typedef vector<int> VI;
typedef vector<VI> VVI;
typedef pair<int, int> PII;

#define FOR(i, a, n) for (ll i = (ll)a; i < (ll)n; ++i)
#define REP(i, n) FOR(i, 0, n)
#define ALL(x) x.begin(), x.end()
#define IN(a, b, x) (a<=x&&x<b)
#define PB push_back

const ll LLINF = (1LL<<60);
const int INF = (1LL<<30);
const int MOD = 1000000007;

template <typename T> T &chmin(T &a, const T &b) { return a = min(a, b); }
template <typename T> T &chmax(T &a, const T &b) { return a = max(a, b); }
template<class S,class T>
ostream &operator <<(ostream& out,const pair<S,T>& a){
  out<<'('<<a.first<<','<<a.second<<')';
  return out;
}

int dx[] = {0, 1, 0, -1}, dy[] = {1, 0, -1, 0};

// ax + by = gcd(a, b) となる {x, y, gcd(a, b)} を返す
// O(log(min(a, b)))
ll extgcd(ll a, ll b, ll &x, ll &y) {
  ll g = a; x = 1, y = 0;
  if(b != 0) g = extgcd(b, a%b, y, x), y -= (a/b) * x;
  return g;
}

// a^-1 mod n を返す　存在しなければ-1
// O(log(n))
ll inv(ll a, ll n) {
  ll s, t;
  extgcd(a, n, s, t);
  return (n+s) % n;
}

// 二分累乗法 x^e % mod O(log(e))
ll binpow(ll x, ll e, ll mod) {
  ll a = 1, p = x;
  while(e > 0) {
    if(e%2 == 0) {p = (p*p) % mod; e /= 2;}
    else {a = (a*p) % mod; e--;}
  }
  return a % mod;
}

// x = a1 mod m1, x2 = a2 mod m2 を解く オーバーフローには注意
// O(log(min(m1, m2)))
pair<ll, ll> crt(ll a1, ll a2, ll m1, ll m2) {
  auto normal = [](ll x, ll m) { return x>=-x ? x%m : m-(-x)%m; };
  auto modmul = [&normal](ll a, ll b, ll m) { return normal(a, m)*normal(b, m)%m; };
  ll k1, k2;
  ll g = extgcd(m1, m2, k1, k2);
  if(normal(a1, g) != normal(a2, g)) return {-1, -1};
  ll l = m1 / g * m2;
  ll x = a1 + modmul(modmul((a2-a1)/g, k1, l), m1, l);
  return {x, l};
}

pair<ll, ll> crt(vector<ll> a, vector<ll> m) {
  ll mod = 1, ans = 0;
  int n = a.size();
  REP(i, n) {
    tie(ans, mod) = crt(ans, a[i], mod, m[i]);
    if(ans == -1) return {-1, -1};
  }
  return {ans, mod};
}

ll fact[1000010], ifact[1000010];
void makeFac(ll p, ll q) {
  ll pr = 1;
  REP(i, q) pr *= p;
  fact[0] = ifact[0] = 1;
  FOR(i, 1, pr+1) {
    if(i%p == 0) {
      fact[i] = fact[i-1];
    } else {
      fact[i] = fact[i-1] * i % pr;
    }
    ifact[i] = inv(fact[i], pr);
  }
}

ll C(ll n, ll r, ll p, ll q) {
  if(n < 0 || r < 0 || r > n) return 0;
  // pr = p^q
  int pr = 1;
  REP(i, q) pr *= p;

  ll z = n-r;
  int e0 = 0;
  for(ll u=n/p;u;u/=p) e0 += u;
  for(ll u=r/p;u;u/=p) e0 -= u;
  for(ll u=z/p;u;u/=p) e0 -= u;
  int em = 0;
  for(ll u=n/pr;u;u/=p) em += u;
  for(ll u=r/pr;u;u/=p) em -= u;
  for(ll u=z/pr;u;u/=p) em -= u;

  ll ret = 1;
  while(n > 0) {
    ret = ret*fact[n%pr]%pr*ifact[r%pr]%pr*ifact[z%pr]%pr;
    n /= p; r /= p; z /= p;
  }
  (ret *= binpow(p, e0, pr)) %= pr;
  if(!(p==2 && q >= 3) && em%2) ret = (pr-ret) % pr;
  return ret;
}

ll func(ll n, ll r, ll mod) {
  ll x = mod;
  vector<ll> a, m;
  FOR(i, 2, mod+1) if(x%i == 0) {
    ll cnt=0, pr=1;
    while(x%i==0) x/=i, cnt++, pr*=i;
    makeFac(i, cnt);
    a.PB(C(n, r, i, cnt));
    m.PB(pr);
  }

  return crt(a, m).first;
}

signed main(void)
{
  ll n, m;
  cin >> n >> m;
  const string ans = to_string(func(n, m, 100000000));
  cout << string(8 - ans.length(), '0') << ans << '\n';

  return 0;
}
0