[笔记] 拉格朗日插值法

拉格朗日插值法

问题:给你 (n+1) 个点值,求这 (n+1) 个点确定的 (n) 次多项式 (f(x))(求出给定点 (x_0) 的值 (f(x_0)) 即可)。

我们可以直接高斯消元,(mathcal{O}(n^3))

一般的拉格朗日插值法

简单来说,拉格朗日插值法可以找出一个恰好经过直角坐标系内 (n) 个给定点的函数。

我们设所求的多项式为 (f(x)) ,点的坐标为 ((x_i,y_i)) ,那么我们有:

[f(x_0)=sum_{i=1}^n y_i prod_{i eq j} frac{x_0-x_i}{x_i-x_j} ]

我们发现,代入任意一个 (x_k,kin[1,n]) ,当 (k eq i) 时,后面的 (prod_{j eq i} frac{x_0-x_i}{x_i-x_j})(0) ;而 (k==i) 时,后面的 (prod_{j eq i} frac{x_0-x_i}{x_i-x_j})(1)

时间复杂度 (mathcal{O}(n^2))

(x) 取值连续时的拉格朗日插值法

很多时候我们做题都是先发现某个函数是多少次的多项式,然后自己取一些值代入插值求 (f(x_0))

这时我们爆算的点的横坐标可以从 (1) 开始连续取,这样我们把上面的式子中的 (x_i) 换成 (i) 有:

[f(x_0)=sum_{i=1}^{n}y_i prod_{i eq j} frac{x_0−j}{i−j} ]

我们发现:

对于分子我们可以与处理前缀积和后缀积:

(pl[i]=prod_{j=1}^i (x_0-j))

(pr[i]=prod_{j=i}^n(x_0-j))

这样可以 (mathcal{O}(1)) 求分子;

对于分母,我们可以预处理阶乘,分母实际上就是:

((-1)^{n-i}(i-1)!cdot (n-i)!)

这样我们可以一个 (mathcal{O}(log n)) 求出分母。(如果你预处理逆元则可以 (mathcal{O}(1))

重心拉格朗日插值法

我们发现,一般的拉格朗日插值法每多加入一个点时就要整个重新计算,考虑如何利用已经计算的信息。

我们对于

(f(x_0)=sum_{i=1}^n y_i prod_{i eq j} frac{x_0−x_j}{x_i−x_j})

把分子提取出来,设为 (g=prod_{i=1}^n x_0-x_i) ,则此时:

(f(x_0)=g imes sum_{i=1}^n prod_{i eq j} frac{y_i}{(x_0-x_i)(x_i-x_j)})

(t_i=frac{y_i}{prod_{i eq j} x_i-x_j}) 则:

(f(x_0)=g imes sum_{i=1}^n frac{t_i}{x_0-x_i})

因此每次多加入一个点只需要重新 (mathcal{O}(n))算它的 (t_i) 就好了。

拉格朗日插值法的应用

暴力算式子:确定题目中的某个函数的次数,就可以暴力代入,然后差值去算出这个式子。

或是说题目给定的式子不确定,你需要在程序中算出来(有回校内考试出现过

拉格朗日插值法的板子(?

Luogu P4781 【模板】拉格朗日插值

#include<cstdio>
#include<iostream>
#define ll long long
#define R register int
using namespace std;
namespace Luitaryi {
template<class I> inline I g(I& x) { x=0; register I f=1;
  register char ch; while(!isdigit(ch=getchar())) f=ch=='-'?-1:f;
  do x=x*10+(ch^48); while(isdigit(ch=getchar())); return x*=f;
} const int N=2010,M=998244353;
int n,k;
int x[N],y[N];
inline ll qpow(ll a,ll b) { register ll ret=1;
  for(;b;b>>=1,(a*=a)%=M) if(b&1) (ret*=a)%=M; return ret;
}
inline int Lagrange(const int& n,const int& k,int* x,int* y) { R ret=0;
  for(R i=1;i<=n;++i) { R up=1,dn=1;
    for(R j=1;j<=n;++j) if(i!=j) 
      up=1ll*up*(k-x[j]+M)%M,dn=1ll*dn*(x[i]-x[j]+M)%M;
    ret=(ret+1ll*y[i]*up%M*qpow(dn,M-2))%M;
  } return ret;
}
inline void main() {
  g(n),g(k); for(R i=1;i<=n;++i) g(x[i]),g(y[i]);
  printf("%d
",Lagrange(n,k,x,y));
}
} signed main() {Luitaryi::main(); return 0;}

CF622F The Sum of the k-th Powers

#include<bits/stdc++.h>
#define R register int
using namespace std;
namespace Luitaryi {
inline int g() { R x=0,f=1;
  register char s; while(!isdigit(s=getchar())) f=s=='-'?-1:f;
  do x=x*10+(s^48); while(isdigit(s=getchar())); return x*f;
} const int N=1000010,M=1000000007;
int n,k,ans;
int fac[N],pl[N],pr[N];
inline int qpow(int a,int b) { R ret=1;
  for(;b;b>>=1,a=1ll*a*a%M) if(b&1) ret=1ll*ret*a%M; return ret;
}
inline void main() { n=g(),k=g();
  fac[0]=pl[0]=pr[k+3]=1;
  for(R i=1;i<=k+2;++i) fac[i]=1ll*fac[i-1]*i%M;
  for(R i=1;i<=k+2;++i) pl[i]=1ll*pl[i-1]*(n-i)%M;
  for(R i=k+2;i;--i) pr[i]=1ll*pr[i+1]*(n-i)%M;
  for(R i=1,y=0,up,dn;i<=k+2;++i) {
    y=(y+qpow(i,k))%M;
    up=1ll*pl[i-1]*pr[i+1]%M;
    dn=((k-i)&1?-1ll:1ll)*fac[i-1]*fac[k+2-i]%M;
    ans=(ans+1ll*y*up%M*qpow(dn,M-2))%M;
  } printf("%d
",(ans+M)%M);
}
} signed main() {Luitaryi::main(); return 0;}
原文地址:https://www.cnblogs.com/Jackpei/p/12051257.html