当给出了n个点的坐标,求经过这 n 个点的多项式,以及某个点的点值,可以在 O(N^2) 的时间内求出,而朴素的高斯消元需要O(N^3)
给出
我们对于每一个 i,构造一个多项式
发现,把 xi 带进去时,函数值是 yi, 而将其他的x带进去时,函数值都是0
于是我们把每个 i 的 li(x)相加,就是最后的多项式了
当要求 N 的函数值时,带进去就可以 n^2 做了
ll lagrange(ll k){ ll ans = 0; for(int i=1; i<=n; i++){ ll s1 = 1, s2 = 1; for(int j=1; j<=n; j++){ if(i != j){ s1 = mul(s1, k - x[j]); s2 = mul(s2, x[i] - x[j]); } } ans = add(ans, mul(y[i], mul(s1, power(s2, Mod-2)))); } return ans; }如果x是连续的话,可以优化到O(n)
也就是求出 所有(N - xj) 的积,然后不要的那项求一个逆元
分母是 (k-i)! *(i-1)!,要注意一下符号问题
ll lagrange(ll n, ll k){ ll pre = 1, ans = 0; if(n <= k) return f[n]; for(int i=1; i<=k; i++) pre = mul(pre, n - i); for(int i=1; i<=k; i++){ ll flag = ((k-i) & 1) ? -1 : 1; ll inv1 = power(n - i, Mod-2); ll inv2 = (flag * power(mul(fac[i-1], fac[k-i]), Mod-2) + Mod) % Mod; ans = add(ans, mul(mul(inv1, inv2), mul(f[i], pre))); } return ans; }经典应用
首先讲一个经典应用:计算
老祖宗告诉我们,这个东西是个以n为自变量的k + 1次多项式
然后直接带入k+1个点后用拉格朗日插值算即可,复杂度O(k)
也就是说,只要能判断要求的是一个多项式, 我们就可以用拉格朗日来插
最后放两道道模板题 P4593 [TJOI2018]教科书般的亵渎
这道题主要就是求 , 直接插就可以了
#include<bits/stdc++.h> #define N 100 using namespace std; const int Mod = 1e9+7; typedef long long ll; ll a[N], f[N], fac[N]; int T; ll add(ll a, ll b){ return (a+b) % Mod;} ll mul(ll a, ll b){ return (a*b) % Mod;} ll dec(ll a, ll b){ return (a+Mod-b) % Mod;} ll power(ll a, ll b){ ll ans = 1; for(;b;b>>=1){ if(b&1) ans = mul(ans, a); a = mul(a, a); } return ans; } ll lagrange(ll n, ll k){ ll pre = 1, ans = 0; if(n <= k) return f[n]; for(int i=1; i<=k; i++) pre = mul(pre, n - i); for(int i=1; i<=k; i++){ ll flag = ((k-i) & 1) ? -1 : 1; ll inv1 = power(n - i, Mod-2); ll inv2 = (flag * power(mul(fac[i-1], fac[k-i]), Mod-2) + Mod) % Mod; ans = add(ans, mul(mul(inv1, inv2), mul(f[i], pre))); } return ans; } int main(){ scanf("%d", &T); fac[1] = fac[0] = 1; for(int i=2; i<=N-5; i++) fac[i] = mul(fac[i-1], i); while(T--){ ll n, m; scanf("%lld%lld", &n, &m); for(int i=1; i<=m; i++) scanf("%lld", &a[i]); sort(a+1, a+m+1); for(int i=1; i<=m+3; i++) f[i] = add(f[i-1], power(i, m+1)); ll ans = 0; for(int i=0; i<=m; i++){ ans = add(ans, lagrange(n-a[i], m+3)); for(int j=i+1; j<=m; j++) ans = dec(ans, power(a[j]-a[i], m+1)); } printf("%lld\n", (ans % Mod + Mod) % Mod); } return 0; }[tyvj1858] xlkxc
求
首先 求
可以知道 f 是k+1项的多项式, g差分一下就是 f, 所以 g 是 k+2 项的多项式, F是k+3项的多项式
于是把 F的k+4项用 g 插出来, 然后暴力插 F的 第n项
#include<bits/stdc++.h> #define N 200 using namespace std; typedef long long ll; const ll Mod = 1234567891; int T, k, a, n, d; ll f[N], g[N], fac[N]; ll add(ll a, ll b){ return (a+b) % Mod;} ll mul(ll a, ll b){ return (a*b) % Mod;} ll power(ll a, ll b){ ll ans = 1; for(;b;b>>=1){ if(b&1) ans = mul(ans, a); a = mul(a, a); } return ans; } ll lagrange(ll *a, ll k, ll n){ ll pre = 1, ans = 0; if(n <= k) return a[n]; for(int i=1; i<=k; i++) pre = mul(pre, n - i); for(int i=1; i<=k; i++){ ll inv1 = power(n - i, Mod-2); ll inv2 = power(mul(fac[i-1], fac[k-i]), Mod-2); ll flag = ((k-i) & 1) ? -1 : 1; ans = add(ans, flag * mul(mul(inv1, inv2), mul(a[i], pre))); ans = (ans % Mod + Mod) % Mod; } return ans; } int main(){ scanf("%d", &T); fac[1] = fac[0] = 1; for(int i=2; i<=N-5; i++) fac[i] = mul(fac[i-1], i); while(T--){ scanf("%d%d%d%d", &k, &a, &n, &d); for(int i=1; i<=k+4; i++) f[i] = add(f[i-1], power(i, k)); for(int i=1; i<=k+4; i++) f[i] = add(f[i], f[i-1]); for(int i=0; i<=k+4; i++) g[i] = lagrange(f, k+4, add(a, mul(i, d))); for(int i=1; i<=k+4; i++) g[i] = add(g[i], g[i-1]); printf("%lld\n", lagrange(g, k+4, n)); } return 0; }