結果
問題 | No.2459 Stampaholic (Hard) |
ユーザー |
👑 ![]() |
提出日時 | 2023-08-25 04:32:11 |
言語 | C++17 (gcc 13.3.0 + boost 1.87.0) |
結果 |
AC
|
実行時間 | 637 ms / 4,000 ms |
コード長 | 2,860 bytes |
コンパイル時間 | 1,989 ms |
コンパイル使用メモリ | 117,348 KB |
最終ジャッジ日時 | 2025-02-16 13:25:53 |
ジャッジサーバーID (参考情報) |
judge5 / judge1 |
(要ログイン)
ファイルパターン | 結果 |
---|---|
sample | AC * 3 |
other | AC * 19 |
ソースコード
#include <iostream>#include <string>#include <vector>#include <algorithm>#include <atcoder/modint>#include <atcoder/convolution>using namespace std;using i32 = int;using u32 = unsigned int;using i64 = long long;using u64 = unsigned long long;#define rep(i,n) for(int i=0; i<(int)(n); i++)const i64 INF = 1001001001001001001;using Modint = atcoder::static_modint<998244353>;vector<Modint> invfps(vector<Modint> x, int n){Modint ix0 = x[0].inv();for(auto& a : x) a *= ix0;vector<Modint> res;res.push_back(1);int z = 1;while(z < n){if((int)x.size() < z*2) x.resize(z*2);auto xx = vector(x.begin(), x.begin() + z*2);auto tmp = atcoder::convolution(res, xx);tmp.resize(z*2);tmp = atcoder::convolution(tmp, res);rep(i,z) res.push_back(-tmp[z+i]);z *= 2;}for(auto& a : res) a *= ix0;return vector(res.begin(), res.begin() + z);}vector<Modint> factorial;vector<Modint> invFactorial;void extend(int n){factorial.resize(1+n);factorial[0] = 1;for(int i=1; i<=n; i++) factorial[i] = factorial[i-1] * i;invFactorial.resize(1+n);invFactorial[n] = factorial[n].inv();for(int i=n; i>=1; i--) invFactorial[i-1] = invFactorial[i] * i;}Modint comb(int n, int k){ return factorial[n] * invFactorial[k] * invFactorial[n-k]; }vector<Modint> powsum(Modint q, int n){vector<Modint> a(n+1), b(n+1);a[0] = q; rep(i,n) a[i+1] = a[i] * q;rep(i,n+1) a[i] *= invFactorial[i+1];rep(i,n+1) b[i] = invFactorial[i+1];auto invb = invfps(move(b), n+1);a = atcoder::convolution(a, invb);rep(i,n+1) a[i] *= factorial[i];return vector(a.begin(), a.begin() + (n+1));}int main(){i64 H, W, N, K; cin >> H >> W >> N >> K;extend(N+2);i64 f = (H-K+1) * (W-K+1);Modint p = (Modint(f)).inv().pow(N);Modint ans = 0;i64 xcent = min(K, W-K+1);i64 xcent_width = W - xcent * 2 + 2;i64 ycent = min(K, H-K+1);i64 ycent_width = H - ycent * 2 + 2;auto xpowsum = powsum(xcent, N); xpowsum[0] -= 1;auto ypowsum = powsum(ycent, N); ypowsum[0] -= 1;vector<Modint> powf(N+1);powf[0] = Modint((N%2==0) ? 1 : -1);rep(i,N) powf[i+1] = powf[i] * (-f);vector<Modint> pow_xcent(N+1);pow_xcent[0] = 1;rep(i,N) pow_xcent[i+1] = pow_xcent[i] * xcent;vector<Modint> pow_ycent(N+1);pow_ycent[0] = 1;rep(i,N) pow_ycent[i+1] = pow_ycent[i] * ycent;rep(k,N+1) ans += comb(N,k) * powf[N-k] * (xpowsum[k] * 2 + pow_xcent[k] * xcent_width) * (ypowsum[k] * 2 + pow_ycent[k] * ycent_width);ans *= p;ans = Modint(H*W) - ans;cout << ans.val() << endl;return 0;}struct ios_do_not_sync{ios_do_not_sync(){ios::sync_with_stdio(false);cin.tie(nullptr);}} ios_do_not_sync_instance;