提交记录 18000


用户 题目 状态 得分 用时 内存 语言 代码长度
cuking 1002. 测测你的多项式乘法 Accepted 100 251.493 ms 101844 KB C++11 5.58 KB
提交时间 评测时间
2022-09-07 21:55:16 2022-09-07 21:55:18
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define U unsigned
#define LL long long
#define UL U LL
#define S int
#define fr(i,l) for(S i=0;i<l;i++)
using std::min;
using std::max;
using std::swap;
constexpr U mod=998244353u;
U pow(U a,U b)
{
    U ans=1;
    for(;b;a=(UL)a*a%mod,b>>=1)
        if(b&1)ans=(UL)ans*a%mod;
    return ans;
}
U mo(U x){return x>=mod?x-mod:x;}
U& mul(U&a,U b){return a=(UL)a*b%mod;}

//memmove简易封装
void mov(U*a,const U*b,S len){memmove(a,b,len*sizeof(U));}
namespace Poly
{
constexpr U g=3u,gi=332748118u;
//ml,多项式最大长度,mn,同时存在的多项式最大个数(内存池数量)
constexpr S ml=1<<21,mn=3;
U mem[(ml+16)*mn],*stk[mn],top=mn,f[(ml+16)*2],wr[ml+16],wi[ml+16],ninv[ml+16];
//统一分配回收内存,每次固定分配ml个u32,注意空间。
U*m(){return stk[--top];}void m(U*p){stk[top++]=p;}
//求最小的大于等于x的2的幂次
S up(S x){S l=1;while(l<x)l<<=1;return l;}
//使用前调用
void init()
{
    fr(i,mn)stk[i]=mem+i*(ml+16);
    U*fp;
    for(S len=1;fp=f+len,len<=ml;len<<=1)
        fr(i,len)fp[i]=(fp[i>>1]>>1)|(i&1?len>>1:0);
    for(S len=1;len<ml;len<<=1)
    {
        U Wr=pow(g,(mod-1)/(len<<1));
        U Wi=pow(gi,(mod-1)/(len<<1));
        U tr=1,ti=1;
        fr(i,len)
        {
            wr[len+i]=tr;mul(tr,Wr);
            wi[len+i]=ti;mul(ti,Wi);
        }
    }
    //预处理逆元,如果只需要乘法和求逆那就不需要写
    ninv[1]=1;
    for(S i=2;i<ml;i++)
        ninv[i]=mo(mod-mod/i*(UL)ninv[mod%i]%mod);
}
#define lst(n,a,x) poly&n(a){x return*this;}
struct poly
{
    U*mem,*a;
    S len;
    //有符号数初始化一个长度为len的多项式,系数全为0
    poly(S len):len(len){a=mem=m();cls(0,len);}
    //初始化一个长度为1的多项式
    poly():poly(1){}
    //初始化一个长度为1的常数项为x的多项式
    poly(U x):poly(){a[0]=x;}
    //用给定的长度为len的多项式初始化
    poly(U*l,S len):len(len){a=mem=m();mov(a,l,len);}
    //用另一个poly类初始化
    poly(const poly&b):poly(b.a,b.len){}
    //回收内存,注意判断不能删
    ~poly(){if(mem)m(mem);}
    //拷贝,如果是截取的多项式直接拷,否则先rsz再拷
    lst(operator=,const poly&b,if(mem)rsz(b.len);mov(a,b.a,b.len);)
    U& operator[](S idx) {return a[idx];}

    //取得多项式的一段,注意范围不要取越界(可能有未初始化的部分),或者手动先rsz,本方法为了方便直接用new生成临时对象,不会释放内存,注意不要调用过多次。
    poly& operator()(S l,S len){poly&t=*new poly;m(t.mem);t.mem=0;t.a=a+l;t.len=len;return t;}
    //cls从l下标开始清空len个数.
    poly& cls(S l,S len){memset(a+l,0,len*sizeof(U));return *this;}
    //修改多项式长度,自动清空
    lst(rsz,S nlen,if(nlen>len)cls(len,nlen-len);len=nlen;)
    //左移,不改变多项式长度
    lst(lmov,S off,mov(a,a+off,len-off);cls(len-off,off);)
    //右移,改变多项式长度
    lst(rmov,S off,rsz(len+off);mov(a+off,a,len-off);cls(0,off);)
    //反转多项式
    lst(rsv,,fr(i,len/2)swap(a[i],a[len-i-1]);)
    //NTT和NTTi,不传参以多项式长度计算,传参则改变多项式长度然后计算
    template<U*wp=wr>
    lst(NTT,S len=-1,static UL a[ml+16];
        if(~len)rsz(len);else len=this->len;
        fr(i,len)a[i]=this->a[f[len+i]];

        for(S i=1;i<len;i<<=1)
        {
            if(i==1u<<18)fr(j,len)a[j]%=mod;
            U*w=wp+i;
            for(S j=0;j<len;j+=i<<1)
                fr(k,i)
                {
                    UL x=a[j+k];
                    UL y=a[i+j+k]*w[k]%mod;
                    a[j+k]=x+y;
                    a[i+j+k]=x-y+mod;
                }
        }
        fr(i,len)this->a[i]=a[i]%mod;
    )
    lst(NTTi,S len=-1,NTT<wi>(len);*this*=pow(this->len,mod-2);)
    
