【學習筆記】拉格朗日插值

學習多項式的第一步。

參考資料：

attack的luogu博客

oi wiki拉格朗日插值

Apocryphal的luogu博客

1.拉格朗日插值法的簡介

問題：

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

解法1：高斯消元

顯然\(\deg \geqslant n\)的多項式有無窮個，因為根據高斯消元，一定會出現自由元。

直接把這\(n\)個點列出方程組，用高斯消元求解。

求出多項式后再求出\(f(k)\)即可。

時間復雜度\(O(n^3).\)

解法2：拉格朗日插值法

給出一個關於點\((x_i,y_i)\)的基函數：

\[g(k)=\prod_{j\not=i}^{1\leqslant j\leqslant n} \dfrac{k-x_j}{x_i-x_j} \]

容易發現：\(\forall j\not=i,g(x_j)=0.\)因為累乘中總有\(k=x_j\),使得\(k-x_j=x_j-x_j=0,g(x_j)=0.\)

對於\(j=i,g(x_j)=1.\)因為累乘中每一項均為\(\dfrac{x_i-x_j}{x_i-x_j}=1,g(x_j)=1.\)

於是我們的多項式就可以表示為：

\[f(k)=\sum_{i=1}^{n}y_i\times\prod_{j\not=i}^{1\leqslant j\leqslant n} \dfrac{k-x_j}{x_i-x_j} \]

這樣\(\forall(x_i,y_i),f(x_i)=y_i\),也可以據此求出\(f(k).\)

大概因為要求逆元，所以時間復雜度為\(O(n^2+n\log n)=O(n^2).\)

\(\tt{code:}\)

#include <bits/stdc++.h>
#define ll long long
using namespace std;
int n;
const ll mod = 998244353;
ll x[2011], y[2011], k, ans;
ll ksm(ll s1, ll s2) {
	if(!s2) return 1;
	if(s2 & 1) return s1 * ksm(s1, s2 - 1) % mod;
	ll ret = ksm(s1, s2 / 2);
	return ret * ret % mod;
}
int main() {
	scanf("%d%lld", &n, &k);
	for(int i = 1; i <= n; i++) scanf("%lld%lld", &x[i], &y[i]);
	for(int i = 1; i <= n; i++) {
		ll sum = 1, mul = 1;
		for(int j = 1; j <= n; j++) {
			if(i == j) continue;
			sum *= (k - x[j]); sum %= mod; sum += mod; sum %= mod;
			mul *= (x[i] - x[j]); mul %= mod; mul += mod; mul %= mod;
		}
		sum *= ksm(mul, mod-2); sum %= mod;
		ans += (sum * y[i]) % mod; ans %= mod;
	}
	printf("%lld\n", ans);
	return 0; 
}

2.拉格朗日插值法的拓展

問題：

同上，加入\(x_i\)的值連續的限制。

解法：

假設\(x_i\in[1,n],\)則公式化為：

\[f(k)=\sum_{i=1}^{n}y_i\times\prod_{j\not=i}^{1\leqslant j\leqslant n} \dfrac{k-j}{i-j} \]

可以維護\(num=\prod_{1\leqslant j\leqslant n}k-j,\)於是分子部分化為\(\dfrac{num}{k-i},\)因為這么做要求逆元，帶個log，很不優，所以可以預處理前綴和\(pre_i\)和后綴和\(suf_i.\)

對於分母，可以化為\(1\times 2\times···\times i-1\times (-1)\times (-2)\times ···\times(i-n)\)的形式。

實質上是一個階乘，|分母|\(=fac_{i-1}\times fac_{n-i},\)然后再處理一下正負即可。

這樣再加上一些\(O(n)\)的預處理\(pre,suf,fac,inv...\)，我們只需要\(O(n)\)即可完成。

3.拉格朗日插值法的具體應用&例題

【例1】CF622F The Sum of the k-th Powers

拉格朗日插值的經典模型，即求\(\sum_{i=1}^n i^k\).

這是個\(k+1\)次的多項式，感性證明一下：

考慮\(n\)次多項式\(g(x)=a_nx^n+a_{n-1}x^{n-1}+...+a_0\).

它的一階差分\(\Delta g(x)=g(x+1)-g(x)=a_n(x+1)^n-a_nx^n+...\).

通過二項式定理，我們發現這時的\(n\)次項會被消掉，即：

每做一次差分，多項式的度便會減一。

令\(f(x)=\sum_{i=1}^ni^k,\)則\(\Delta f(x)=f(x+1)-f(x)=(x+1)^k,\)為\(k\)次多項式。

所以\(f(x)\)為\(k+1\)次多項式。

至此，我們只需確定\(k+2\)個點，做一個插值，就可以求出\(f(n)\)了。

因為\(i\)連續，所以可以用到2.2中的拉格朗日插值拓展做到\(O(k\log k).\)

\(i^k\)也可以線性求出，這里不多作說明。

\(\tt{code:}\)

#include <bits/stdc++.h>

#define ll long long

using namespace std;

const int N = 1e6;
const ll mod = 1e9 + 7;

ll n, k, fac[N+11], pre[N+11], suf[N+11], ans;

ll ksm(ll s1, ll s2) {
	if(!s2) return 1ll;
	if(s2 % 2) return s1 * ksm(s1, s2-1) % mod;
	ll ret = ksm(s1, s2/2);
	return ret * ret % mod;
}

int main() {
	scanf("%lld%lld", &n, &k);
	ll sum = 0; 
	fac[0] = 1;
	for(int i = 1; i <= k+2; i++) 
		fac[i] = fac[i-1] * (ll)i % mod;
	pre[0] = suf[k+3] = 1;
	for(int i = 1; i <= k+2; i++)
		pre[i] = pre[i-1] * (n - i) % mod;
	for(int i = k+2; i >= 1; i--)
		suf[i] = suf[i+1] * (n - i) % mod;
	for(int i = 1; i <= k+2; i++) {
		sum += ksm(i, k); sum %= mod;
		ll num = pre[i-1] * suf[i+1] % mod;
		ll fm = fac[i-1] * fac[k+2-i] % mod;
		num *= sum; num %= mod;
		num *= ksm(fm, mod-2); num %= mod;
		if((k+2-i) % 2 == 0) ans += num;
		else ans -= num; 
		ans %= mod; ans = (ans + mod) % mod;
	}
	printf("%lld\n", ans);
	return 0;
}

【例2】luogu P4593 [TJOI2018]教科書般的褻瀆

問題轉化：先求出\(k,\)即“褻瀆”的次數。忽略大於\(n\)的\(a_i\)（雖然不知道有沒有），容易發現一次“褻瀆”可以使一段血量連續的怪物全部死掉，所以\(k=m+1.\)

又因為怪物的血量互不相同，所以我們可以把求的值轉化為\(\sum_{i=1}^n i^k,\)於是就和【例1】做法相似了。

注意如果有一段連續的\([x,n]\)的空位，不要計入答案。

代碼就不放了。