ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

Luogu4512 【模板】多项式除法

2021-05-05 14:33:49  阅读:132  来源: 互联网

标签:Luogu4512 include frac int void return sizeof 除法 模板


Luogu4512 【模板】多项式除法

\(NTT\)

\[F(x)=Q(x) \times G(x) + R(x) \]

对于一个\(n\)次多项式\(f(x)=\sum_{i=0}^n a_i x^i\),定义\(f^r(x)=\sum_{i=0}^n a_{n-i} x^i\)。

\[F(\frac{1}{x})=Q(\frac{1}{x}) \times G(\frac{1}{x}) + R(\frac{1}{x})\\ x^n F(\frac{1}{x})=x^n (Q(\frac{1}{x}) \times G(\frac{1}{x}) + R(\frac{1}{x}))\\ F^r(x)=Q^r(x)G^r(x)+x^{n-m+1}R^r(x)\\ F^r(x)=Q^r(x)G^r(x) \mod(x^{n-m+1})\\ Q^r(x)=\frac{F^r(x)}{G^r(x)} \mod(x^{n-m+1}) \\ R(x)=F(x)-G(x)Q(x) \]

\(Code:\)

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#define N 400005
#define ll long long
using namespace std;
const int p=998244353;
int n,m,f[N],g[N],F[N],G[N],Q[N],R[N],c[N],d[N];
int s,l,E[2][35],rev[N];
void Add(int &x,int y)
{
    x=(x+y>=p)?(x+y-p):(x+y);
}
void Del(int &x,int y)
{
    x=(x>=y)?(x-y):(x-y+p);
}
void Mul(int &x,int y)
{
    x=(ll)x*y%p;
}
int add(int x,int y)
{
    return (x+y>=p)?(x+y-p):(x+y);
}
int del(int x,int y)
{
    return (x>=y)?(x-y):(x-y+p);
}
int mul(int x,int y)
{
    return (ll)x*y%p;
}
int ksm(int x,int y)
{
    int ans(1);
    while (y)
    {
        if (y & 1)
            Mul(ans,x);
        Mul(x,x);
        y >>=1;
    }
    return ans;
}
void Pre()
{
    E[0][23]=ksm(3,(p-1)/(1 << 23));
    E[1][23]=ksm(E[0][23],p-2);
    for (int i=22;i;--i)
    {
        E[0][i]=mul(E[0][i+1],E[0][i+1]);
        E[1][i]=mul(E[1][i+1],E[1][i+1]);
    }
}
void solve(int n)
{
    s=1,l=0;
    while (s<n)
        s <<=1,++l;
    for (int i=0;i<s;++i)
        rev[i]=(rev[i >> 1] >> 1) | ((i & 1) << l-1);
}
void Cut(int *a,int n,int s)
{
    if (n>s)
        return;
    memset(a+n,0,(s-n)*sizeof(int));
}
void NTT(int *a,int t)
{
    for (int i=0;i<s;++i)
        if (i<rev[i])
            swap(a[i],a[rev[i]]);
    for (int mid=1,o=1;mid<s;mid <<=1,++o)
        for (int j=0;j<s;j+=(mid << 1))
        {
            int g(1);
            for (int k=0;k<mid;++k,Mul(g,E[t][o]))
            {
                int x(a[j+k]),y(mul(g,a[j+k+mid]));
                a[j+k]=add(x,y),a[j+k+mid]=del(x,y);
            }
        }
}
void GetInv(int *f,int *g,int R)
{
    if (R==2)
    {
        g[0]=ksm(f[0],p-2);
        return;
    }
    GetInv(f,g,R >> 1);
    memcpy(c,g,(R >> 2)*sizeof(int)),memcpy(d,f,(R >> 1)*sizeof(int));
    solve(R),NTT(c,0),NTT(d,0);
    for (int i=0;i<s;++i)
        c[i]=del(add(c[i],c[i]),mul(d[i],mul(c[i],c[i])));
    NTT(c,1);
    int iv(ksm(s,p-2));
    for (int i=0;i<s;++i)
        Mul(c[i],iv);
    memcpy(g,c,(R >> 1)*sizeof(int)),memset(c,0,s*sizeof(int)),memset(d,0,s*sizeof(int));
}
int main()
{
    Pre();
    scanf("%d%d",&n,&m);
    for (int i=0;i<=n;++i)
        scanf("%d",&f[i]);
    for (int i=0;i<=m;++i)
        scanf("%d",&g[i]);
    for (int i=0;i<=n;++i)
        F[i]=f[n-i];
    for (int i=0;i<=m;++i)
        G[i]=g[m-i];
    Cut(F,n-m+1,n+1);
    Cut(G,n-m+1,m+1);
    s=1;
    while (s<n-m+1)
        s <<=1;
    s <<=1;
    GetInv(G,Q,s);
    Cut(Q,n-m+1,s);
    solve((n-m+1) << 1);
    NTT(Q,0),NTT(F,0);
    for (int i=0;i<s;++i)
        Mul(Q[i],F[i]);
    NTT(Q,1);
    int iv(ksm(s,p-2));
    for (int i=0;i<s;++i)
        Mul(Q[i],iv);
    Cut(Q,n-m+1,s);
    reverse(Q,Q+n-m+1);
    for (int i=0;i<n-m+1;++i)
        printf("%d ",Q[i]);
    putchar('\n');
    solve(n-m+1+m+1);
    NTT(Q,0),NTT(g,0);
    for (int i=0;i<s;++i)
        Mul(Q[i],g[i]);
    NTT(Q,1);
    iv=ksm(s,p-2);
    for (int i=0;i<s;++i)
        Mul(Q[i],iv);
    for (int i=0;i<n;++i)
        Del(f[i],Q[i]);
    for (int i=0;i<m;++i)
        printf("%d ",f[i]);
    putchar('\n');
    return 0;
}

标签:Luogu4512,include,frac,int,void,return,sizeof,除法,模板
来源: https://www.cnblogs.com/GK0328/p/14731637.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有