\(BM\) 算法

用处

它可以用来求常系数线性递推的系数，并且可以求出最短的

求出来有什么用呢？

你可以闷声Cayley-Hamilton定理优化递推矩阵快速幂

算法简介

首先设一个数列 \(f\)，我们想要试出其中满足

\(f_n=\sum_{i=1}^{m}a_if_{n-i}(n>m)\)

的最小的 \(m\) 以及对应的系数 \(a\)

考虑增量法构造

首先因为要求 \(n>m\)，所以 \(m=n\) 且 \(a\) 都为 \(0\) 显然是满足条件的，所以初始可以就是全 \(0\)
假设有一个长度为 \(m\) 的 \(a\) 在 \(f_{1...n-1}\) 都满足条件，并且 \(f_n\) 不满足了

设 \(delta_n=\sum_{i=1}^{m}f_{n-i}a_i-f_n\)

我们只要构造出一个长度为 \(m'\) 最短的 \(a'\)

使得 \(\sum_{i=1}^{m'}f_{n-i}a'_i=-delta_n\) 然后 \(a,a'\) 按位相加就好了

怎么找到呢，实际上我们之前已经存在有一些不满足条件的情况

假设有个 \(x\)

\(delta_x=\sum_{i=1}^{m'}f_{x-i}a'_i-f_x\)

把 \(a'\) 向后移动 \(n-x\) 位，前面补 \(n-x-1\) 个 \(0\)，第 \(n-x\) 位搞个 \(-1\)

这样得到的长度为 \(m'+n-x\) 的 \(b\) 再搞个 \(\frac{-delta_i}{delta_x}\) 乘起来就好了

搞出来的 \(b\) 显然就是我们要求的，但是可能不是最短的

万物皆可持久化把之前所有求过的 \(a\) 全部记录下来

(其实记录那个最短的系数就好了)

然后又搞个 \(fail_i\) 表示第 \(i\) 个 \(a\) 挂了的位置

最后弄个变量记录一下最短的就好了

代码

可能是对的

可以去 zzq 的博客里面搞个数据测一下正确性

# include <bits/stdc++.h>

using namespace std;

typedef long long ll;

const int maxn(3005);

const int mod(1e9 + 7);

inline void Inc(int &x, int y) {

    x = x + y >= mod ? x + y - mod : x + y;

}

inline void Dec(int &x, int y) {

    x = x - y < 0 ? x - y + mod : x - y;

}

inline int Add(int x, int y) {

    return x + y >= mod ? x + y - mod : x + y;

}

inline int Sub(int x, int y) {

    return x - y < 0 ? x - y + mod : x - y;

}

inline int Pow(ll x, int y) {

    register ll ret = 1;

    for (; y; y >>= 1, x = x * x % mod)

        if (y & 1) ret = ret * x % mod;

    return ret;

}

int n, f[maxn], dt[maxn], fail[maxn], cnt, inv, mn;

vector <int> cur, now, mncoef;

int main() {

    freopen("BM-in.txt", "r", stdin);

    register int i, j, l;

    scanf("%d", &n), mn = 0;

    for (i = 1; i <= n; ++i) scanf("%d", &f[i]);

    for (i = 1; i <= n; ++i) {

        dt[i] = mod - f[i], l = now.size();

        for (j = 0; j < l; ++j) Inc(dt[i], (ll)f[i - j - 1] * now[j] % mod);

        if (!dt[i]) continue;

        fail[cnt] = i;

        if (!cnt) {

			now.clear(), now.resize(i), ++cnt;

            continue;

        }

        inv = mod - (ll)dt[i] * Pow(dt[fail[mn]], mod - 2) % mod, l = mncoef.size();

        cur.clear(), cur.resize(i - fail[mn] - 1), cur.push_back(mod - inv);

        for (j = 0; j < l; ++j) cur.push_back((ll)inv * mncoef[j] % mod);

        if (now.size() > cur.size()) cur.resize(now.size());

        for (l = now.size(), j = 0; j < l; ++j) Inc(cur[j], now[j]);

        if (now.size() - i < mncoef.size() - fail[mn]) mn = cnt, mncoef = now;

		now = cur, ++cnt;

    }

    cout << cur.size() << endl;

    return 0;

}

巴特西

Berlekamp-Massey算法

用处

算法简介

代码

最新文章

热门文章