結果
| 問題 |
No.2120 場合の数の下8桁
|
| コンテスト | |
| ユーザー |
kwm_t
|
| 提出日時 | 2022-11-04 22:29:27 |
| 言語 | C++17(gcc12) (gcc 12.3.0 + boost 1.87.0) |
| 結果 |
AC
|
| 実行時間 | 10 ms / 2,000 ms |
| コード長 | 8,228 bytes |
| コンパイル時間 | 5,621 ms |
| コンパイル使用メモリ | 125,268 KB |
| 最終ジャッジ日時 | 2025-02-08 18:04:08 |
|
ジャッジサーバーID (参考情報) |
judge5 / judge1 |
(要ログイン)
| ファイルパターン | 結果 |
|---|---|
| other | AC * 20 |
ソースコード
/**
* date : 2021-04-27 12:43:05
*/
#define NDEBUG
#include <cassert>
#include <tuple>
#include <utility>
#include <vector>
#include <string>
//
using namespace std;
constexpr pair<long long, long long> inv_gcd(long long a, long long b) {
a = a % b;
if (a < 0) a += b;
if (a == 0) return { b, 0 };
long long s = b, t = a;
long long m0 = 0, m1 = 1;
while (t) {
long long u = s / t;
s -= t * u;
m0 -= m1 * u;
auto tmp = s;
s = t;
t = tmp;
tmp = m0;
m0 = m1;
m1 = tmp;
}
if (m0 < 0) m0 += b / s;
return { s, m0 };
}
vector<long long> pre_im(vector<long long> m) {
vector<long long> ims;
long long m0 = 1;
int n = m.size();
for (int i = 0; i < n; i++) {
long long m1 = m[i];
if (m0 < m1) swap(m0, m1);
long long _, im;
tie(_, im) = inv_gcd(m0, m1);
ims.push_back(im);
m0 *= m1;
}
return ims;
}
pair<long long, long long> crt(const vector<long long>& r,
const vector<long long>& m,
const vector<long long>& ims) {
assert(r.size() == m.size());
int n = int(r.size());
long long r0 = 0, m0 = 1;
for (int i = 0; i < n; i++) {
long long r1 = r[i], m1 = m[i], im = ims[i];
if (m0 < m1) swap(r0, r1), swap(m0, m1);
long long x = (r1 - r0) * im % m1;
r0 += x * m0;
m0 *= m1;
if (r0 < 0) r0 += m0;
}
return { r0, m0 };
}
using namespace std;
#define PRIME_POWER_BINOMIAL_M_MAX ((1LL << 30) - 1)
#define PRIME_POWER_BINOMIAL_N_MAX 20000000
struct prime_power_binomial {
int p, q, M;
vector<int> fac, ifac, inv;
int delta;
using i64 = long long;
using u64 = unsigned long long;
u64 iM, ip;
prime_power_binomial(int _p, int _q) : p(_p), q(_q) {
assert(p <= PRIME_POWER_BINOMIAL_M_MAX);
assert(_q > 0);
long long m = 1;
while (_q--) {
m *= p;
assert(m <= PRIME_POWER_BINOMIAL_M_MAX);
}
M = m;
iM = u64(-1) / M + 1;
ip = u64(-1) / p + 1;
enumerate();
delta = (p == 2 && q >= 3) ? 1 : M - 1;
}
inline i64 modulo_M(u64 n) {
u64 x = u64((__uint128_t(n) * iM) >> 64);
i64 r = i64(n - x * M);
if (r < 0) r += M;
return r;
}
int modpow(int a, long long e) {
int r = 1;
while (e) {
if (e & 1) r = modulo_M(1LL * r * a);
a = modulo_M(1LL * a * a);
e >>= 1;
}
return r;
}
inline i64 divide_p(u64 n) {
u64 x = u64((__uint128_t(n) * ip) >> 64);
i64 r = i64(n - x * p);
if (r < 0) x--;
return i64(x);
}
inline pair<i64, int> quorem_p(u64 n) {
u64 x = u64((__uint128_t(n) * ip) >> 64);
i64 r = i64(n - x * p);
if (r < 0) r += M, x--;
return make_pair(i64(x), r);
}
void enumerate() {
int MX = min<int>(M, PRIME_POWER_BINOMIAL_N_MAX + 10);
fac.resize(MX);
ifac.resize(MX);
inv.resize(MX);
fac[0] = ifac[0] = inv[0] = 1;
fac[1] = ifac[1] = inv[1] = 1;
for (int i = 2; i < MX; i++) {
if (i % p == 0) {
fac[i] = fac[i - 1];
fac[i + 1] = modulo_M(1LL * fac[i - 1] * (i + 1));
i++;
}
else {
fac[i] = modulo_M(1LL * fac[i - 1] * i);
}
}
ifac[MX - 1] = modpow(fac[MX - 1], M / p * (p - 1) - 1);
for (int i = MX - 2; i > 1; --i) {
if (i % p == 0) {
ifac[i] = modulo_M(1LL * ifac[i + 1] * (i + 1));
ifac[i - 1] = ifac[i];
i--;
}
else {
ifac[i] = modulo_M(1LL * ifac[i + 1] * (i + 1));
}
}
}
long long Lucas(long long n, long long m) {
int res = 1;
while (n) {
int n0, m0;
tie(n, n0) = quorem_p(n);
tie(m, m0) = quorem_p(m);
if (n0 < m0) return 0;
res = modulo_M(1LL * res * fac[n0]);
res = modulo_M(1LL * res * ifac[m0]);
res = modulo_M(1LL * res * ifac[n0 - m0]);
}
return res;
}
long long C(long long n, long long m) {
if (n < m || n < 0 || m < 0) return 0;
if (q == 1) return Lucas(n, m);
long long r = n - m;
int e0 = 0, eq = 0, i = 0;
int res = 1;
while (n) {
res = modulo_M(1LL * res * fac[modulo_M(n)]);
res = modulo_M(1LL * res * ifac[modulo_M(m)]);
res = modulo_M(1LL * res * ifac[modulo_M(r)]);
n = divide_p(n);
m = divide_p(m);
r = divide_p(r);
int eps = n - m - r;
e0 += eps;
if (e0 >= q) return 0;
if (++i >= q) eq += eps;
}
res = modulo_M(1LL * res * modpow(delta, eq));
res = modulo_M(1LL * res * modpow(p, e0));
return res;
}
};
// constraints:
// (M <= 1e7 and max(N) <= 1e18) or (M < 2^30 and max(N) <= 2e7)
struct arbitrary_mod_binomial {
int mod;
vector<int> M;
vector<long long> ims;
vector<prime_power_binomial> cs;
arbitrary_mod_binomial(long long md) : mod(md) {
assert(1 <= md);
assert(md <= PRIME_POWER_BINOMIAL_M_MAX);
for (int i = 2; i * i <= md; i++) {
if (md % i == 0) {
int j = 0, k = 1;
while (md % i == 0) md /= i, j++, k *= i;
M.push_back(k);
cs.emplace_back(i, j);
assert(M.back() == cs.back().M);
}
}
if (md != 1) {
M.push_back(md);
cs.emplace_back(md, 1);
}
assert(M.size() == cs.size());
vector<long long> ms;
for (auto& c : cs) ms.push_back(c.M);
ims = pre_im(ms);
}
long long C(long long n, long long m) {
if (mod == 1) return 0;
vector<long long> rem, d;
for (int i = 0; i < (int)cs.size(); i++) {
rem.push_back(cs[i].C(n, m));
d.push_back(M[i]);
}
return crt(rem, d, ims).first;
}
};
#undef PRIME_POWER_BINOMIAL_M_MAX
#undef PRIME_POWER_BINOMIAL_N_MAX
#include <cstring>
#include <type_traits>
using namespace std;
namespace fastio {
static constexpr int SZ = 1 << 17;
char ibuf[SZ], obuf[SZ];
int pil = 0, pir = 0, por = 0;
struct Pre {
char num[40000];
constexpr Pre() : num() {
for (int i = 0; i < 10000; i++) {
int n = i;
for (int j = 3; j >= 0; j--) {
num[i * 4 + j] = n % 10 + '0';
n /= 10;
}
}
}
} constexpr pre;
inline void load() {
memcpy(ibuf, ibuf + pil, pir - pil);
pir = pir - pil + fread(ibuf + pir - pil, 1, SZ - pir + pil, stdin);
pil = 0;
}
inline void flush() {
fwrite(obuf, 1, por, stdout);
por = 0;
}
inline void skip_space() {
if (pil + 32 > pir) load();
while (ibuf[pil] <= ' ') pil++;
}
inline void rd(char& c) {
if (pil + 32 > pir) load();
c = ibuf[pil++];
}
template <typename T>
inline void rd(T& x) {
if (pil + 32 > pir) load();
char c;
do c = ibuf[pil++];
while (c < '-');
[[maybe_unused]] bool minus = false;
if constexpr (is_signed<T>::value == true) {
if (c == '-') minus = true, c = ibuf[pil++];
}
x = 0;
while (c >= '0') {
x = x * 10 + (c & 15);
c = ibuf[pil++];
}
if constexpr (is_signed<T>::value == true) {
if (minus) x = -x;
}
}
inline void rd() {}
template <typename Head, typename... Tail>
inline void rd(Head& head, Tail&... tail) {
rd(head);
rd(tail...);
}
inline void wt(char c) {
if (por > SZ - 32) flush();
obuf[por++] = c;
}
inline void wt(bool b) {
if (por > SZ - 32) flush();
obuf[por++] = b ? '1' : '0';
}
template <typename T>
inline void wt(T x) {
if (por > SZ - 32) flush();
if (!x) {
obuf[por++] = '0';
return;
}
if constexpr (is_signed<T>::value == true) {
if (x < 0) obuf[por++] = '-', x = -x;
}
int i = 12;
char buf[16];
while (x >= 10000) {
memcpy(buf + i, pre.num + (x % 10000) * 4, 4);
x /= 10000;
i -= 4;
}
if (x < 100) {
if (x < 10) {
obuf[por] = '0' + x;
++por;
}
else {
uint32_t q = (uint32_t(x) * 205) >> 11;
uint32_t r = uint32_t(x) - q * 10;
obuf[por] = '0' + q;
obuf[por + 1] = '0' + r;
por += 2;
}
}
else {
if (x < 1000) {
memcpy(obuf + por, pre.num + (x << 2) + 1, 3);
por += 3;
}
else {
memcpy(obuf + por, pre.num + (x << 2), 4);
por += 4;
}
}
memcpy(obuf + por, buf + i + 4, 12 - i);
por += 12 - i;
}
inline void wt() {}
template <typename Head, typename... Tail>
inline void wt(Head&& head, Tail&&... tail) {
wt(head);
wt(forward<Tail>(tail)...);
}
template <typename... Args>
inline void wtn(Args&&... x) {
wt(forward<Args>(x)...);
wt('\n');
}
struct Dummy {
Dummy() { atexit(flush); }
} dummy;
} // namespace fastio
using fastio::rd;
using fastio::skip_space;
using fastio::wt;
using fastio::wtn;
#include <iostream>
int main() {
unsigned int M, N;
rd(M, N);
arbitrary_mod_binomial C(100000000);
auto val = C.C(M, N);
string str = to_string(val);
string ans;
for (int i = 0; i < 8 - str.size(); ++i)ans += '0';
ans += str;
std::cout << ans << endl;
}
kwm_t