明明有最厲害最好想的插值做法,怎麼沒有人寫呢。
思路
考慮 \(n\) 個點可以確定一個 \(n-1\) 次多項式。
如何確定。
令 \(l_i(x)=\prod_{j\not =i}\frac{(x-x_j)}{(x_i-x_j)}\)。
可以發現這個多項式在 \(x=x_i\) 時值為一,在 \(x=x_j(j\not = i)\) 時值為零。
那麼就有:
\[F(x)=\sum_{i=0}^{i<n}y_il_i(x)
\]
容易發現這個多項式恰好滿足上面的條件,當然,這就是拉格朗日插值。
如何得到這個多項式?
可以先求出:
\[G(x)=\prod(x-x_i)
\]
發現:
\[l_i(x)=\frac{G(x)}{(x-x_i)k_i}
\]
其它的是一個常數所以和起來寫成 \(k_i\) 即可。
那麼就可以 \(O(n^2)\) 求解了。
思路
#include <bits/stdc++.h>
using namespace std;
#define int long long
int mod;
inline int power(int x, int y) {
int res = 1;
while (y) {
if (y & 1) res = res * x % mod;
x = x * x % mod, y /= 2;
}
return res;
}
inline vector<int> lagrange(const vector<int> &x, const vector<int> &y) {
int n = x.size();
vector<int> a(n + 1, 0), f(n, 0);
a[0] = 1;
auto add = [&](int x) {
for (int j = n; j >= 1; j--)
a[j] = (a[j - 1] - a[j] * x) % mod;
a[0] = -a[0] * x % mod;
};
for (int i = 0; i < n; i++) add(x[i]);
for (int i = 0; i < n; i++) {
if (x[i] == 0) {
for (int j = 0; j <= n; j++) a[j] = a[j + 1];
a[n] = 0;
} else {
int iv = power(x[i], mod - 2);
a[0] = -a[0] * iv % mod;
for (int j = 1; j <= n; j++) {
a[j] = a[j] - a[j - 1];
a[j] = -a[j] * iv % mod;
}
}
int s = 1;
for (int j = 0; j < n; j++)
if (i != j) s = s * (x[i] - x[j]) % mod;
s = power(s, mod - 2) * y[i] % mod;
for (int j = 0; j < n; j++)
f[j] = (f[j] + a[j] * s) % mod;
add(x[i]);
}
for (int i = 0; i < n; i++) f[i] = (f[i] % mod + mod) % mod;
return f;
}
signed main() {
ios::sync_with_stdio(0), cin.tie(0);
vector<int> a, b, f;
cin >> mod;
for (int i = 1, x; i <= mod; i++) {
cin >> x;
a.push_back(i - 1);
b.push_back(x);
}
f = lagrange(a, b);
for (int i = 0; i < mod; i++)
cout << f[i] << " \n"[i == mod - 1];
return 0;
}