拉格朗日插值(连续插值)
\(\text{Problem}:\)拉格朗日插值2
\(\text{Solution}:\)
前置知识:\(O(n^2)\) 拉格朗日插值。
要对 \(k\in[0,n]\) 求出:
\[f(m+k)=\sum\limits_{i=0}^{n}f(n)\prod\limits_{j\not=i}\frac{m+k-j}{i-j}
\]
发现 \(\prod\) 内分子和分母值的变化是连续的,有:
\[\begin{aligned}
f(m+k)&=\sum\limits_{i=0}^{n}f(i)\frac{(m+k)^{\underline{i}}(m+k-n)^{\overline{n-i}}}{i!(n-i)!}\\
&=\sum\limits_{i=0}^{n}f(i)(-1)^{n-i}\frac{(m+k)^{\underline{n+1}}}{(m+k-i)i!(n-i)!}\\
&=(m+k)^{\underline{n+1}}\sum\limits_{i=0}^{n}\frac{(-1)^{n-i}f(i)}{i!(n-i)!}\cdot \frac{1}{m+k-i}
\end{aligned}
\]
设 \(A_{i}=\frac{(-1)^{n-i}f(i)}{i!(n-i)!}\),\(B_{i}=\frac{1}{m-n+i}\),有:
\[C_{k}=\sum\limits_{i=0}^{k}\frac{(-1)^{n-i}f(i)}{i!(n-i)!}\cdot \frac{1}{m-n+k-i}=A_{i}B_{k-i}\\
C_{n+k}=\sum\limits_{i=0}^{n}\frac{(-1)^{n-i}f(i)}{i!(n-i)!}\cdot\frac{1}{m+k-i}=A_{i}B_{k+n-i}
\]
即求出 \(C_{n},C_{n+1},...,C_{n+n}\),最后向左平移 \(n\) 位即可求得答案。利用 \(\text{NTT}\) 可以在 \(O(n\log n)\) 的时间复杂度内解决。
\(\text{Code}:\)
#include <bits/stdc++.h>
#pragma GCC optimize(3)
//#define int long long
#define ri register
#define mk make_pair
#define fi first
#define se second
#define pb push_back
#define eb emplace_back
#define is insert
#define es erase
#define vi vector<int>
#define vpi vector<pair<int,int>>
using namespace std; const int N=550010, Mod=998244353;
inline int read()
{
int s=0, w=1; ri char ch=getchar();
while(ch<'0'||ch>'9') { if(ch=='-') w=-1; ch=getchar(); }
while(ch>='0'&&ch<='9') s=(s<<3)+(s<<1)+(ch^48), ch=getchar();
return s*w;
}
int n,m,f[N];
vector<int> A,B;
int rev[N],r[24][2],fac[N],inv[N];
inline int ksc(int x,int p) { int res=1; for(;p;p>>=1, x=1ll*x*x%Mod) if(p&1) res=1ll*res*x%Mod; return res; }
inline void Get_Rev(int T) { for(ri int i=0;i<T;i++) rev[i]=(rev[i>>1]>>1)|((i&1)?(T>>1):0); }
inline void DFT(int T,vector<int> &s,int type)
{
for(ri int i=0;i<T;i++) if(rev[i]<i) swap(s[i],s[rev[i]]);
for(ri int i=2,cnt=1;i<=T;i<<=1,cnt++)
{
int wn=r[cnt][type];
for(ri int j=0,mid=(i>>1);j<T;j+=i)
{
for(ri int k=0,w=1;k<mid;k++,w=1ll*w*wn%Mod)
{
int x=s[j+k], y=1ll*w*s[j+mid+k]%Mod;
s[j+k]=(x+y)%Mod;
s[j+mid+k]=x-y;
if(s[j+mid+k]<0) s[j+mid+k]+=Mod;
}
}
}
if(!type) for(ri int i=0,inv=ksc(T,Mod-2);i<T;i++) s[i]=1ll*s[i]*inv%Mod;
}
inline void NTT(int n,int m,vector<int> &A,vector<int> B)
{
int len=n+m;
int T=1;
while(T<=len) T<<=1;
Get_Rev(T);
A.resize(T), B.resize(T);
DFT(T,A,1), DFT(T,B,1);
for(ri int i=0;i<T;i++) A[i]=1ll*A[i]*B[i]%Mod;
DFT(T,A,0);
}
signed main()
{
r[23][1]=ksc(3,119), r[23][0]=ksc(ksc(3,Mod-2),119);
for(ri int i=22;~i;i--) r[i][0]=1ll*r[i+1][0]*r[i+1][0]%Mod, r[i][1]=1ll*r[i+1][1]*r[i+1][1]%Mod;
n=read(), m=read();
for(ri int i=0;i<=n;i++) f[i]=read();
fac[0]=1;
for(ri int i=1;i<=n;i++) fac[i]=1ll*fac[i-1]*i%Mod;
inv[n]=ksc(fac[n],Mod-2);
for(ri int i=n;i;i--) inv[i-1]=1ll*inv[i]*i%Mod;
A.resize(n+1);
for(ri int i=0;i<=n;i++)
{
A[i]=1ll*f[i]*inv[i]%Mod*inv[n-i]%Mod;
if((n-i)&1) A[i]=Mod-A[i];
}
B.resize(n+n+1);
for(ri int i=0;i<=(n<<1);i++) B[i]=ksc(m-n+i,Mod-2);
NTT(n,n+n,A,B);
int now=1;
for(ri int i=m;i>=m-n;i--) now=1ll*now*i%Mod;
for(ri int i=0;i<=n;i++)
{
printf("%d ",1ll*now*A[i+n]%Mod);
now=1ll*now*(m+i+1)%Mod*ksc(m-n+i,Mod-2)%Mod;
}
puts("");
return 0;
}
夜畔流离回,暗叹永无殿。
独隐万花翠,空寂亦难迁。
千秋孰能为,明灭常久见。
但得心未碎,踏遍九重天。