    poly operator+(const poly&b){return poly(*this)+=b;}
    poly operator-(const poly&b){return poly(*this)-=b;}
    poly operator*(const poly&b){return poly(*this)*=b;}
    poly operator/(const poly&b){return poly(*this)/=b;}
    poly operator*(U x){return poly(*this)*=x;}

    //加减法长度变为两者长度最大值
    lst(operator+=,const poly&b,rsz(max(len,b.len));fr(i,b.len)a[i]=mo(a[i]+b.a[i]);)
    lst(operator-=,const poly&b,rsz(max(len,b.len));fr(i,b.len)a[i]=mo(a[i]-b.a[i]+mod);)
    //乘法变为两者长度之和再减一再取到2的幂次
    lst(operator*=,poly b,S l=len+b.len-1;NTT(up(l)).vmul(b.NTT(up(l))).NTTi().rsz(l);)
    
    //求逆再乘,长度变为被除数长度向上取到2的幂次
    lst(operator/=,poly b,S l=up(len);*this*=(b.rsz(l).inv()).rsz(l);)
    /*
    //高效版求逆再乘,有需要时使用
    lst(operator/=,poly b,
        S l=up(len);
        if(l==1)return*this*=pow(b[0],mod-2);
        poly g(b.a,l/2);
        poly q=(g.inv()*poly(a,l/2)).rsz(l/2);
        *this=q-b.NTT(l).vmul(poly(q).NTT(l)).NTTi().operator-=(*this).lmov(l/2).NTT().vmul(g.NTT(l)).NTTi().rsz(l/2).rmov(l/2);
    )
    */
    //乘数
    lst(operator*=,U x,fr(i,len)mul(a[i],x);)
    //以本多项式长度,进行点值乘法
    lst(vmul,const poly&b,fr(i,len)mul(a[i],b.a[i]);)
    //求逆,长度变为原长度向上取到2的幂次
    lst(inv,,poly f(pow(a[0],mod-2));
        rsz(up(len));
        for(S n=2;n<=len;n<<=1)(f=f*2u-f*f*poly(a,n)).rsz(n);
        *this=f;
    )
    /*
    //求逆的高效版本
    lst(inv,,poly f(rsz(up(len)));
        rsz(1);a[0]=pow(a[0],mod-2);
        for(S n=2;n<=f.len;n<<=1)
        {
            poly t(*this);t.NTT(n);
            *this-=poly(f.a,n).NTT().vmul(t).NTTi().operator-=(1u).lmov(n/2).NTT().vmul(t).NTTi().rsz(n/2).rmov(n/2);
        }
    )
    */
    //多项式除法
    lst(div,poly b,S nl=len-b.len+1;rsv().rsz(nl).operator/=(b.rsv()).rsz(nl).rsv();)
    //多项式取模
    lst(Mod,poly b,*this-=poly(*this).div(b)*b;rsz(b.len-1);)

    lst(D,,fr(i,len)mul(a[i],i);lmov(1);)
    lst(I,,rmov(1).rsz(len-1);fr(i,len)mul(a[i],ninv[i]);)
    lst(ln,,poly t(*this);D().operator/=(t).I();)
    lst(exp,,poly f(rsz(up(len)));rsz(1);a[0]=1;for(S n=2;n<=f.len;n<<=1)(*this*=poly(1u)-poly(*this).rsz(n).ln()+f(0,n)).rsz(n);)
    lst(Pow,U n,U x=a[0];operator*=(pow(x,mod-2)).ln().operator*=(n).exp().operator*=(pow(x,n));)

    //打印输出,可不写
    void cprint() const {fr(i,len)printf("%u ",a[i]);puts("");}
    lst(print,,fr(i,len)printf("%u ",a[i]);puts("");)
};};
using Poly::poly;
void poly_multiply(unsigned *a, int n, unsigned *b, int m, unsigned *c)
{
    Poly::init();
    poly f(a,n+1);
    poly g(b,m+1);
    f*=g;
    mov(c,f.a,n+m+1);
}

CompilationN/AN/ACompile OKScore: N/A

Testcase #1251.493 ms99 MB + 468 KBAcceptedScore: 100


Judge Duck Online | 评测鸭在线
Server Time: 2026-03-16 15:57:41 | Loaded in 1 ms | Server Status
个人娱乐项目,仅供学习交流使用 | 捐赠