#ifndef LOCAL #pragma GCC optimize ("O3") #pragma GCC optimize ("unroll-loops") #endif #include using namespace std; using ll=long long; #define int ll bool dbg=false; #define rng(i,a,b) for(int i=int(a);i=int(a);i--) #define per(i,b) gnr(i,0,b) #define pb push_back #define eb emplace_back #define a first #define b second #define bg begin() #define ed end() #define all(x) x.bg,x.ed #define si(x) int(x.size()) #ifdef LOCAL #define dmp(x) cerr<<__LINE__<<" "<<#x<<" "< bool chmax(t&a,u b){if(a bool chmin(t&a,u b){if(b using vc=vector; template using vvc=vc>; using pi=pair; using vi=vc; using vvi=vc>; template ostream& operator<<(ostream& os,const pair& p){ return os<<"{"< ostream& operator<<(ostream& os,const vc& v){ os<<"{"; for(auto e:v)os< void dmpr(ostream&os,const T&t,const Args&... args){ os< ostream& operator<<(ostream&os,const array&a){ return os<(all(a)); } ll rand_int(ll l, ll r) { //[l, r] //#ifdef LOCAL static mt19937_64 gen; /*#else static mt19937_64 gen(chrono::steady_clock::now().time_since_epoch().count()); #endif*/ return uniform_int_distribution(l, r)(gen); } ll rand_int(ll k){ //[0,k) return rand_int(0,k-1); } string rand_string(int n,char lw,char up){ string s(n,'?'); rep(i,n)s[i]=rand_int(lw,up); return s; } int current_run_id,run_batch_size=1000; int calc_random_limit(){ return current_run_id/run_batch_size+1; } template void generate_single(t&a){ a=rand_int(1,calc_random_limit()); } void generate_single(string&a){ int n;generate_single(n); a=rand_string(n,'a','b'); } template void generate_single(pair&a){ generate_single(a.a); generate_single(a.b); } //https://trap.jp/post/1224/ template void input(Args&... a){ if(dbg){ (generate_single(a),...); }else{ #ifdef USE_FAST_IO sc.read(a...); #else (cin >> ... >> a); #endif } } #define INT(...) int __VA_ARGS__;input(__VA_ARGS__) #define LL(...) ll __VA_ARGS__;input(__VA_ARGS__) #define ULL(...) ull __VA_ARGS__;input(__VA_ARGS__) #define STR(...) string __VA_ARGS__;input(__VA_ARGS__) #define CHR(...) char __VA_ARGS__;input(__VA_ARGS__) #define DBL(...) double __VA_ARGS__;input(__VA_ARGS__) #define LD(...) ld __VA_ARGS__;input(__VA_ARGS__) #define overload3(a,b,c,d,...) d #define VI2(name,size) vi name(size);rep(i_##name,size)input(name[i_##name]); #define VI3(name,size,offset) vi name(size);rep(i_##name,size)input(name[i_##name]),name[i_##name]+=offset; #define VI(...) overload3(__VA_ARGS__,VI3,VI2)(__VA_ARGS__) #define VPI(name,size) vc name(size);rep(i_##name,size)input(name[i_##name].a,name[i_##name].b); #define VVI(name,sizeN,sizeM) vvi name(sizeN,vi(sizeM));\ rep(i_##name,sizeN)rep(j_##name,sizeM)input(name[i_##name][j_##name]); #define VS(name,size) vc name(size);rep(i_##name,size)input(name[i_##name]); #define VMI(name,size) vc name(size);rep(i_##name,size){INT(tmp_##name);name[i_##name]=tmp_##name;}; #define overload5(a,b,c,d,e,f,...) f #define VVC4(type,name,sizeN,sizeM) vvc name(sizeN,vc(sizeM)); #define VVC5(type,name,sizeN,sizeM,ini) vvc name(sizeN,vc(sizeM,ini)); #define VVC(...) overload5(__VA_ARGS__,VVC5,VVC4)(__VA_ARGS__) template T vvvc(T v){ return v; } template auto vvvc(int n,T v,Args...args){ return vector(n,vvvc(v,args...)); } template void print_tuple(ostream&,const T&){ } template void print_tuple(ostream&os,const T&t){ if(i)os<<","; os<(t); print_tuple(os,t); } template ostream& operator<<(ostream&os,const tuple&t){ os<<"{"; print_tuple<0,tuple,Args...>(os,t); return os<<"}"; } void printsuc(int suc){ #ifdef USE_FAST_IO if(suc==1)pr.write('\n'); if(suc==2)pr.write(' '); #else if(suc==1){ if(dbg)cout< void print_single(t x,int suc=1){ #ifdef USE_FAST_IO pr.write(x); #else cout< void print_single(const pair&p,int suc=1){ print_single(p.a,2); print_single(p.b,suc); } template void print_single(const vector&v,int suc=1){ rep(i,v.size()) print_single(v[i],i==int(v.size())-1?3:2); printsuc(suc); } template void print_single(const array&v,int suc=1){ rep(i,N) print_single(v[i],i==int(N)-1?3:2); printsuc(suc); } template void print(const T&t){ print_single(t); } template void print(const T&t,const Args&...args){ print_single(t,2); print(args...); } template void printvv(const vvc&vs){ for(const auto&row:vs)print(row); } string readString(){ string s; cin>>s; return s; } template T sq(const T& t){ return t*t; } void YES(bool ex=true){ cout<<"YES\n"; if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } void NO(bool ex=true){ cout<<"NO\n"; if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } void Yes(bool ex=true){ cout<<"Yes\n"; if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } void No(bool ex=true){ cout<<"No\n"; if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } //#define CAPITAL /* void yes(bool ex=true){ #ifdef CAPITAL cout<<"YES"<<"\n"; #else cout<<"Yes"<<"\n"; #endif if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } void no(bool ex=true){ #ifdef CAPITAL cout<<"NO"<<"\n"; #else cout<<"No"<<"\n"; #endif if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif }*/ void possible(bool ex=true){ #ifdef CAPITAL cout<<"POSSIBLE"<<"\n"; #else cout<<"Possible"<<"\n"; #endif if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } void impossible(bool ex=true){ #ifdef CAPITAL cout<<"IMPOSSIBLE"<<"\n"; #else cout<<"Impossible"<<"\n"; #endif if(ex)exit(0); #ifdef LOCAL cout.flush(); #endif } constexpr ll ten(int n){ return n==0?1:ten(n-1)*10; } const ll infLL=LLONG_MAX/3; #ifdef int const int inf=infLL; #else const int inf=INT_MAX/2-100; #endif int topbit(signed t){ return t==0?-1:31-__builtin_clz(t); } int topbit(ll t){ return t==0?-1:63-__builtin_clzll(t); } int topbit(ull t){ return t==0?-1:63-__builtin_clzll(t); } int botbit(signed a){ return a==0?32:__builtin_ctz(a); } int botbit(ll a){ return a==0?64:__builtin_ctzll(a); } int botbit(ull a){ return a==0?64:__builtin_ctzll(a); } int popcount(signed t){ return __builtin_popcount(t); } int popcount(ll t){ return __builtin_popcountll(t); } int popcount(ull t){ return __builtin_popcountll(t); } int bitparity(ll t){ return __builtin_parityll(t); } bool ispow2(int i){ return i&&(i&-i)==i; } ll mask(int i){ return (ll(1)< void mkuni(S&v){ sort(all(v)); v.erase(unique(all(v)),v.ed); } template bool isuni(vc v){ int s=si(v); mkuni(v); return si(v)==s; } bool isperm(const vi&p){ int n=si(p); vc used(n); for(auto v:p){ if(!inc(0,v,n-1))return false; if(used[v])return false; used[v]=true; } return true; } template void myshuffle(vc&a){ rep(i,si(a))swap(a[i],a[rand_int(0,i)]); } template int lwb(const S&v,const u&a){ return lower_bound(all(v),a)-v.bg; } template bool bis(const vc&v,const u&a){ return binary_search(all(v),a); } //VERIFY: yosupo //KUPC2017J //AOJDSL1A //without rank struct unionfind{ vi p,s; int c; unionfind(int n):p(n,-1),s(n,1),c(n){} void clear(){ fill(all(p),-1); fill(all(s),1); c=si(p); } int find(int a){ return p[a]==-1?a:(p[a]=find(p[a])); } //set b to a child of a bool unite(int a,int b){ a=find(a); b=find(b); if(a==b)return false; p[b]=a; s[a]+=s[b]; c--; return true; } bool same(int a,int b){ return find(a)==find(b); } int sz(int a){ return s[find(a)]; } }; vvc readGraph(int n,int m){ vvc g(n); rep(i,m){ int a,b; cin>>a>>b; //sc.read(a,b); a--;b--; g[a].pb(b); g[b].pb(a); } return g; } vvc initUG(int n,const vc&es){ vvc g(n); for(auto [a,b]:es){ g[a].pb(b); g[b].pb(a); } return g; } vvc rand_tree(int n){ vvc t(n); unionfind uf(n); while(uf.c>1){ int a=rand_int(n); int b=rand_int(n); if(uf.unite(a,b)){ t[a].pb(b); t[b].pb(a); } } return t; } vvc readTree(int n){ if(dbg){ return rand_tree(n); }else{ return readGraph(n,n-1); } } vi readRooted(int n){ assert(!dbg); vi par(n,-1); rng(i,1,n){ input(par[i]); par[i]--; assert(inc(0,par[i],i-1)); } return par; } void printTree(const vvc t){ int n=si(t); int degsum=0; rep(i,n)degsum+=si(t[i]); if(degsum==n-1){ //directed rep(i,si(t))for(auto j:t[i]){ print(i+1,j+1); } }else if(degsum==2*(n-1)){ //undirected rep(i,si(t))for(auto j:t[i])if(i vc presum(const vc&a){ vc s(si(a)+1); rep(i,si(a))s[i+1]=s[i]+a[i]; return s; } vc presum(const vi&a){ vc s(si(a)+1); rep(i,si(a))s[i+1]=s[i]+a[i]; return s; } //BIT で数列を管理するときに使う (CF850C) template vc predif(vc a){ gnr(i,1,si(a))a[i]-=a[i-1]; return a; } template vvc imos(const vvc&a){ int n=si(a),m=si(a[0]); vvc b(n+1,vc(m+1)); rep(i,n)rep(j,m) b[i+1][j+1]=b[i+1][j]+b[i][j+1]-b[i][j]+a[i][j]; return b; } //verify してないや void transvvc(int&n,int&m){ swap(n,m); } template void transvvc(int&n,int&m,vvc&a,Args&...args){ assert(si(a)==n); vvc b(m,vi(n)); rep(i,n){ assert(si(a[i])==m); rep(j,m)b[j][i]=a[i][j]; } a.swap(b); transvvc(n,m,args...); } //CF854E void rotvvc(int&n,int&m){ swap(n,m); } template void rotvvc(int&n,int&m,vvc&a,Args&...args){ assert(si(a)==n); vvc b(m,vi(n)); rep(i,n){ assert(si(a[i])==m); rep(j,m)b[m-1-j][i]=a[i][j]; } a.swap(b); rotvvc(n,m,args...); } //ソートして i 番目が idx[i] //CF850C template vi sortidx(const vc&a){ int n=si(a); vi idx(n);iota(all(idx),0); sort(all(idx),[&](int i,int j){return a[i] vc a_idx(const vc&a,const vi&idx){ int n=si(a); assert(si(idx)==n); vc vs(n); rep(i,n)vs[i]=a[idx[i]]; return vs; } //CF850C vi invperm(const vi&p){ int n=si(p); vi q(n); rep(i,n)q[p[i]]=i; return q; } template s SUM(const vc&a){ return accumulate(all(a),s(0)); } ll SUM(const vi&a){ return accumulate(all(a),ll(0)); } template s SUM(const array&a){ return accumulate(all(a),s(0)); } template t MAX(const vc&a){ return *max_element(all(a)); } template pair MAXi(const vc&a){ auto itr=max_element(all(a)); return mp(*itr,itr-a.bg); } template auto MIN(const A&a){ return *min_element(all(a)); } template pair MINi(const vc&a){ auto itr=min_element(all(a)); return mp(*itr,itr-a.bg); } vi vid(int n){ vi res(n);iota(all(res),0); return res; } template void soin(S&s){ sort(all(s)); } template void soin(S&s,F&&f){ sort(all(s),forward(f)); } template S soout(S s){ soin(s); return s; } template void rein(S&s){ reverse(all(s)); } template S reout(S s){ rein(s); return s; } template pair&operator+=(pair&a,pair b){ a.a+=b.a;a.b+=b.b;return a;} template pair&operator-=(pair&a,pair b){ a.a-=b.a;a.b-=b.b;return a;} template pair operator+(pair a,pair b){return mp(a.a+b.a,a.b+b.b);} template pair operator-(pair a,pair b){return mp(a.a-b.a,a.b-b.b);} template pair&operator*=(pair&a,v b){ a.a*=b;a.b*=b;return a;} template pair operator*(pair a,v b){return a*=b;} template pair operator-(pair a){return mp(-a.a,-a.b);} namespace std{ template istream&operator>>(istream&is,pair&a){ return is>>a.a>>a.b; } } template array&operator+=(array&a,const array&b){ rep(i,n)a[i]+=b[i]; return a; } template array&operator-=(array&a,const array&b){ rep(i,n)a[i]-=b[i]; return a; } template array&operator*=(array&a,v b){ rep(i,n)a[i]*=b; return a; } template array operator+(array a,const array&b){return a+=b;} template array operator-(array a,const array&b){return a-=b;} template array operator*(array a,v b){return a*=b;} template t gpp(vc&vs){ assert(si(vs)); t res=move(vs.back()); vs.pop_back(); return res; } template void pb(vc&a,const vc&b){ a.insert(a.ed,all(b)); } template vc cat(vc a,Args&&...b){ (pb(a,forward(b)),...); return a; } template vc& operator+=(vc&a,u x){ for(auto&v:a)v+=x; return a; } template vc operator+(vc a,u x){ return a+=x; } template vc& operator+=(vc&a,const vc&b){ a.resize(max(si(a),si(b))); rep(i,si(b))a[i]+=b[i]; return a; } template vc operator+(const vc&a,const vc&b){ vc c(max(si(a),si(b))); rep(i,si(a))c[i]+=a[i]; rep(i,si(b))c[i]+=b[i]; return c; } template vc& operator-=(vc&a,u x){ for(auto&v:a)v-=x; return a; } template vc operator-(vc a,u x){ return a-=x; } template vc& operator-=(vc&a,const vc&b){ a.resize(max(si(a),si(b))); rep(i,si(b))a[i]-=b[i]; return a; } template vc operator-(const vc&a,const vc&b){ vc c(max(si(a),si(b))); rep(i,si(a))c[i]+=a[i]; rep(i,si(b))c[i]-=b[i]; return c; } template vc& operator*=(vc&a,u x){ for(auto&v:a)v*=x; return a; } template vc operator*(vc a,u x){ return a*=x; } template vc& operator/=(vc&a,u x){ for(auto&v:a)v/=x; return a; } template vc operator/(vc a,u x){ return a/=x; } template vc& operator<<=(vc&a,int k){ assert(k>=0); a.insert(a.bg,k,t(0)); return a; } template vc operator<<(vc a,int k){ return a<<=k; } template vc& operator>>=(vc&a,int k){ if(si(a)<=k)a.clear(); else a.erase(a.bg,a.bg+k); return a; } template vc operator>>(vc a,int k){ return a>>=k; } //消した要素の個数を返してくれる //not verified template int remval(vc&a,const u&v){ auto itr=remove(all(a),v); int res=a.ed-itr; a.erase(itr,a.ed); return res; } //消した要素の個数を返してくれる //UCUP 2-8-F template int remif(vc&a,F f){ auto itr=remove_if(all(a),f); int res=a.ed-itr; a.erase(itr,a.ed); return res; } template void rempos(vc&a,int i){ assert(inc(0,i,si(a)-1)); a.erase(a.bg+i); } template void fila(VS&vs,const u&a){ fill(all(vs),a); } template int findid(const vc&vs,const u&a){ auto itr=find(all(vs),a); if(itr==vs.ed)return -1; else return itr-vs.bg; } template void rtt(vc&vs,int i){ rotate(vs.bg,vs.bg+i,vs.ed); } //Multiuni2023-8 C //f(lw)=false,...,f(n-1)=false,f(n)=true,...,f(up)=true, //のときに n を返す template int find_min_true(int lw,int up,F f){ while(up-lw>1){ const int mid=(lw+up)/2; if(f(mid))up=mid; else lw=mid; } return up; } //f(lw)=true,f(up)=false template int find_max_true(int lw,int up,F f){ while(up-lw>1){ const int mid=(lw+up)/2; if(f(mid))lw=mid; else up=mid; } return lw; } template using pqmin=priority_queue,greater>; template using pqmax=priority_queue; using T=tuple; //mint107 は verify してねえ //#define DYNAMIC_MOD struct modinfo{uint mod,root; #ifdef DYNAMIC_MOD constexpr modinfo(uint m,uint r):mod(m),root(r),im(0){set_mod(m);} ull im; constexpr void set_mod(uint m){ mod=m; im=ull(-1)/m+1; } uint product(uint a,uint b)const{ ull z=ull(a)*b; uint x=((unsigned __int128)z*im)>>64; uint v=uint(z)-x*mod; return v struct modular{ static constexpr uint const &mod=ref.mod; static modular root(){return modular(ref.root);} uint v; //modular(initializer_listls):v(*ls.bg){} modular(ll vv=0){s(vv%mod+mod);} modular& s(uint vv){ v=vv>=1; } return res; } modular inv()const{return pow(mod-2);} /*modular inv()const{ int x,y; int g=extgcd(v,mod,x,y); assert(g==1); if(x<0)x+=mod; return modular(x); }*/ friend modular operator+(ll x,const modular&y){ return modular(x)+y; } friend modular operator-(ll x,const modular&y){ return modular(x)-y; } friend modular operator*(ll x,const modular&y){ return modular(x)*y; } friend modular operator/(ll x,const modular&y){ return modular(x)/y; } friend ostream& operator<<(ostream&os,const modular&m){ return os<>(istream&is,modular&m){ ll x;is>>x; m=modular(x); return is; } bool operator<(const modular&r)const{return v void inplace_fmt(const int n,mint*const f,bool inv){ static constexpr uint mod=mint::mod; static constexpr uint mod2=mod*2; static constexpr int L=30; static mint g[L],ig[L],p2[L]; if(g[0].v==0){ rep(i,L){ mint w=-mint::root().pow(((mod-1)>>(i+2))*3); g[i]=w; ig[i]=w.inv(); p2[i]=mint(1<>=1){//input:[0,mod) rep(i,b){ uint x=f[i+b].v; f[i+b].v=f[i].v+mod-x; f[i].v+=x; } } if(b>>=1){//input:[0,mod*2) mint p=1; for(int i=0,k=0;i>=1){//input:[0,mod*3) mint p=1; for(int i=0,k=0;i>=1){//input:[0,mod*4) mint p=1; for(int i=0,k=0;i void inplace_fmt(vector&f,bool inv){ inplace_fmt(si(f),f.data(),inv); } //size of input must be a power of 2 //output elements are in the range [0,mod*4) template void half_fmt(const int n,mint*const f){ static constexpr uint mod=mint::mod; static constexpr uint mod2=mod*2; static const int L=30; static mint g[L],h[L]; if(g[0].v==0){ rep(i,L){ g[i]=-mint::root().pow(((mod-1)>>(i+2))*3); h[i]=mint::root().pow((mod-1)>>(i+2)); } } int b=n; int lv=0; if(b>>=1){//input:[0,mod) mint p=h[lv++]; for(int i=0,k=0;i>=1){//input:[0,mod*2) mint p=h[lv++]; for(int i=0,k=0;i>=1){//input:[0,mod*3) mint p=h[lv++]; for(int i=0,k=0;i>=1){//input:[0,mod*4) mint p=h[lv++]; for(int i=0,k=0;i void half_fmt(vector&f){ half_fmt(si(f),f.data()); } #ifdef USE_GOOD_MOD template vc multiply(vc x,const vc&y,bool same=false){ if(x.empty()||y.empty())return {}; int n=si(x)+si(y)-1; int s=1; while(s z; z.clear();z.resize(s); rep(i,si(y))z[i]=y[i]; inplace_fmt(z,false); rep(i,s)x[i]*=z[i]; }else{ rep(i,s)x[i]*=x[i]; } inplace_fmt(x,true);x.resize(n); return x; } template vc multiply_givenlength(vc x,const vc&y,bool same=false){ int s=si(x); assert(ispow2(s)); assert(si(y)); x.resize(s);inplace_fmt(x,false); if(!same){ static vc z; z.clear();z.resize(s); rep(i,si(y))z[i]=y[i]; inplace_fmt(z,false); rep(i,s)x[i]*=z[i]; }else{ rep(i,s)x[i]*=x[i]; } inplace_fmt(x,true); return x; } #else #ifndef USE_MOD64 //59501818244292734739283969-1=5.95*10^25 までの値を正しく計算 //最終的な列の大きさが 2^24 までなら動く //最終的な列の大きさが 2^20 以下のときは,下の 3 つの素数を使ったほうが速い(は?) //VERIFY: yosupo //Yukicoder No980 (same=true) namespace arbitrary_convolution{ constexpr modinfo base0{167772161,3};//2^25 * 5 + 1 constexpr modinfo base1{469762049,3};//2^26 * 7 + 1 constexpr modinfo base2{754974721,11};//2^24 * 45 + 1 //extern constexpr modinfo base0{1045430273,3};//2^20 * 997 + 1 //extern constexpr modinfo base1{1051721729,6};//2^20 * 1003 + 1 //extern constexpr modinfo base2{1053818881,7};//2^20 * 1005 + 1 using mint0=modular; using mint1=modular; using mint2=modular; template vc sub(const vc&x,const vc&y,bool same=false){ int n=si(x)+si(y)-1; int s=1; while(s z(s);rep(i,si(x))z[i]=x[i].v; inplace_fmt(z,false); if(!same){ vc w(s);rep(i,si(y))w[i]=y[i].v; inplace_fmt(w,false); rep(i,s)z[i]*=w[i]; }else{ rep(i,s)z[i]*=z[i]; } inplace_fmt(z,true);z.resize(n); return z; } template vc multiply(const vc&x,const vc&y,bool same=false){ auto d0=sub(x,y,same); auto d1=sub(x,y,same); auto d2=sub(x,y,same); int n=si(d0); vc res(n); static const mint1 r01=mint1(mint0::mod).inv(); static const mint2 r02=mint2(mint0::mod).inv(); static const mint2 r12=mint2(mint1::mod).inv(); static const mint2 r02r12=r02*r12; static const mint w1=mint(mint0::mod); static const mint w2=w1*mint(mint1::mod); rep(i,n){ ull a=d0[i].v; ull b=(d1[i].v+mint1::mod-a)*r01.v%mint1::mod; ull c=((d2[i].v+mint2::mod-a)*r02r12.v+(mint2::mod-b)*r12.v)%mint2::mod; res[i].v=(a+b*w1.v+c*w2.v)%mint::mod; } return res; } template vc&sub_givenlength(const vc&x,const vc&y,bool same=false){ int s=si(x); assert(ispow2(s)); assert(si(y)==s); static vc z; z.clear();z.resize(s); rep(i,si(x))z[i]=x[i].v; inplace_fmt(z,false); if(!same){ static vc w; w.clear();w.resize(s); rep(i,si(y))w[i]=y[i].v; inplace_fmt(w,false); rep(i,s)z[i]*=w[i]; }else{ rep(i,s)z[i]*=z[i]; } inplace_fmt(z,true); return z; } template vc multiply_givenlength(vc x,const vc&y,bool same=false){ auto&d0=sub_givenlength(x,y,same); auto&d1=sub_givenlength(x,y,same); auto&d2=sub_givenlength(x,y,same); int n=si(d0); x.resize(n); static const mint1 r01=mint1(mint0::mod).inv(); static const mint2 r02=mint2(mint0::mod).inv(); static const mint2 r12=mint2(mint1::mod).inv(); static const mint2 r02r12=r02*r12; static const mint w1=mint(mint0::mod); static const mint w2=w1*mint(mint1::mod); rep(i,n){ ull a=d0[i].v; ull b=(d1[i].v+mint1::mod-a)*r01.v%mint1::mod; ull c=((d2[i].v+mint2::mod-a)*r02r12.v+(mint2::mod-b)*r12.v)%mint2::mod; x[i].v=(a+b*w1.v+c*w2.v)%mint::mod; } return x; } } using arbitrary_convolution::multiply; using arbitrary_convolution::multiply_givenlength; #else //multiuni-2025-1-D //mod 10^10 くらいで使用 //列の長さは 10^5 くらいかな namespace arbitrary_convolution_64{ constexpr modinfo base0{1045430273,3};//2^20 * 997 + 1 constexpr modinfo base1{1051721729,6};//2^20 * 1003 + 1 constexpr modinfo base2{1053818881,7};//2^20 * 1005 + 1 using mint0=modular; using mint1=modular; using mint2=modular; template vc sub(const vc&x,const vc&y,bool same=false){ int n=si(x)+si(y)-1; int s=1; while(s z(s);rep(i,si(x))z[i]=x[i].val(); inplace_fmt(z,false); if(!same){ vc w(s);rep(i,si(y))w[i]=y[i].val(); inplace_fmt(w,false); rep(i,s)z[i]*=w[i]; }else{ rep(i,s)z[i]*=z[i]; } inplace_fmt(z,true);z.resize(n); return z; } template vc multiply(const vc&x,const vc&y,bool same=false){ auto d0=sub(x,y,same); auto d1=sub(x,y,same); auto d2=sub(x,y,same); int n=si(d0); vc res(n); static const mint1 r01=mint1(mint0::mod).inv(); static const mint2 r02=mint2(mint0::mod).inv(); static const mint2 r12=mint2(mint1::mod).inv(); static const mint2 r02r12=r02*r12; using D=__int128; static const D w1=mint0::mod; static const D w2=(D)mint0::mod*mint1::mod%mint::mod; rep(i,n){ ull a=d0[i].v; ull b=(d1[i].v+mint1::mod-a)*r01.v%mint1::mod; ull c=((d2[i].v+mint2::mod-a)*r02r12.v+(mint2::mod-b)*r12.v)%mint2::mod; res[i]=(a+b*w1+c*w2)%mint::mod; } return res; } } using arbitrary_convolution_64::multiply; #endif #endif //UTPC2021 C namespace integer_convolution{ extern constexpr modinfo base0{1045430273,3};//2^20 * 997 + 1 extern constexpr modinfo base1{1051721729,6};//2^20 * 1003 + 1 //extern constexpr modinfo base0{469762049,3};//2^26 * 7 + 1 //extern constexpr modinfo base1{754974721,11};//2^24 * 45 + 1 using mint0=modular; using mint1=modular; template vc sub(const vi&x,const vi&y,bool same=false){ int n=si(x)+si(y)-1; int s=1; while(s z(s);rep(i,si(x))z[i]=x[i]; inplace_fmt(z,false); if(!same){ vc w(s);rep(i,si(y))w[i]=y[i]; inplace_fmt(w,false); rep(i,s)z[i]*=w[i]; }else{ rep(i,s)z[i]*=z[i]; } inplace_fmt(z,true);z.resize(n); return z; } vi multiply(const vi&x,const vi&y,bool same=false){ auto d0=sub(x,y,same); auto d1=sub(x,y,same); const mint1 r=mint1(mint0::mod).inv(); const ll v=ll(mint0::mod)*mint1::mod; int n=si(d0); vi res(n); rep(i,n){ res[i]=d0[i].v+(r*(d1[i]-d0[i].v)).v*(ull)mint0::mod; if(res[i]>v/2)res[i]-=v; } return res; } } //最大で 1< vc large_convolution(const vc&a,const vc&b,int mx){ int n=si(a),m=si(b); vc c(n+m-1); int len=1<<(mx-1); for(int i=0;i(a.bg+i,a.bg+i+x),vc(b.bg+j,b.bg+j+y)); rep(k,si(d)) c[i+j+k]+=d[k]; } } return c; } //input A: N 次,B ?,M //output D: M 次多項式 //C を M 次多項式として //[x^N] A*B*C = [x^M] D*C //となるような D を返す //CF796F template vc transpose_advance(const vc&a,const vc&b,int m){ int n=si(a)-1; auto d=multiply(a,b); vc res(m+1); if(n>=m){ rep(i,m+1)res[i]=d[i+n-m]; }else{ rng(i,m-n,m+1)res[i]=d[i+n-m]; } return res; } //Yukicoder 2166 template void chmult(vc&x,const vc&y,int s){ x=multiply(move(x),y); x.resize(s); } //Poly というのは常にサイズ 1 以上であることにしよう //low のあたりをかならずサイズ s のものを返すようにいじった //その影響で何かが起きているかも知れないし,起きていないかも知れない template struct Poly:public vc{ /*Poly(){} //Poly(const vc&val):vc(val){} //Poly(vc&&val):vc(val){} explicit Poly(int n):vc(n){} explicit Poly(int n,mint v):vc(n,v){} //Poly(initializer_listinit):vc(all(init)){} template Poly(Args&&...args):vc(forward(args)...){}*/ Poly(const vc&val):vc(val){} Poly(vc&&val):vc(val){} using vc::vector; using vc::operator=; //Poly(const Poly&rhs):vc((vc)rhs){} //Poly(Poly&&rhs):vc(forward>(rhs)){} int size()const{ return vc::size(); } void ups(int s){ if(size()resize(s,0); } Poly low(int s)const{ assert(s); Poly res(s); rep(i,min(s,size()))res[i]=(*this)[i]; return res; } Poly rev()const{ auto r=*this; reverse(all(r)); return r; } Poly& operator>>=(int x){ if(x>=size()){ this->resize(1); (*this)[0]=0; }else{ this->erase(this->bg,this->bg+x); } return *this; } Poly operator>>(int x)const{return Poly(*this)>>=x;} Poly&operator<<=(int x){ this->insert(this->bg,x,0); return *this; } Poly operator<<(int x)const{return Poly(*this)<<=x;} mint freq(int i)const{ return i&r){ return operator+=((const Poly&)r); } template Poly& operator+=(T t){ ups(1); (*this)[0]+=t; return *this; } Poly& operator-=(const Poly&r){ ups(r.size()); rep(i,r.size()) (*this)[i]-=r[i]; return *this; } Poly& operator-=(const vc&r){ return operator-=((const Poly&)r); } template Poly& operator-=(T t){ ups(1); (*this)[0]-=t; return *this; } template Poly& operator*=(T t){ for(auto&v:*this) v*=t; return *this; } Poly& operator*=(const Poly&r){ *this=multiply(*this,r); return *this; } Poly& operator*=(const vc&r){ return operator*=((const Poly&)r); } Poly square()const{ return multiply(*this,*this,true); } #ifndef USE_GOOD_MOD Poly inv(int s)const{ Poly r{mint(1)/(*this)[0]}; for(int n=1;n f=low(2*n); f.resize(2*n); inplace_fmt(f,false); vc g=r.low(2*n); g.resize(2*n); inplace_fmt(g,false); rep(i,2*n)f[i]*=g[i]; inplace_fmt(f,true); rep(i,n)f[i]=0; inplace_fmt(f,false); rep(i,2*n)f[i]*=g[i]; inplace_fmt(f,true); rng(i,n,min(2*n,s))r[i]=-f[i]; } return r; } #endif template Poly& operator/=(T t){ return *this*=mint(1)/mint(t); } Poly quotient(const Poly&r,const Poly&rri)const{ int m=r.size(); assert(r[m-1].v); int n=size(); int s=n-m+1; if(s<=0) return {0}; return (rev().low(s)*rri.low(s)).low(s).rev(); } Poly& operator/=(const Poly&r){ return *this=quotient(r,r.rev().inv(max(size()-r.size(),int(0))+1)); } Poly& operator%=(const Poly&r){ *this-=*this/r*r; return *this=low(r.size()-1); } Poly operator+(const Poly&r)const{return Poly(*this)+=r;} template Poly operator+(T t)const{return Poly(*this)+=t;} template friend Poly operator+(T t,Poly r){ r[0]+=t; return r; } Poly operator-(const Poly&r)const{return Poly(*this)-=r;} template Poly operator-(T t)const{return Poly(*this)-=t;} template friend Poly operator-(T t,Poly r){ for(auto&v:r)v=-v; r[0]+=t; return r; } template Poly operator*(T t)const{return Poly(*this)*=t;} template friend Poly operator*(const T&t,const Poly&f){return f*t;} Poly operator*(const Poly&r)const{return Poly(*this)*=r;} template Poly operator/(T t)const{return Poly(*this)/=t;} Poly operator/(const Poly&r)const{return Poly(*this)/=r;} Poly operator%(const Poly&r)const{return Poly(*this)%=r;} Poly dif()const{ assert(size()); if(size()==1){ return {0}; }else{ Poly r(size()-1); rep(i,r.size()) r[i]=(*this)[i+1]*(i+1); return r; } } Poly inte(const mint invs[])const{ Poly r(size()+1,0); rep(i,size()) r[i+1]=(*this)[i]*invs[i+1]; return r; } //VERIFY: yosupo //opencupXIII GP of Peterhof H Poly log(int s,const mint invs[])const{ assert((*this)[0]==1); if(s==1)return {0}; return (low(s).dif()*inv(s-1)).low(s-1).inte(invs); } //Petrozavodsk 2019w mintay1 G //yosupo judge //UOJ Round23 C Poly exp(int s,const mint invs[])const{ assert((*this)[0]==mint(0)); Poly f{1},g{1}; for(int n=1;;n*=2){ if(n>=s)break; g=g*2-(g.square()*f).low(n); //if(n>=s)break; Poly q=low(n).dif(); q=q+g*(f.dif()-f*q).low(2*n-1); f=f+(f*(low(2*n)-q.inte(invs))).low(2*n); } return f.low(s); } //exp(x),exp(-x) のペアを返す //UOJ Round23 C pair exp2(int s,const mint invs[])const{ assert((*this)[0]==mint(0)); Poly f{1},g{1}; for(int n=1;;n*=2){ //if(n>=s)break; g=g*2-(g.square()*f).low(n); if(n>=s)break; Poly q=low(n).dif(); q=q+g*(f.dif()-f*q).low(2*n-1); f=f+(f*(low(2*n)-q.inte(invs))).low(2*n); } return make_pair(f.low(s),g.low(s)); } #ifndef USE_GOOD_MOD //CF250 E Poly sqrt(int s)const{ assert((*this)[0]==1); static const mint half=mint(1)/mint(2); Poly r{1}; for(int n=1;n f{1},g{1},z{1}; for(int n=1;n delta(2*n); rep(i,n)delta[n+i]=z[i]-freq(i)-freq(n+i); inplace_fmt(delta,false); vc gbuf(2*n); rep(i,n)gbuf[i]=g[i]; inplace_fmt(gbuf,false); rep(i,2*n)delta[i]*=gbuf[i]; inplace_fmt(delta,true); f.resize(2*n); rng(i,n,2*n)f[i]=-half*delta[i]; if(2*n>=s)break; z=f; inplace_fmt(z,false); vc eps=gbuf; rep(i,2*n)eps[i]*=z[i]; inplace_fmt(eps,true); rep(i,n)eps[i]=0; inplace_fmt(eps,false); rep(i,2*n)eps[i]*=gbuf[i]; inplace_fmt(eps,true); g.resize(2*n); rng(i,n,2*n)g[i]=-eps[i]; } f.resize(s); return f; } #endif pair divide(const Poly&r,const Poly&rri)const{ Poly a=quotient(r,rri); Poly b=*this-a*r; return make_pair(a,b.low(r.size()-1)); } //Yukicoder No.215 Poly pow_mod(int n,const Poly&r)const{ Poly rri=r.rev().inv(r.size()); Poly cur{1},x=*this%r; while(n){ if(n%2) cur=(cur*x).divide(r,rri).b; x=(x*x).divide(r,rri).b; n/=2; } return cur; } int lowzero()const{ rep(i,size())if((*this)[i]!=0)return i; return size(); } //VERIFY: yosupo //UOJ Round23 C (z=0,p<0) //Multiuni2023-4-B Poly pow(int s,int p,const mint invs[])const{ assert(s>0); if(p==0){ Poly res(s,0); res[0]=1; return res; } int n=size(),z=0; for(;z0); if(z*p>=s)return Poly(s,0); mint c=(*this)[z],cinv=c.inv(); mint d=c.pow(p); int t=s-z*p; Poly x(t); rng(i,z,min(z+t,n))x[i-z]=(*this)[i]*cinv; x=x.log(t,invs); rep(i,t)x[i]*=p; x=x.exp(t,invs); rep(i,t)x[i]*=d; Poly y(s); rep(i,t)y[z*p+i]=x[i]; return y; } //verified in compositional inverse Poly pow_const1(mint p,const mint invs[])const{ assert((*this)[0]==1); Poly x=log(size(),invs); rep(i,size())x[i]*=p; return x.exp(size(),invs); } mint eval(mint x)const{ mint r=0,w=1; for(auto v:*this){ r+=w*v; w*=x; } return r; } }; template void print_single(const Poly&v,int suc=1){ rep(i,v.size()) print_single(v[i],i==int(v.size())-1?3:2); if(suc==1){ if(dbg)cout< struct Laurent{ Poly f; int a; Laurent(const Poly&num,const Poly&den,int s){ a=den.lowzero(); assert(a>a).inv(s)).low(s); } Laurent(const Poly&ff,int aa):f(ff),a(aa){} Laurent dif()const{ return Laurent(f*(-a)+(f.dif()<<1),a+1); } mint&operator[](int i){ assert(inc(0,i+a,si(f)-1)); return f[i+a]; } }; #ifndef DYNAMIC_MOD extern constexpr modinfo base{998244353,3}; //extern constexpr modinfo base{1000000007,0}; //extern constexpr modinfo base{2147483579,1689685080};//2^31 未満の最大の安全素数 //modinfo base{1,0}; #ifdef USE_GOOD_MOD static_assert(base.mod==998244353); #endif #else modinfo base(1,0); extern constexpr modinfo base107(1000000007,0); using mint107=modular; #endif using mint=modular; mint parity(int i){ return i%2==0?1:-1; } #ifdef LOCAL const int vmax=10010; #else const int vmax=(1<<21)+10; #endif mint fact[vmax],finv[vmax],invs[vmax]; void initfact(){ const int s=min(vmax,base.mod); fact[0]=1; rng(i,1,s){ fact[i]=fact[i-1]*i; } finv[s-1]=fact[s-1].inv(); for(int i=s-2;i>=0;i--){ finv[i]=finv[i+1]*(i+1); } for(int i=s-1;i>=1;i--){ invs[i]=finv[i]*fact[i-1]; } } mint choose(int n,int k){ return inc(0,k,n)?fact[n]*finv[n-k]*finv[k]:0; } mint binom(int a,int b){ return 0<=a&&0<=b?fact[a+b]*finv[a]*finv[b]:0; } mint catalan(int n){ return binom(n,n)-(n-1>=0?binom(n-1,n+1):0); } //対角線を超えず (x,y) に至る方法の数 mint catalan(int x,int y){ assert(y<=x); return binom(x,y)-binom(x+1,y-1); } //y=x+c を超えず (x,y) に至る方法の数 mint catalan(int x,int y,int c){ assert(y<=x+c); return binom(x,y)-binom(x+c+1,y-c-1); } /* const int vmax=610; mint fact[vmax+1],binbuf[vmax+1][vmax+1]; mint choose(int n,int k){ return 0<=k&&k<=n?binbuf[n-k][k]:0; } mint binom(int a,int b){ return 0<=a&&0<=b?binbuf[a][b]:0; } void initfact(int n){ fact[0]=1; rep(i,n)fact[i+1]=fact[i]*(i+1); rep(i,n+1)rep(j,n+1){ if(i==0&&j==0){ binbuf[i][j]=1; }else{ binbuf[i][j]=0; if(i)binbuf[i][j]+=binbuf[i-1][j]; if(j)binbuf[i][j]+=binbuf[i][j-1]; } } } */ mint p2[vmax],p2inv[vmax]; void initp2(){ p2[0]=1; rep(i,vmax-1)p2[i+1]=p2[i]*2; p2inv[vmax-1]=p2[vmax-1].inv(); per(i,vmax-1)p2inv[i]=p2inv[i+1]*2; } int m2i(mint a){ uint v=a.v; return vabs(a)+b){ res=pi(a,b); } }; rng(b,1,1000){ mint y=x*b; upd(y.v,b); upd(-int((-y).v),b); } return res; } ll m2l(mint a){ return a.v&a){ vi tmp(si(a)); rep(i,si(a)){ tmp[i]=m2l(a[i]); } cerr< a,Poly b,int k){ int s=1;while(s c(s*2); while(s<=k){ fill(all(c),mint(0)); rep(i,s){ c[i]=b[i]; if(i&1)c[i]=-c[i]; } inplace_fmt(a,false); inplace_fmt(b,false); inplace_fmt(c,false); rep(i,2*s)a[i]*=c[i]; rep(i,2*s)b[i]*=c[i]; inplace_fmt(a,true); inplace_fmt(b,true); if(k%2==0){ rep(i,s)a[i]=a[i*2]; rng(i,s,2*s)a[i]=0; }else{ rep(i,s)a[i]=a[i*2+1]; rng(i,s,2*s)a[i]=0; } rep(i,s)b[i]=b[i*2]; rng(i,s,s*2)b[i]=0; k/=2; } a.resize(s); b.resize(s); auto d=a*b.inv(s); return d[k]; } using P=Poly; P extend(const P&a,int k){ assert(a[0]==1); int n=si(a)-1,s=n*k; auto b=a.log(s+1,invs); rep(i,s+1){ if(i%k==0)b[i]*=k; else b[i]*=0; } auto c=b.exp(s+1,invs); P res(n+1); rep(i,n+1)res[i]=c[i*k]; return res; } P magnify(const P&a,int k){ int n=si(a)-1; P res(n*k+1); rep(i,n+1)res[i*k]=a[i]; return res; } P select_mod(const P&a,int k,int tar){ P res; for(int i=tar;i ns{{N}}; while(1){ int v=ns.back()[0]/2; if(v==0)break; ns.pb(vi{v}); } for(auto&ls:ns){ while(1){ int v=ls.back()/3; if(v==0)break; ls.pb(v); } } vvc

