Berlekamp-Massey算法学习笔记
Berlekamp-Massey算法
很久之前就听说过这个算法，当时六校联考的时候Day1T1是一道很有意思的递推，神仙zzx不会做于是就拿BM算法艹出了递推式Orzzzzzzzzzzx
我就不详细说了，只记一下自己感觉比较难理解的地方
设\(r(m)\)表示序列的递推式且长度为\(m\)
\(f(r, i)\)表示\(\sum_{j = 1}^m r_j * a[i - j]\)
\(\delta(r, i)\)表示\(a[i] - f(r, i)\)
\(fail_i\)表示第\(i\)个递推式出错的位置
对于某一个位置\(i\)，如果我们求出的\(\delta(r, i) \not = 0\)，这时候我们需要构造一个递推式\(r'(m')\)，满足\(\forall j \in [m' + 1, i - 1] f(r', j) = 0\)且\(f(r, i) = \delta(r, i)\)
这样我们令\(r = r + r'\)就得到新位置的递推式了
\(r'\)可以这么构造
设\(mul = \frac{\delta(r, i)}{\delta(r, fail_{cnt - 1})}\)
那么\(r' = \{0, 0, 0 \dots, 0, mul, -mul * R_{cnt - 1} \}\)
\(0\)的个数为\(i - fail_{cnt - 1} - 1\)
至于为什么这么构造是对的，我思考了挺长时间，简单的证明一下
首先对于\(\forall j \in [m' + 1, i - 1]\), \(\delta(r', j) = 0\)
仔细想了想，，发现自己并不会证。。如果哪位大佬会的话可以教教本蒟蒻
感性理解就是因为\(r\)在\([1, M]\)处满足任意位置为\(0\)，然后右移一下还满足？。。
至于为什么\(f(r', i) = \delta(r, i)\)
可以这么考虑，前\(i - fail_{cnt - 1} - 1\)个位置产生的贡献为\(0\)
\(mul\)产生的贡献为\(mul * a_{fail_{cnt - 1}}\)
\(-mul * R_{cnt - 1}\)产生的贡献为\(-mul * (a[fail_{cnt - 1}] - \delta(r, fail_{cnt - 1]})\)
合并同类项后可以得到\(mul * \delta(r, fail_{cnt - 1}) = \delta(r, i)\)
代码如下
#include<bits/stdc++.h>
using namespace std;
const int MAXN = 2005;
const double eps = 1e-8;
int cnt, fail[MAXN];
double val[MAXN], delta[MAXN];
vector <double> ans[MAXN];
int main() {
int N; scanf("%d", &N);
for (int i = 1; i <= N; i++) scanf("%lf", &val[i]);
for (int i = 1; i <= N; i++) {
double tmp = val[i];
for (int j = 0; j < ans[cnt].size(); j++)
tmp -= ans[cnt][j] * val[i - j - 1];
delta[i] = tmp;
if (fabs(tmp) <= eps) continue;
fail[cnt] = i;
if (cnt == 0) {
ans[++cnt].resize(i);
continue;
}
double mul = delta[i] / delta[fail[cnt - 1]];
cnt++; ans[cnt].resize(i - fail[cnt - 2] - 1);
ans[cnt].push_back(mul);
for (int j = 0; j < ans[cnt - 2].size(); j++)
ans[cnt].push_back(ans[cnt - 2][j] * -mul);
if (ans[cnt].size() < ans[cnt - 1].size()) ans[cnt].resize(ans[cnt - 1].size());
for (int j = 0; j < ans[cnt - 1].size(); j++)
ans[cnt][j] += ans[cnt - 1][j];
}
for (int i = 0; i < ans[cnt].size(); i++)
cout << ans[cnt][i] << ' ';
return 0;
}