src/math/convolution.hpp

• View this file on GitHub
• Last update: 2025-12-30 15:03:15+09:00
• Include: #include "src/math/convolution.hpp"

Depends on

• src/math/modint.hpp
• src/template.hpp

Code

#pragma once
#include "./modint.hpp"

//acl
mint g=3;
void fft(vector<mint>& a,bool inv=false){
    int n=si(a),s=__lg(n);
    static vector<mint> z,iz;
    while(si(z)<=s){
        z.emplace_back(g.pow(mint(-1).val()/(1<<si(z))));
        iz.emplace_back(z.back().inv());
    }
    vector<mint> b(n);
    for(int i=1;i<=s;i++){
        int w=1<<s-i;
        mint base=inv?iz[i]:z[i],now=1;
        for(int y=0;y<n/2;y+=w){
            rep(x,w){
                auto l=a[y<<1|x],r=now*a[y<<1|x|w];
                b[y|x]=l+r,b[y|x|n>>1]=l-r;
            }
            now*=base;
        }
        swap(a,b);
    }
}

vector<mint> convolution(vector<mint> a,vector<mint> b){
    int n=si(a),m=si(b);
    if(!n or !m) return {};
    if(min(n,m)<=30){
        vector<mint> ans(n+m-1);
        rep(i,n) rep(j,m) ans[i+j]+=a[i]*b[j];
        return ans;
    }
    int N=n+m-1;
    int z=__bit_ceil(unsigned(N));
    a.resize(z),b.resize(z);
    fft(a),fft(b);
    rep(i,z) a[i]*=b[i];
    fft(a,true);
    a.resize(n+m-1);
    mint iz=mint(z).inv();
    for(auto&& e : a) e*=iz;
    return a;
}

#line 2 "src/template.hpp"
#include <bits/stdc++.h>
using namespace std;
#define si(a) (long)a.size()
#define fi first
#define se second
#define all(x) x.begin(),x.end()
#define rep(i,n) for(int i=0;i<(int)(n);++i)
template<typename S,typename F> bool chmin(S&a,F b){return b<a?(a=b,1):0;}
template<typename S,typename F> bool chmax(S&a,F b){return b>a?(a=b,1):0;}
bool _=(ios::sync_with_stdio(0),cin.tie(0),cout<<fixed<<setprecision(16),0);
#line 3 "src/math/modint.hpp"

template<int mod=998244353> struct modint {
    int x;
    constexpr modint(long long x_=0):x(((x_%mod)+mod)%mod){}
    constexpr modint operator-(){
        auto res=*this;
        res.x=(x?mod-x:0);
        return res;
    }
    constexpr modint& operator+=(modint r){
        if((x+=r.x)>=mod) x-=mod;
        return *this;
    }
    constexpr modint& operator-=(modint r){
        if((x-=r.x)<0) x+=mod;
        return *this;
    }
    constexpr modint& operator*=(modint r){
        x=1ll*x*r.x%mod;
        return *this;
    }
    constexpr modint& operator/=(modint r){return *this*=r.inv();}
    constexpr friend modint operator+(modint a,modint b){return a+=b;}
    constexpr friend modint operator-(modint a,modint b){return a-=b;}
    constexpr friend modint operator*(modint a,modint b){return a*=b;}
    constexpr friend modint operator/(modint a,modint b){return a/=b;}
    constexpr modint inv() const {return pow(mod-2);}
    constexpr modint pow(long long b) const {
        assert(0<=b);
        modint a=*this,c=1;
        while(b){
            if(b&1) c*=a;
            a*=a;
            b>>=1;
        }
        return c;
    }
    constexpr int val() const {return x;}
    constexpr friend ostream& operator<<(ostream& os,const modint& m){return os<<m.val();}
    constexpr friend istream& operator>>(istream& is,modint& m){
        long long v;
        is>>v;
        m=modint(v);
        return is;
    }
};
#line 3 "src/math/convolution.hpp"

//acl
mint g=3;
void fft(vector<mint>& a,bool inv=false){
    int n=si(a),s=__lg(n);
    static vector<mint> z,iz;
    while(si(z)<=s){
        z.emplace_back(g.pow(mint(-1).val()/(1<<si(z))));
        iz.emplace_back(z.back().inv());
    }
    vector<mint> b(n);
    for(int i=1;i<=s;i++){
        int w=1<<s-i;
        mint base=inv?iz[i]:z[i],now=1;
        for(int y=0;y<n/2;y+=w){
            rep(x,w){
                auto l=a[y<<1|x],r=now*a[y<<1|x|w];
                b[y|x]=l+r,b[y|x|n>>1]=l-r;
            }
            now*=base;
        }
        swap(a,b);
    }
}

vector<mint> convolution(vector<mint> a,vector<mint> b){
    int n=si(a),m=si(b);
    if(!n or !m) return {};
    if(min(n,m)<=30){
        vector<mint> ans(n+m-1);
        rep(i,n) rep(j,m) ans[i+j]+=a[i]*b[j];
        return ans;
    }
    int N=n+m-1;
    int z=__bit_ceil(unsigned(N));
    a.resize(z),b.resize(z);
    fft(a),fft(b);
    rep(i,z) a[i]*=b[i];
    fft(a,true);
    a.resize(n+m-1);
    mint iz=mint(z).inv();
    for(auto&& e : a) e*=iz;
    return a;
}

Back to top page