fs(si(ns)); rep(i,si(ns))fs[i].resize(si(ns[i])); rep(i,si(fs))rep(j,si(fs[i])){ if(i==0&&j==0){ fs[i][j]=P{1,-1,-1,-1}; }else if(i>0){ fs[i][j]=extend(fs[i-1][j],2); }else if(j>0){ fs[i][j]=extend(fs[i][j-1],3); } } auto den=fs; rep(i,si(den))rep(j,si(den[i])){ if(i>0)den[i][j]*=den[i-1][j]; } rep(i,si(den))rep(j,si(den[i])){ if(j>0)den[i][j]*=den[i][j-1]; } vvc

dp(si(ns)); rep(i,si(ns))dp[i].resize(si(ns[i])); rep(i,si(dp))rep(j,si(dp[i])){ if(i==0&&j==0)dp[i][j]=P{1}; else{ if(i>0){ auto f=dp[i-1][j]*(magnify(den[i][j],2)/den[i-1][j]); dp[i][j]+=select_mod(f,2,ns[i-1][j]%2); } if(j>0){ auto f=dp[i][j-1]*(magnify(den[i][j],3)/den[i][j-1]); dp[i][j]+=select_mod(f,3,ns[i][j-1]%3); } dp[i][j]/=P{1,-1,-1,-1}; } } mint ans=0; rep(i,si(dp))rep(j,si(dp[i])){ ans+=get_kth_div_fast(dp[i][j],den[i][j],ns[i][j]-1); } print(ans); } signed main(signed argc,char*argv[]){ if(argc>1&&strcmp(argv[1],"D")==0)dbg=true; cin.tie(0); ios::sync_with_stdio(0); cout<