結果
問題 |
No.1145 Sums of Powers
|
ユーザー |
![]() |
提出日時 | 2025-07-29 15:27:11 |
言語 | C++23 (gcc 13.3.0 + boost 1.87.0) |
結果 |
AC
|
実行時間 | 247 ms / 2,000 ms |
コード長 | 18,121 bytes |
コンパイル時間 | 1,626 ms |
コンパイル使用メモリ | 120,224 KB |
実行使用メモリ | 8,940 KB |
最終ジャッジ日時 | 2025-07-29 15:27:15 |
合計ジャッジ時間 | 3,467 ms |
ジャッジサーバーID (参考情報) |
judge5 / judge1 |
(要ログイン)
ファイルパターン | 結果 |
---|---|
other | AC * 6 |
ソースコード
#include<iostream> #include<vector> #include<cassert> #include<optional> #include<algorithm> #include<array> #include<type_traits> #include<numeric> template<typename T> constexpr std::enable_if_t<(std::numeric_limits<T>::digits<=32),T>pow_mod(T a,T n,T mod){ using u64=unsigned long long; u64 res=1; while(n>0){ if(n&1)res=((u64)res*a)%mod; a=((u64)a*a)%mod; n>>=1; } return T(res); } template<typename T> constexpr std::enable_if_t<(std::numeric_limits<T>::digits>32),T>pow_mod(T a,T n,T mod){ using u128=__uint128_t; u128 res=1; while(n>0){ if(n&1)res=((u128)res*a)%mod; a=((u128)a*a)%mod; n>>=1; } return T(res); } constexpr int primitive_root_constexpr(int x){ if(x==167772161)return 3; if(x==469762049)return 3; if(x==754974721)return 11; if(x==880803841)return 26; if(x==998244353)return 3; if(x==2)return 1; int x2=x; int p[20]={}; int c=0; x--; for(int i=2;i*i<=x;i++){ if(x%i==0){ p[c++]=i; while(x%i==0)x/=i; } } if(x!=1)p[c++]=x; x=x2; for(int g=2;;g++){ bool ok=true; for(int i=0;i<c;i++)if(pow_mod(g,(x-1)/p[i],x)==1){ ok=false; break; } if(ok)return g; } } #include<limits> #include<concepts> template<typename T> constexpr std::enable_if_t<std::numeric_limits<T>::digits<=32,int>msb(T n){return n==0?-1:31-__builtin_clz(n);} template<typename T> constexpr std::enable_if_t<(std::numeric_limits<T>::digits>32),int>msb(T n){return n==0?-1:63-__builtin_clzll(n);} template<typename T> constexpr std::enable_if_t<std::numeric_limits<T>::digits<=32,int>lsb(T n){return n==0?-1:__builtin_ctz(n);} template<typename T> constexpr std::enable_if_t<(std::numeric_limits<T>::digits>32),int>lsb(T n){return n==0?-1:__builtin_ctzll(n);} template<typename T> constexpr std::enable_if_t<std::is_integral_v<T>,T>floor_pow2(T n){return n==0?0:T(1)<<msb(n);} template<typename T> constexpr std::enable_if_t<std::is_integral_v<T>,T>ceil_pow2(T n){return n<=1?1:T(1)<<(msb(n-1)+1);} template<std::integral T> constexpr T safe_div(T a,T b){return a/b-(a%b&&(a^b)<0);} template<std::integral T> constexpr T safe_ceil(T a,T b){return a/b+(a%b&&(a^b)>0);} template<int m> struct ntt_root{ static constexpr int rank2=lsb(m-1); static constexpr int g=primitive_root_constexpr(m); std::array<int,rank2+1>root,invroot; std::array<int,std::max(0,rank2-1)>rate2,invrate2; std::array<int,std::max(0,rank2-2)>rate3,invrate3; constexpr ntt_root(){ root[rank2]=pow_mod(g,m>>rank2,m); invroot[rank2]=pow_mod(root[rank2],m-2,m); for(int i=rank2-1;i>=0;i--){ root[i]=(long long)root[i+1]*root[i+1]%m; invroot[i]=(long long)invroot[i+1]*invroot[i+1]%m; } int prod=1,invprod=1; for(int i=0;i<rank2-1;i++){ rate2[i]=(long long)root[i+2]*prod%m; invrate2[i]=(long long)invroot[i+2]*invprod%m; prod=(long long)prod*invroot[i+2]%m; invprod=(long long)invprod*root[i+2]%m; } prod=invprod=1; for(int i=0;i<rank2-2;i++){ rate3[i]=(long long)root[i+3]*prod%m; invrate3[i]=(long long)invroot[i+3]*invprod%m; prod=(long long)prod*invroot[i+3]%m; invprod=(long long)invprod*root[i+3]%m; } } }; template<typename T> void dft(std::vector<T>&a){ static constexpr ntt_root<T::mod()>r; static constexpr unsigned long long mod2=(unsigned long long)T::mod()*T::mod(); int n=a.size(); int h=lsb(n); int len=0; while(len<h){ if(h-len==1){ T rot=T::raw(1); for(int s=0;s<(1<<len);s++){ int of=s*2; T u=a[of],v=a[of+1]*rot; a[of]=u+v; a[of+1]=u-v; rot*=T::raw(r.rate2[lsb(~(unsigned int)s)]); } len++; } else{ int p=1<<(h-len-2); T rot=T::raw(1),imag=T::raw(r.root[2]); for(int s=0;s<(1<<len);s++){ const unsigned long long rot1=rot.val(),rot2=rot1*rot1%T::mod(),rot3=rot1*rot2%T::mod(); int of=s<<(h-len); for(int i=0;i<p;i++){ const unsigned long long a0=a[i+of].val(),a1=(unsigned long long)a[i+of+p].val()*rot1,a2=(unsigned long long)a[i+of+p*2].val()*rot2,a3=(unsigned long long)a[i+of+p*3].val()*rot3; const unsigned long long m=(unsigned long long)T(a1+mod2-a3).val()*imag.val(); const unsigned long long k=mod2-a2; a[i+of]=a0+a2+a1+a3; a[i+of+p]=a0+a2+(mod2*2-a1-a3); a[i+of+p*2]=a0+k+m; a[i+of+p*3]=a0+k+(mod2-m); } rot*=T::raw(r.rate3[lsb(~(unsigned int)s)]); } len+=2; } } } template<typename T> void idft(std::vector<T>&a){ static constexpr ntt_root<T::mod()>r; int n=a.size(); int h=lsb(n); int len=h; while(len){ if(len==1){ int p=1<<(h-1); for(int i=0;i<p;i++){ T u=a[i],v=a[i+p]; a[i]=u+v; a[i+p]=u-v; } len--; } else{ int p=1<<(h-len); T rot=T::raw(1),imag=T::raw(r.invroot[2]); for(int s=0;s<(1<<(len-2));s++){ const unsigned long long rot1=rot.val(),rot2=rot1*rot1%T::mod(),rot3=rot1*rot2%T::mod(); int of=s<<(h-len+2); for(int i=0;i<p;i++){ const unsigned long long a0=a[i+of].val(),a1=a[i+of+p].val(),a2=a[i+of+p*2].val(),a3=a[i+of+p*3].val(); const unsigned long long k=T((T::mod()+a2-a3)*imag.val()).val(); a[i+of]=a0+a1+a2+a3; a[i+of+p]=(a0+T::mod()-a1+k)*rot1; a[i+of+p*2]=(a0+a1+T::mod()*2-a2-a3)*rot2; a[i+of+p*3]=(a0+T::mod()*2-a1-k)*rot3; } rot*=T::raw(r.invrate3[lsb(~(unsigned int)s)]); } len-=2; } } } template<typename T> std::vector<T>ntt_convolution(std::vector<T> a,std::vector<T> b){ int n=a.size(),m=b.size(),s=n+m-1; if(std::min(n,m)<60){ std::vector<T>ret(s,0); if(n<m)for(int i=0;i<m;i++)for(int j=0;j<n;j++)ret[i+j]+=a[j]*b[i]; else for(int i=0;i<n;i++)for(int j=0;j<m;j++)ret[i+j]+=a[i]*b[j]; return ret; } int z=ceil_pow2(s); a.resize(z,0); b.resize(z,0); dft(a),dft(b); std::vector<T>c(z); for(int i=0;i<z;i++)c[i]=a[i]*b[i]; idft(c); T g=T::raw(z).inv(); for(int i=0;i<s;i++)c[i]*=g; return {c.begin(),c.begin()+s}; } template<typename T> std::vector<std::vector<T>>convolution2d(const std::vector<std::vector<T>>&a,const std::vector<std::vector<T>>&b){ int x=a.size()+b.size()-1,y=a[0].size()+b[0].size()-1; std::vector<T>f(a.size()*y),g(b.size()*y); for(int i=0;i<a.size();i++)for(int j=0;j<a[i].size();j++)f[i*y+j]=a[i][j]; for(int i=0;i<b.size();i++)for(int j=0;j<b[i].size();j++)g[i*y+j]=b[i][j]; auto c=ntt_convolution(f,g); std::vector<std::vector<T>>ret(x,std::vector<T>(y)); for(int i=0;i<x;i++)for(int j=0;j<y;j++)ret[i][j]=c[i*y+j]; return ret; } template<typename T> std::vector<T> fps_inv(const std::vector<T> &a,int deg=-1){ int n=a.size(); if(deg==-1)deg=n; const T zero=T::raw(0); assert(a[0]!=zero); std::vector<T> ret(ceil_pow2(deg)); ret[0]=a[0].inv(); for(int m=1;m<deg;m<<=1){ std::vector<T> f(a.begin(),a.begin()+std::min(n,m*2)); if(f.size()<m*2)f.resize(m*2,0); std::vector<T> g(ret); f.resize(m*2); dft(f); g.resize(m*2); dft(g); for(int i=0;i<m*2;i++)f[i]*=g[i]; idft(f); T inv=T::raw(m*2).inv(); for(int i=0;i<m;i++)f[i]=zero; for(int i=m;i<m*2;i++)f[i]*=inv; dft(f); for(int i=0;i<m*2;i++)f[i]*=g[i]; idft(f); for(int i=0;i<m*2;i++)f[i]*=inv; for(int i=m;i<m*2;i++)ret[i]-=f[i]; } ret.resize(deg); return ret; } template<typename T> std::vector<T> fps_diff(const std::vector<T>&a){ int n=a.size(); std::vector<T>b(std::max(0,n-1)); for(int i=1;i<n;i++)b[i-1]=a[i]*T::raw(i); return b; } template<typename T> std::vector<T> fps_integral(const std::vector<T>&a){ int n=a.size(); std::vector<T>b(n+1); b[0]=0; if(n)b[1]=1; for(int i=2;i<=n;i++)b[i]=-b[T::mod()%i]*(T::mod()/i); for(int i=0;i<n;i++)b[i+1]*=a[i]; return b; } template<typename T> std::vector<T> fps_log(const std::vector<T>&a,int deg=-1){ int n=a.size(); if(deg==-1)deg=n; std::vector<T>b=fps_integral(ntt_convolution(fps_diff(a),fps_inv(a,deg))); return {b.begin(),b.begin()+deg}; } template<typename T> std::vector<T> fps_exp(const std::vector<T>&a,int deg=-1){ assert(a.empty()||a[0]==0); int n=a.size(); if(deg==-1)deg=n; std::vector<T>inv(ceil_pow2(deg)); inv[0]=0,inv[1]=1; for(int i=2;i<inv.size();i++)inv[i]=-inv[T::mod()%i]*(T::mod()/i); std::vector<T>b{1,1<n?a[1]:T()},c{T::raw(1)},f1,f2(2,T::raw(1)); T sinv=T::raw(2).inv(); for(int m=2;m<deg;m*=2){ std::vector<T>y(b); y.resize(m*2); dft(y); f1=f2; std::vector<T>z(m); for(int i=0;i<m;i++)z[i]=y[i]*f1[i]; idft(z); for(int i=m/2;i<m;i++)z[i]*=sinv; std::fill(z.begin(),z.begin()+m/2,T()); dft(z); for(int i=0;i<m;i++)z[i]*=-f1[i]; idft(z); for(int i=m/2;i<m;i++)z[i]*=sinv; c.insert(c.end(),z.begin()+m/2,z.end()); f2=c; f2.resize(m*2); dft(f2); std::vector<T>x(m); std::copy(a.begin(),a.begin()+std::min(n,m),x.begin()); for(int i=0;i<m-1;i++)x[i]=x[i+1]*T::raw(i+1); x.back()=T(); dft(x); for(int i=0;i<m;i++)x[i]*=y[i]; idft(x); for(int i=0;i<m;i++)x[i]*=sinv; x.resize(m*2); for(int i=0;i<m-1;i++)x[m+i]=x[i]-b[i+1]*T::raw(i+1),x[i]=T(); dft(x); for(int i=0;i<m*2;i++)x[i]*=f2[i]; idft(x); sinv/=2; for(int i=0;i<m*2;i++)x[i]*=sinv; for(int i=(int)x.size()-1;i>=1;i--)x[i]=x[i-1]*inv[i]; x[0]=T(); for(int i=m;i<std::min(n,m*2);i++)x[i]+=a[i]; std::fill(x.begin(),x.begin()+m,T()); dft(x); for(int i=0;i<m*2;i++)x[i]*=y[i]; idft(x); for(int i=m;i<m*2;i++)x[i]*=sinv; b.insert(b.end(),x.begin()+m,x.end()); } b.resize(deg); return b; } template<typename T> std::vector<T> fps_pow(std::vector<T> a,unsigned long long k,int deg=-1){ int n=a.size(); if(deg==-1)deg=n; if(k==0){ std::vector<T>ret(deg,0); ret[0]=1; return ret; } int of=0; while(a[of]==0&&of!=n)of++; if(of==n)return std::vector<T>(deg,0); if(of!=0&&k>=deg)return std::vector<T>(deg,0); if(of*k>=deg)return std::vector<T>(deg,0); a.erase(a.begin(),a.begin()+of); n=a.size(); T a0=a[0]; T inv=a[0].inv(); for(int i=0;i<n;i++)a[i]*=inv; std::vector<T>lg=fps_log(a,deg); T tk=T(k); for(int i=0;i<deg;i++)lg[i]*=tk; std::vector<T>ep=fps_exp(lg,deg); T pw=a0.pow(k); for(int i=0;i<deg;i++)ep[i]*=pw; std::vector<T>ret(deg,0); for(int i=of*k;i<deg;i++)ret[i]=ep[i-of*k]; return ret; } template<typename T> std::vector<T>fps_divmod(std::vector<T>f,std::vector<T>g){ if(f.empty())return f; int n=f.size(),m=g.size(); std::vector<T>r(f); std::reverse(f.begin(),f.end()),std::reverse(g.begin(),g.end()); std::vector<T>q=ntt_convolution(f,fps_inv(g,n)); q.resize(std::max(0,n-m+1)); auto p=ntt_convolution(g,q); std::reverse(p.begin(),p.end()); for(int i=0;i<std::min<int>(p.size(),m);i++)r[i]-=p[i]; r.resize(m-1); while(!r.empty()&&r.back().val()==0)r.pop_back(); return r; } template<typename T> std::optional<std::vector<T>>fps_sqrt(std::vector<T>f,int deg=-1){ if(deg==-1)deg=f.size(); int prefix_zero=0; while(prefix_zero<f.size()&&f[prefix_zero].val()==0)prefix_zero++; if(prefix_zero==f.size())return std::make_optional(std::vector<T>(deg,0)); if(prefix_zero&1)return std::nullopt; f.erase(f.begin(),f.begin()+prefix_zero); prefix_zero/=2; auto opt_sq=f[0].sqrt(); if(!opt_sq)return std::nullopt; T sq=*opt_sq; T inv0=f[0].inv(); for(int i=0;i<f.size();i++)f[i]*=inv0; std::vector<T>g{1}; T inv2=T::raw(2).inv(); while(g.size()<deg-prefix_zero){ std::vector<T>fp(f.begin(),f.begin()+std::min(f.size(),g.size()*2)); fp=ntt_convolution(fp,fps_inv(g,g.size()*2)); std::vector<T>nxtg(g); nxtg.resize(g.size()*2); for(int i=0;i<nxtg.size();i++){ if(i<fp.size())nxtg[i]+=fp[i]; nxtg[i]*=inv2; } std::swap(g,nxtg); } g.insert(g.begin(),prefix_zero,0); g.resize(deg); for(int i=0;i<g.size();i++)g[i]*=sq; return g; } template<typename T> void ntt_doubling(std::vector<T>&a){ static constexpr ntt_root<T::mod()>r; int n=a.size()/2; std::vector<T>b(a.begin(),a.begin()+n); idft(b); T now=T::raw(n).inv(),zeta=T::raw(r.root[msb(n)+1]); for(int i=0;i<n;i++){ b[i]*=now; now*=zeta; } dft(b); std::copy(b.begin(),b.end(),a.begin()+n); } //0≤j<mに対するsum(f_i*g_i^j)の列挙 template<typename T> std::vector<T>power_prod_sum(std::vector<T>f,std::vector<T>g,int m){ assert(f.size()==g.size()); int s=f.size(); int n=ceil_pow2(f.size()),log2n=msb(n); for(int i=0;i<s;i++)g[i]=-g[i]; f.resize(n),g.resize(n); for(int i=0;i<log2n;i++){ int b=1<<i; std::vector<T>lf(b*2),lg(b*2),rf(b*2),rg(b*2); int j; for(j=0;j+b<s;j+=b*2){ std::copy(f.begin()+j,f.begin()+j+b,lf.begin()); std::copy(f.begin()+j+b,f.begin()+j+b*2,rf.begin()); std::copy(g.begin()+j,g.begin()+j+b,lg.begin()); std::copy(g.begin()+j+b,g.begin()+j+b*2,rg.begin()); ntt_doubling(lf),ntt_doubling(rf),ntt_doubling(lg),ntt_doubling(rg); for(int k=0;k<b;k++)lg[k]++,rg[k]++; for(int k=b;k<b*2;k++)lg[k]--,rg[k]--; for(int k=0;k<b*2;k++){ f[j+k]=lf[k]*rg[k]+lg[k]*rf[k]; g[j+k]=lg[k]*rg[k]; g[j+k]--; } } if(j<s){ std::copy(f.begin()+j,f.begin()+j+b,lf.begin()); std::copy(g.begin()+j,g.begin()+j+b,lg.begin()); ntt_doubling(lf),ntt_doubling(lg); for(int k=0;k<b;k++){ f[j+k]=lf[k]; g[j+k]=lg[k]; } for(int k=b;k<b*2;k++){ f[j+k]=-lf[k]; g[j+k]=-lg[k]; } } } idft(f),idft(g); g.push_back(n); std::reverse(f.begin(),f.end()); std::reverse(g.begin(),g.end()); f.resize(m),g.resize(m); f=ntt_convolution(f,fps_inv(g)); f.resize(m); return f; } constexpr int carmichael_constexpr(int n){ if(n==998244353)return 998244352; if(n==1000000007)return 1000000006; if(n<=1)return n; int res=1; int t=0; while(n%2==0){ n/=2; t++; } if(t==2)res=2; else if(t>=3)res=1<<(t-2); for(int i=3;i*i<=n;i++)if(n%i==0){ int c=0; while(n%i==0){ n/=i; c++; } int prod=i-1; for(int j=0;j<c-1;j++)prod*=i; res=std::lcm(res,prod); } if(n!=1)res=std::lcm(res,n-1); return res; } template<int m> struct mod_int{ private: static constexpr unsigned int umod=static_cast<unsigned int>(m); static constexpr unsigned int car=carmichael_constexpr(m); using uint=unsigned int; using mint=mod_int; uint v; static_assert(m<uint(1)<<31); mint sqrt_impl()const{ if(this->val()<=1)return *this; if constexpr(m%8==1){ mint b=2; while(b.pow((m-1)/2).val()==1)b++; int m2=m-1,e=0; while(m2%2==0)m2>>=1,e++; mint x=this->pow((m2-1)/2); mint y=(*this)*x*x; x*=*this; mint z=b.pow(m2); while(y.val()!=1){ int j=0; mint t=y; while(t.val()!=1)t*=t,j++; z=z.pow(1<<(e-j-1)); x*=z; z*=z; y*=z;e=j; } return x; } else if constexpr(m%8==5){ mint ret=this->pow((m+3)/8); if((ret*ret).val()==this->val())return ret; else return ret*mint::raw(2).pow((m-1)/4); } else{ return this->pow((m+1)/4); } } public: using value_type=uint; mod_int():v(0){} template<typename T,std::enable_if_t<std::is_signed_v<T>,std::nullptr_t> =nullptr> mod_int(T a){ a%=m; if(a<0)v=a+umod; else v=a; } template<typename T,std::enable_if_t<std::is_unsigned_v<T>,std::nullptr_t> =nullptr> mod_int(T a):v(a%umod){} static constexpr mint raw(int a){ mint ret; ret.v=a; return ret; } inline uint val()const{return this->v;} static constexpr int mod(){return m;} inline mint &operator+=(const mint &b){ this->v+=b.v; if(this->v>=umod)this->v-=umod; return *this; } inline mint &operator-=(const mint &b){ this->v-=b.v; if(this->v>=umod)this->v+=umod; return *this; } inline mint &operator*=(const mint &b){ this->v=((unsigned long long)this->v*b.v)%umod; return *this; } inline mint &operator/=(const mint &b){ *this*=b.inv(); return *this; } inline mint operator+()const{return *this;} inline mint operator-()const{return mint()-*this;} friend inline mint operator+(const mint &a,const mint &b){return mint(a)+=b;} friend inline mint operator-(const mint &a,const mint &b){return mint(a)-=b;} friend inline mint operator*(const mint &a,const mint &b){return mint(a)*=b;} friend inline mint operator/(const mint &a,const mint &b){return mint(a)/=b;} friend inline bool operator==(const mint &a,const mint &b){return a.val()==b.val();} friend inline bool operator!=(const mint &a,const mint &b){return !(a==b);} inline mint operator++(int){ mint ret=*this; *this+=mint::raw(1); return ret; } inline mint operator--(int){ mint ret=*this; *this-=mint::raw(1); return ret; } mint pow(long long n)const{ mint ret=mint::raw(1),a(*this); while(n){ if(n&1)ret*=a; a*=a; n>>=1; } return ret; } inline mint inv()const{ assert(this->v!=0); return pow(car-1); } std::optional<mint>sqrt()const{ if(this->val()<=1||this->pow((m-1)/2)==1)return std::make_optional(this->sqrt_impl()); else return std::nullopt; } static constexpr unsigned int order(){return car;} friend std::istream &operator>>(std::istream &is,mint &b){ long long a; is>>a; b=mint(a); return is; } friend std::ostream &operator<<(std::ostream &os,const mint &b){ os<<b.val(); return os; } }; template<int m> struct std::hash<mod_int<m>>{ std::size_t operator()(mod_int<m>x)const{ return std::hash<unsigned int>()(x.val()); } }; using mint998=mod_int<998244353>; using mint107=mod_int<1000000007>; using mint=mint998; using namespace std; int main(){ int n,m; cin>>n>>m; vector<mint>a(n),b(n); for(int i=0;i<n;i++){ a[i]=1; cin>>b[i]; } a=power_prod_sum(a,b,m+1); for(int i=1;i<=m;i++)cout<<a[i]<<" \n"[i==m]; }