多项式多点求值

给定一个\(n\)次多项式\(A(x)\)和\(m\)个值\(a_i\),求出对于\(任意i\in [0,m-1],A(a_i)\)的值

前置知识:

分治FFT

多项式除法

一般优化多项式要么倍增要么分治……

然而这题看上去不像能倍增的亚子,所以就分治吧

考虑先将要求的点分为两部分

\(x[0]=\{x_0,x_1,……,x_{\frac{m}{2}}\},x[1]=\{x_{\frac{m}{2}+1},x_{\frac{m}{2}+2},……x_{m-1}\}\)

我们记\(p[0]=\prod\limits_{i=1}^{\frac{m}{2}}(x-x_i),p[1]=\prod\limits_{i=\frac{m}{2}+1}^{m-1}(x-x_i)\)

显然\(p\)可以用类似线段树建树的方法求出

考虑对\(A(x)\)进行分治

\(A(x)=D(x)p[0](x)+A[0](x)\)

当\(x\in x[0]\)的时候,\(A(x)≡A[0](x)\, (mod\, p[0])\)

\(A[0]\)的次数是小于\(p[0]\)的

这里\(A[0]\)是可以用多项式除法求出来的

\(A[1]\)同理

当\(A\)的次数小于\(100\)的时候其实就可以暴力求解了

时间复杂度\(O(nlog^2n)\)

……其实我自己也觉得没理解透彻,有锅欢迎指出

