結果
問題 | No.1510 Simple Integral |
ユーザー |
![]() |
提出日時 | 2021-05-15 01:18:44 |
言語 | C++17 (gcc 13.3.0 + boost 1.87.0) |
結果 |
AC
|
実行時間 | 36 ms / 2,000 ms |
コード長 | 9,258 bytes |
コンパイル時間 | 3,966 ms |
コンパイル使用メモリ | 239,276 KB |
最終ジャッジ日時 | 2025-01-21 12:56:20 |
ジャッジサーバーID (参考情報) |
judge1 / judge1 |
(要ログイン)
ファイルパターン | 結果 |
---|---|
sample | AC * 3 |
other | AC * 43 |
ソースコード
#include<bits/stdc++.h>using namespace std;#pragma region atcoder#include <atcoder/modint>#include <atcoder/convolution>using namespace atcoder;//using mint = modint998244353;//using mint = modint1000000007;#pragma endregion#pragma region debug for var, v, vv#define debug(var) do{std::cerr << #var << " : ";view(var);}while(0)template<typename T> void view(T e){std::cerr << e << std::endl;}template<typename T> void view(const std::vector<T>& v){for(const auto& e : v){ std::cerr << e << " "; } std::cerr << std::endl;}template<typename T> void view(const std::vector<std::vector<T> >& vv){cerr << endl;int cnt = 0;for(const auto& v : vv){cerr << cnt << "th : "; view(v); cnt++;} cerr << endl;}#pragma endregionusing ll = long long;const ll mod = 1000000007;const int inf = 1001001001;const ll INF = 1001001001001001001ll;int dx[]={1,0,-1,0};int dy[]={0,1,0,-1};template<class T, class K>bool chmax(T &a, const K b) { if (a<b) { a=b; return 1; } return 0; }template<class T, class K>bool chmin(T &a, const K b) { if (b<a) { a=b; return 1; } return 0; }ll rudiv(ll a, ll b) { return a/b+((a^b)>0&&a%b); } // 20 / 3 == 7ll rddiv(ll a, ll b) { return a/b-((a^b)<0&&a%b); } // -20 / 3 == -7ll power(ll a, ll p){ll ret = 1; while(p){if(p & 1){ret = ret * a;} a = a * a; p >>= 1;} return ret;}ll modpow(ll a, ll p, ll m){ll ret = 1; while(p){if(p & 1){ret = ret * a % m;} a = a * a % m; p >>= 1;} return ret;}/*---------------------------------------------------------------------------------------------------------------------------------*/template<class T>struct FormalPowerSeries : vector<T>{using vector<T>::vector;using vector<T>::operator=;using F = FormalPowerSeries;// operations in O(N)F operator-()const{F res(*this);for (auto &e : res) e = -e;return res;}F &operator*=(const T &g) {for (auto &e : *this) e *= g;return *this;}F &operator/=(const T &g) {assert(g != T(0));*this *= g.inv();return *this;}F &operator+=(const F &g) {int n = (*this).size(), m = g.size();for(int i = 0; i < min(n, m); i++) (*this)[i] += g[i];return *this;}F &operator-=(const F &g) {int n = (*this).size(), m = g.size();for(int i = 0; i < min(n, m); i++) (*this)[i] -= g[i];return *this;}// f -> f*z^dF &operator<<=(const int d) {int n = (*this).size();(*this).insert((*this).begin(), d, 0);(*this).resize(n);return *this;}// f -> f*z^(-d)F &operator>>=(const int d) {int n = (*this).size();(*this).erase((*this).begin(), (*this).begin() + min(n, d));(*this).resize(n);return *this;}// return 1/f, deg(1/f) = d - 1// using atcoder::convolutionF inv(int d = -1) const {int n = (*this).size();assert(n != 0 && (*this)[0] != 0);if (d == -1) d = n;assert(d > 0);F res{(*this)[0].inv()};while (res.size() < d){int m = size(res);F f(begin(*this), begin(*this) + min(n, 2*m));F r(res);f.resize(2*m), internal::butterfly(f);r.resize(2*m), internal::butterfly(r);for(int i = 0; i < 2*m; i++) f[i] *= r[i];internal::butterfly_inv(f);f.erase(f.begin(), f.begin() + m);f.resize(2*m), internal::butterfly(f);for(int i = 0; i < 2*m; i++) f[i] *= r[i];internal::butterfly_inv(f);T iz = T(2*m).inv(); iz *= -iz;for(int i = 0; i < m; i++) f[i] *= iz;res.insert(res.end(), f.begin(), f.begin() + m);}return {res.begin(), res.begin() + d};}// computation for f*g & f/g, use either one of them//// fast: FMT-friendly modulus only -> O((N + M)log(N + M))/* F &operator*=(const F &g) {int n = (*this).size();*this = convolution(*this, g);(*this).resize(n);return *this;}F &operator/=(const F &g) {int n = (*this).size();*this = convolution(*this, g.inv(n));(*this).resize(n);return *this;}*///// naive -> O(NM)/* F &operator*=(const F &g) {int n = (*this).size(), m = g.size();for(int i = n - 1; i >= 0; i--) {(*this)[i] *= g[0];for(int j = 1; j < min(i+1, m); j++) (*this)[i] += (*this)[i-j] * g[j];}return *this;}F &operator/=(const F &g) {assert(g[0] != T(0));T ig0 = g[0].inv();int n = (*this).size(), m = g.size();for(int i = 0; i < n; i++) {for(int j = 1; j < min(i+1, m); j++) (*this)[i] -= (*this)[i-j] * g[j];(*this)[i] *= ig0;}return *this;}*/// sparse(the coefficiets in g have K non-0s) O(NK)// give the pair of (degree, coefficient) in this function// f -> f*gF &operator*=(vector<pair<int, T>> g) {int n = (*this).size();auto [d, c] = g.front();if (d == 0) g.erase(g.begin());else c = 0;for(int i = n - 1; i >= 0; i--){(*this)[i] *= c;for (auto &[j, b] : g) {if (j > i) break;(*this)[i] += (*this)[i-j] * b;}}return *this;}// f -> f/gF &operator/=(vector<pair<int, T>> g) {int n = (*this).size();auto [d, c] = g.front();assert(d == 0 && c != T(0));T ic = c.inv();g.erase(g.begin());for(int i = 0; i < n; i++) {for (auto &[j, b] : g) {if (j > i) break;(*this)[i] -= (*this)[i-j] * b;}(*this)[i] *= ic;}return *this;}// f -> f*(1 + c*z^d) in O(N)void multiply(const int d, const T c) {int n = (*this).size();if (c == T(1)) for(int i = n-d-1; i >= 0; i--) (*this)[i+d] += (*this)[i];else if (c == T(-1)) for(int i = n-d-1; i >= 0; i--) (*this)[i+d] -= (*this)[i];else for(int i = n-d-1; i >= 0; i--) (*this)[i+d] += (*this)[i] * c;}// f -> f/(1 + c*z^d) in O(N)void divide(const int d, const T c) {int n = (*this).size();if (c == T(1)) for(int i = 0; i < n-d; i++) (*this)[i+d] -= (*this)[i];else if (c == T(-1)) for(int i = 0; i < n-d; i++) (*this)[i+d] += (*this)[i];else for(int i = 0; i < n-d; i++) (*this)[i+d] -= (*this)[i] * c;}// return f(a)T eval(const T &a) const {T x(1), res(0);for (auto e : *this) res += e * x, x *= a;return res;}F operator*(const T &g)const{return F(*this) *= g;}F operator/(const T &g)const{return F(*this) /= g;}F operator+(const F &g)const{return F(*this) += g;}F operator-(const F &g)const{return F(*this) -= g;}F operator<<(const int d)const{return F(*this) <<= d;}F operator>>(const int d)const{return F(*this) >>= d;}F operator*(const F &g)const{return F(*this) *= g;}F operator/(const F &g)const{return F(*this) /= g;}F operator*(vector<pair<int, T>> g)const{ return F(*this) *= g;}F operator/(vector<pair<int, T>> g)const{ return F(*this) /= g;}};using mint = modint998244353;using fps = FormalPowerSeries<mint>;using sfps = vector<pair<int, mint>>;mint a[110];typedef long long ll;ll MOD = 998244353;mint pw(mint a,ll x){mint ret = 1;while(x){if(x&1) (ret *= a);(a *= a); x /= 2;}return ret;}vector<ll> d;mint A[110];ll cnt[1000010] = {};int main(){cin.tie(nullptr);ios::sync_with_stdio(false);//cout << fixed << setprecision(20);ll i,j,n;mint pr = -1,u = -1; cin >> n;MOD--;for(i=1;i*i<MOD;i++){if(MOD%i==0){if(i!=1 && i!=MOD) d.push_back(i);if(i!=1 && i!=MOD) d.push_back(MOD/i);}}MOD++;for(i=2;i<MOD;i++){bool f = true;for(int x:d){if(pw(i,x)==1) f = false;}if(f){pr = i;break;}}u = pw(pr,(MOD - 1)/4);if(pr==-1 || u==-1) exit(1);if(pw(u,2)!=MOD - 1) exit(1);for(i=0;i<n;i++){int x; cin >> x; A[i] = x; cnt[x]++;}mint ans = 0;for(i=1;i<=1000000;i++){if(!cnt[i]) continue;fps f(201);f[0] = 1;ll num = 0;for(j=0;j<n;j++){if(i==A[j]){num++;fps g(2);g[1] = (mint)1; g[0] = 2*(mint)u*i;f = convolution(f,g);continue;}fps g(3);g[2] = 1; g[1] = (mint)2*u*i; g[0] = A[j]*A[j] - (mint)i*i;f = convolution(f,g);}f = f.inv();ans += f[num - 1];}ans *= (mint)2*u;cout << ans.val() << endl;}/** review your code when you get WA (typo? index?)* int overflow, array bounds·* special cases (n=1?)*/