多点求值
给出 $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)$
就是对于左半部分来说, $F(x)=L(x)q(x)+r(x)$ , 把 $x_i$ 带入后只剩下了 $r(x)$ 的部分, 右半部分也同理
于是我们在 $O(nlog^2n)$ 完成了多点求值
快速插值
考虑拉格朗日插值: $F(x)=\sum_{i=1}^n\frac{\prod_{j \ne i}(x-x_j)}{\prod_{j \ne i}(x_i-x_j)}y_i$
如果暴力做的话只能做到 $O(n^2)$ , 考虑怎样优化
下面的部分是个常数, 设 $M(x)=\prod_{i=1}^n(x-x_i)$ , 考虑对于每个 $i$ , 下面的式子相当于 $\frac{M(x)}{x-x_i}$ 带入 $x_i$ 之后的值, 但是出现了 $\frac{0}{0}$ , 根据洛必达法则, 它的值相当于 $M'(x)$ 带入 $x_i$ 后的值, 于是可以分治 + $Ntt$ 后求导, 再多点求值求出每个 $x_i$ 对应的值
于是 $F(x)=\sum_{i=1}^n\prod_{j \ne i}(x-x_j)v_i$ , 这一部分就可以递归处理
设 $L(x)=\prod_{i=1}^{\frac{n}{2}}(x-x_i)$ , $R(x)=\prod_{i=\frac{n}{2}+1}^n(x-x_i)$,
于是 $F(x)=R(x)\sum_{i=1}^{\frac{n}{2}}v_i\prod_{j=1}^{\frac{n}{2}}[j \ne i](x-x_j)+L(x)\sum_{i=\frac{n}{2}+1}^{n}v_i\prod_{j=\frac{n}{2}+1}^{n}[j \ne i](x-x_j)$
于是我们在 $O(nlog^2n)$ 完成了快速插值
代码
- (代码本机跑是挂的但交洛谷是 OK 的不知道为什么, 哪位好心人帮我解答感谢您嘞)
- #include <bits/stdc++.h>
- #define _(d) while(d(isdigit(c=getchar())))
- using namespace std;
- int Rd(){
- char c;_(!);int x=c^48;_()x=(x<<3)+(x<<1)+(c^48);return x;
- }
- int vv[20];
- void Wr(int x){
- if (!x){putchar(48);return;}
- int v=0;while(x) vv[++v]=x%10,x/=10;
- for (int i=v;i;i--) putchar(vv[i]^48);
- }
- const int N=4e5+5,P=998244353;
- int A[N],B[N],t,p,re[N],C[N],D[N],E[N],G[2][19][N],c[N],*d[N],e[N],*f[N],g[N],*h[N];
- int X(int x){return x>=P?x-P:x;}
- int K(int x,int y){
- int z=1;
- for (;y;y>>=1,x=1ll*x*x%P)
- if (y&1) z=1ll*z*x%P;
- return z;
- }
- void init(){
- for (int i=0;i<19;i++){
- G[0][i][0]=G[1][i][0]=1;
- G[0][i][1]=K(3,(P-1)/(1<<(i+1)));
- G[1][i][1]=K((P+1)/3,(P-1)/(1<<(i+1)));
- for (int j=2;j<(1<<i);j++)
- G[0][i][j]=1ll*G[0][i][j-1]*G[0][i][1]%P,
- G[1][i][j]=1ll*G[1][i][j-1]*G[1][i][1]%P;
- }
- }
- void pre(int l){
- for (t=1,p=0;t<l;t<<=1,p++);
- for (int i=0;i<t;i++)
- re[i]=(re[i>>1]>>1)|((i&1)<<(p-1));
- }
- void Ntt(int *a,int o){
- for (int i=0;i<t;i++)
- if (i<re[i]) swap(a[i],a[re[i]]);
- for (int v=0,i=1;i<t;i<<=1,v++){
- for (int x,y,j=0;j<t;j+=(i<<1))
- for (int k=0;k<i;k++)
- x=a[j+k],y=1ll*G[o][v][k]*a[i+j+k]%P,
- a[j+k]=X(x+y),a[i+j+k]=X(x-y+P);
- }
- if (o)
- for (int i=0,v=K(t,P-2);i<t;i++)
- a[i]=1ll*a[i]*v%P;
- }
- void dao(int *a,int l){
- for (int i=1;i<l;i++)
- a[i-1]=1ll*i*a[i]%P;
- a[l-1]=0;
- }
- void inv(int *a,int *b,int l){
- if (l==1){
- b[0]=K(a[0],P-2);
- return;
- }
- inv(a,b,(l+1)>>1);
- for (int i=0;i<l;i++)
- A[i]=a[i],B[i]=b[i];
- pre(l<<1);Ntt(A,0);Ntt(B,0);
- for (int i=0;i<t;i++)
- A[i]=1ll*A[i]*B[i]%P*B[i]%P;
- Ntt(A,1);
- for (int i=0;i<l;i++)
- b[i]=X(X(b[i]<<1)+P-A[i]);
- for (int i=0;i<t;i++) A[i]=B[i]=0;
- }
- void dvs(int *f,int *g,int *q,int *d,int n,int m){
- reverse(f,f+n+1);reverse(g,g+m+1);
- for (int i=min(n-m,m);~i;i--) E[i]=g[i];
- inv(E,q,n-m+1);
- for (int i=0;i<=n-m;i++)
- A[i]=q[i],B[i]=f[i];
- pre((n-m+1)<<1);Ntt(A,0);Ntt(B,0);
- for (int i=0;i<t;i++)
- A[i]=1ll*A[i]*B[i]%P;
- Ntt(A,1);
- for (int i=0;i<=n-m;i++) q[i]=A[i];
- for (int i=0;i<t;i++) A[i]=B[i]=0;
- reverse(q,q+n-m+1);
- reverse(f,f+n+1);reverse(g,g+m+1);
- for (int i=0;i<=m;i++) A[i]=g[i];
- for (int i=0;i<=n-m;i++) B[i]=q[i];
- pre(n+2);Ntt(A,0);Ntt(B,0);
- for (int i=0;i<t;i++)
- A[i]=1ll*A[i]*B[i]%P;
- Ntt(A,1);
- for (int i=0;i<m;i++)
- d[i]=X(f[i]-A[i]+P);
- for (int i=0;i<t;i++) A[i]=B[i]=E[i]=0;
- }
- #define Ls k<<1
- #define Rs k<<1|1
- #define mid ((l+r)>>1)
- void build(int k,int l,int r){
- if (l==r){
- c[k]=1;d[k]=new int[2];
- d[k][0]=X(P-C[l]);
- d[k][1]=1; return;
- }
- build(Ls,l,mid);build(Rs,mid+1,r);
- c[k]=c[Ls]+c[Rs];
- d[k]=new int[c[k]+1];
- for (int i=0;i<=c[Ls];i++) A[i]=d[Ls][i];
- for (int i=0;i<=c[Rs];i++) B[i]=d[Rs][i];
- pre(c[k]+2);Ntt(A,0);Ntt(B,0);
- for (int i=0;i<t;i++) A[i]=1ll*A[i]*B[i]%P;
- Ntt(A,1);
- for (int i=0;i<=c[k];i++) d[k][i]=A[i];
- for (int i=0;i<t;i++) A[i]=B[i]=0;
- }
- void push(int k,int l){
- e[k]=l;f[k]=new int[l+1];
- }
- void qry(int k,int l,int r){
- if (l==r){D[l]=f[k][0];return;}
- push(Ls,c[Ls]-1);push(Rs,c[Rs]-1);
- dvs(f[k],d[Ls],C,f[Ls],e[k],c[Ls]);
- for (int i=0;i<=e[k]-c[Ls];i++) C[i]=0;
- dvs(f[k],d[Rs],C,f[Rs],e[k],c[Rs]);
- for (int i=0;i<=e[k]-c[Rs];i++) C[i]=0;
- qry(Ls,l,mid);qry(Rs,mid+1,r);
- }
- void multip(int *a,int n,int *b,int m){
- push(1,m-1);
- for (int i=0;i<m;i++) f[1][i]=a[i];
- qry(1,1,m);
- }
- void work(int k,int x,int y,int l,int r){
- for (int i=0;i<=g[x];i++) A[i]=h[x][i];
- for (int i=0;i<=c[y];i++) B[i]=d[y][i];
- pre(r-l+2);Ntt(A,0);Ntt(B,0);
- for (int i=0;i<t;i++) A[i]=1ll*A[i]*B[i]%P;
- Ntt(A,1);
- for (int i=0;i<=g[k];i++) h[k][i]=X(h[k][i]+A[i]);
- for (int i=0;i<t;i++) A[i]=B[i]=0;
- }
- void calcf(int k,int l,int r){
- if (l==r){
- h[k]=new int[1];
- h[k][0]=D[l];
- return;
- }
- calcf(Ls,l,mid);calcf(Rs,mid+1,r);
- g[k]=r-l;h[k]=new int[r-l+1];
- work(k,Ls,Rs,l,r);work(k,Rs,Ls,l,r);
- }
- void insv(int *a,int *b,int n){
- for (int i=1;i<=n;i++) C[i]=a[i];
- build(1,1,n);
- dao(d[1],n+1);
- for (int i=0;i<=n;i++) C[i]=0;
- multip(d[1],n-1,a,n);
- for (int i=1;i<=n;i++)
- D[i]=1ll*b[i]*K(D[i],P-2)%P;
- calcf(1,1,n);
- for (int i=0;i<n;i++)
- Wr(h[1][i]),putchar(i<n-1?'':'\n');
- }
- int n,a[N],b[N];
- int main(){
- init();n=Rd();
- for (int i=1;i<=n;i++) a[i]=Rd(),b[i]=Rd();
- insv(a,b,n);
- return 0;
- }
来源: http://www.bubuko.com/infodetail-3398809.html