#include<bits/stdc++.h>
using namespace std;
namespace red{
#define int long long
#define ls(p) (p<<1)
#define rs(p) (p<<1|1)
#define eps (1e-8)
    inline int read()
    {
        int x=0;char ch,f=1;
        for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
        if(ch=='-') f=0,ch=getchar();
        while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
        return f?x:-x;
    }
    const int N=64444,mod=998244353;
    int n,m,limit,len;
    vector<int> poly[266666],a;
    int pos[266666],b[N],ret[N];
    int g[21][266666];
    inline int fast(int x,int k)
    {
        int ret=1;
        while(k)
        {
            if(k&1) ret=ret*x%mod;
            x=x*x%mod;
            k>>=1;
        }
        return ret;
    }
    inline int add(int x,const int &y)//卡常
    {
        x+=y;
        return x>mod?x-mod:x;
    }
    inline int del(int x,const int &y)
    {
        x-=y;
        return x<0?x+mod:x;
    }
    inline void init(int x)//封装
    {
        limit=1,len=0;
        while(limit<(x<<1)) limit<<=1,++len;
        for(int i=0;i<limit;++i) pos[i]=(pos[i>>1]>>1)|((i&1)<<(len-1));
    }
    inline void ntt(vector<int> &a,int inv)
    {
        while(a.size()<limit) a.push_back(0);
        for(int i=0;i<limit;++i)
            if(i<pos[i]) swap(a[i],a[pos[i]]);
        for(int mid=1,t=1;mid<limit;mid<<=1,++t)
        {
            for(int r=mid<<1,j=0;j<limit;j+=r)
            {
                for(int k=0;k<mid;++k)
                {
                    int x=a[j+k],y=g[t][k]*a[j+k+mid]%mod;
                    a[j+k]=add(x,y);
                    a[j+k+mid]=del(x,y);
                }
            }
        }
        if(inv) return;
        inv=fast(limit,mod-2);reverse(a.begin()+1,a.begin()+limit);
        for(int i=0;i<limit;++i) a[i]=a[i]*inv%mod;
    }
    inline void NTT(vector<int> a,vector<int> b,vector<int> &c)//封装一下短一点
    {
        c.clear();
        ntt(a,1);ntt(b,1);
        for(int i=0;i<limit;++i) c.push_back(a[i]*b[i]%mod);
        ntt(c,0);
    }
    inline void poly_inv(int pw,vector<int> a,vector<int> &B)//多项式乘法逆
    {
        if(pw==1){B.push_back(fast(a[0],mod-2));return;}
        poly_inv((pw+1)>>1,a,B);
        init(pw);
        while(a.size()<limit) a.push_back(0);
        for(int i=pw;i<limit;++i) a[i]=0;
        ntt(a,1);ntt(B,1);
        for(int i=0;i<limit;++i) B[i]=del(2,a[i]*B[i]%mod)*B[i]%mod;
        ntt(B,0);
        for(int i=pw;i<limit;++i) B[i]=0;
    }
    inline void get_poly(int l,int r,int p)//求出p数组
    {
        if(l==r)
        {
            poly[p].push_back(b[l]?mod-b[l]:0);
            poly[p].push_back(1);
            return;
        }
        int mid=(l+r)>>1;
        get_poly(l,mid,ls(p));get_poly(mid+1,r,rs(p));
        init(r-l+1);
        NTT(poly[ls(p)],poly[rs(p)],poly[p]);
    }
    inline void poly_mod(vector<int> a,vector<int> b,vector<int> &d,int n,int m)//多项式取模(除法)
    {
        while(a.size()<=n) a.push_back(0);
        while(b.size()<=m) b.push_back(0);
        if(n<m) return (void)(d=a);
        vector<int> apos,bpos,bposinv,c,cpos;
        d.clear();
        for(int i=0;i<=n;++i) apos.push_back(a[n-i]);
        for(int i=0;i<=m;++i) bpos.push_back(b[m-i]);
        for(int i=n-m+1;i<apos.size();++i) apos[i]=0;
        for(int i=n-m+1;i<bpos.size();++i) bpos[i]=0;
        poly_inv(n-m+1,bpos,bposinv);
        init(n-m+1);
        NTT(apos,bposinv,cpos);
        for(int i=0;i<=n-m;++i) c.push_back(cpos[n-m-i]);
        init(n);
        NTT(b,c,d);
        for(int i=0;i<m;++i) d[i]=del(a[i],d[i]);
        for(int i=m;i<limit;++i) d[i]=0;
    }
    inline void solve(vector<int> a,int p,int l,int r)
    {
        if(r-l<=100)
        {
            for(int i=l;i<=r;++i)
            {
                int s=0;
                for(int j=a.size()-1;~j;--j)
                {
                    s=add(s*b[i]%mod,a[j]);
                }
                ret[i]=s;
            }
            return;
        }
        vector<int> b;
        int mid=(l+r)>>1;
        poly_mod(a,poly[ls(p)],b,r-l,mid-l+1);
        solve(b,ls(p),l,mid);
        poly_mod(a,poly[rs(p)],b,r-l,r-mid);
        solve(b,rs(p),mid+1,r);
    }
    inline void main()
    {
        n=read(),m=read();
        for(int mid=1,t=1;mid<266666;mid<<=1,++t)//预处理原根,稍微快一点
        {
            g[t][0]=1;int Wn=fast(3,(mod-1)/(mid<<1));
            for(int k=1;k<mid;++k)
            {
                g[t][k]=g[t][k-1]*Wn%mod;
            }
        }
        for(int i=0;i<=n;++i) a.push_back(read());
        for(int i=1;i<=m;++i) b[i]=read();
        get_poly(1,m,1);
        poly_mod(a,poly[1],a,n,m);
        solve(a,1,1,m);
        for(int i=1;i<=m;++i) printf("%lld\n",ret[i]);
    }
}
signed main()
{
    red::main();
    return 0;
}

原文地址:https://www.cnblogs.com/knife-rose/p/12139633.html

时间: 2024-10-29 09:07:40

多项式多点求值的相关文章

多项式多点求值和插值

本文以存板子为主= = 对于比较一般的情况,n次多项式在n个点求值和用n个点插值可以做到,并且这也是下界. 多项式多点求值 给一个多项式F和一堆值,求出. 设,. 那么对于,,对于,.递归即可. 多项式多点插值 给一堆值.,要求求出一个n-1次多项式满足. 考虑拉格朗日插值:. 我们先考虑对于每个i,如何求出.设,那么我们就是要求. 取的时候这个式子分子分母都为0,那么我们可以用洛必达法则,这个式子就等于.那么我们可以用多点求值求出每个. 设为,现在我们就是要求,显然可以分治FFT. 具体地,还

