ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

【笔记】拉格朗日插值还原系数表达式

2022-02-23 13:32:14  阅读:202  来源: 互联网

标签:拉格朗 插值 res void int inline prod ll 表达式


大家都知道拉格朗日插值的公式(已知 \(n\) 个点 \((x_i,y_i)\),求唯一确定的经过这 \(n\) 个点的\(n-1\) 次多项式):

\[f(x)=\sum_{i=1}^ny_i\left(\prod_{j\ne i}\frac{x-x_j}{x_i-x_j}\right) \]

但是洛谷模板只让求了一个点的点值,就不用把系数表达式算出来。本蒻稽不会 \(O(n\log n)\) 的多项式快速插值又想在省选骗分,只好学一下 \(O(n^2)\) 的插值方法。

上面公式分成两部分,一部分是分子 \(y_i \prod_{j\ne i} \frac 1{x_i-x_j}\),\(O(n)\) 可以计算,但分母部分 \(\prod_{j\ne i}(x-x_j)\) 需要 \(O(n^2)\)。

优化方案是先把 \(F(x)=\prod_{i=1}^n (x-x_i)\) 算出来,然后用 \(F(x)/(x-x_j)\) 来得到 \(f_j(x)\).

设 \(F(x)_i\) 为 \(F(x)\) 第 \(i\) 项的系数,\(f_j(x)_i\) 同理,有:

\[F(x)_0=-x_jf_j(x)_0 \]

\[F(x)_i=f_j(x)_{i-1}-x_jf_j(x)_{i},\quad (i>1) \]

于是得到:

\[f_j(x)_0=-F(x)_0/x_j \]

\[f_j(x)_i=(f_j(x)_{i-1}-F(x)_i)/x_j,\quad (i>1) \]

但是注意这样手动模拟的前提是 \(x_j\ne 0\)。\(x_j=0\) 的情况要特判掉(\(f_j(x)=F(x)/x\)).

算出 \(f_j(x)\) 后,把之前的分子部分乘到系数上。最后 \(f(x)=\sum_{j=1}^n f_j(x)\) 就做完了。

然后把代码写出来就好了(洛谷模板):

#include <bits/stdc++.h>

using namespace std;
typedef long long ll;
const int N = 2005, D = 998244353;
int fpm(int a, int b = D - 2, int c = D) {
	int res = 1 % c;
	while(b) {
		if(b & 1) res = (ll)res * a % c;
		b >>= 1, a = (ll)a * a % c;
	}
	return res;
}
inline int _(int a) { return a + (a >> 31 & D); }
inline void Add(int &a, ll b) { a = (a + b) % D; }
inline void Sub(int &a, ll b) { a = _((a - b) % D); }
inline void Div(int &a, int b) { a = (ll)a * fpm(b) % D; }
inline void Mul(int &a, int b) { a = (ll)a * b % D; }
int n, k, x[N], y[N], f[N];
void trans(int x[], int y[], int f[]) {
	static int c[N], F[N], t[N], ix[N];
	for(int i = 1; i <= n; i++) c[i] = 1, f[i] = 0;
	for(int i = 1; i <= n; i++) {
		ix[i] = fpm(x[i]);
		for(int j = 1; j <= n; j++)
			if(j != i) Mul(c[i], _(x[i] - x[j]));
		c[i] = (ll)fpm(_(c[i])) * y[i] % D;
	}
	for(int i = 0; i <= n; i++) F[i] = 0;
	F[0] = 1;
	for(int i = 1; i <= n; i++)
		for(int j = i; j >= 0; j--) 
			F[j] = _((F[j - 1] - (ll)x[i] * F[j]) % D);
	for(int i = 1; i <= n; i++) {
		if(!x[i]) for(int j = 1; j <= n; j++) t[j - 1] = F[j];
		else {
			t[0] = _((ll)-F[0] * ix[i] % D);
			for(int j = 1; j < n; j++)
				t[j] = _((ll)(t[j - 1] - F[j]) * ix[i] % D);
		}
		for(int j = 0; j < n; j++) Add(f[j], (ll)t[j] * c[i] % D);
	}
}

int main() {
	ios::sync_with_stdio(false); cin.tie(nullptr);
	cin >> n >> k;
	for(int i = 1; i <= n; i++) cin >> x[i] >> y[i];
	trans(x, y, f);
	int tmp = 1, ans = f[0];
	for(int i = 1; i < n; i++)
		Mul(tmp, k), Add(ans, (ll)tmp * f[i]);
	cout << ans << '\n';

	return 0;
}

标签:拉格朗,插值,res,void,int,inline,prod,ll,表达式
来源: https://www.cnblogs.com/huaruoji/p/15926926.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有