@算法 - [email&#160;protected] 多项式的多点求值与快速插值

目录 @0 - 参考资料@ @1 - 多点求值@ @理论推导@ @参考代码@ @例题与应用@ @2 - 快速插值@ @理论推导@ @(不建议参考的)代码@ @例题与应用@(暂无) @0 - 参考资料@ Cyhlnj 的博客 @1 - 多点求值@ @理论推导@ 假设已知多项式 \(A(x)\),使用 FFT 可以将 \(A(w_n^0)\),\(A(w_n^1)\),...,\(A(w_n^{n-1})\) 的值在 \(O(n\log n)\) 的时间内快速求出. 那么问题来了,假如我现在要求解任

多点求值与暴力插值

暴力插值 给你 \(n\) 个点 \((x_i,y_i)\) ,要求求出这个 \(n\) 次多项式 \(F(x)\) 我们有 \[ F(x)=\sum_{i=1}^ny_i\frac{\!\prod\limits_{1\le j\le n,i\not=j}\!\!\!\!(x-x_j)}{\!\prod\limits_{1\le j\le n,i\not=j}\!\!\!\!(x_i-x_j)} \] 感性愉悦认识一下这个拉格朗日插值 复杂度 \(O(n^2)\) 多点求值 给你 \(X=\{x

霍纳规则解决多项式的求值问题

霍纳规则用来简化朴素多项式的求值.它将一元n次多项式的求值问题转化为n个一次式. 霍纳规则是采用最少的乘法运算策略,求多项式A(x) = anxn+ an-1xn-1+...+ a1x + a0在x0处的值,该规则是A(x0)=(...((anx0+ an-1)x0+...+ a1)x0+ a0) 1 #include<stdio.h> 2 3 int horner(int a[], int x, int n); 4 5 int main() 6 { 7 int n; 8 scanf(&quo

多点求值与快速插值

多点求值 给出 $n$ 次多项式 $A(x)$ ,求出 $m$ 个 $x_i$ 对应的 $A(x_i)$ 考虑分治,设 $L(x)=\prod_{i=1}^{\frac{n}{2}}(x-x_i)$ , $R(x)=\prod_{i=\frac{n}{2}+1}^n(x-x_i)$ 对于 $i \in [1,\frac{n}{2}],F(x_i)=(F \mod L)(x_i)$ , 对于 $i \in (\frac{n}{2},n],F(x_i)=(F \mod R)(x_i)$ 就是对于左半

多项式求值

一.一维多项式求值: P(x)=3x^6+7x^5+3x^4+3x^3+8x^2+5x+23 一个通用的计算多项式的值的算法可以采用递推的方式.首先可以将上面多项式变形为如下的等价方式: P(x)=(...((an-1x+an-2)x+an-3)x+...+a1)x+a0 通过以上表达式可以看出,只要从里往外逐层按照如下的方式递推,便可以计算得到一个一维多项式的值: Rk=Rk+1*x+ak 具体示例代码如下: 1 /* 2 * @author 3 * 一维多项式求值,从最高一项开始迭代计算 4

算法研究:一维多项式求值(秦九韶算法)

问题通式: 求指定点X处的函数值. 思路: 运用秦九韶算法 将多项式表述成如下嵌套模式: 于是可以转换成程序语言思想: 可通过for循环来实现. 代码实现: C语言代码如下: #include <stdio.h> /* *函数名:double odpf(double array[], int n, double x) *参数:double array[] 系数数组 * int n 表示n项 * double x x的值 *返回值:计算结果sum *说明:One dimensional poly

基于数值分析思想对多项式求值的原理和应用进行探究

摘要:多项式是由多个单项(符号项如:5x或者常数项4)通过四则运算组合起来的式子,如P(x)=2x^4+3x^3-3x^2+5x-1 一般的求解会将特定的x代入到上式中,一个一个的计算,共需要花费10次的乘法和4次加法运算,但是我们可以通过霍纳方法对多项式进行组合计算,在需要每秒对多个不同输入的x多次计算多项式对应的值时,该方法可以极大的提高计算效率. 原理:采用霍纳方法对上式进行分解步骤如下: P(x) = -1+x(5-3x+3x^2+2x^3) = -1+x(5+x(-3+3x+2x^2)

多项式函数插值:多项式形式函数求值的Horner嵌套算法

设代数式序列 $q_1(t), q_2(t), ..., q_{n-1}(t)$ ,由它们生成的多项式形式的表达式(不一定是多项式): $$p(t)=x_1+x_2q_1(t)+...x_nq_1(t)q_2(t)..q_{n-1}(t)=\sum\limits_{i=1}^n(x_i\prod\limits_{j=1}^{i-1}q_j(t))$$ 一般来讲,按照这个形式计算函数在 $t_0$ 点的取值的复杂度为:n-1次 $q_i(t)$ 求值,n-1次浮点数乘法(生成n个不同的乘积),n-