跳转至

Berlekamp-Massey算法

模版

求解递推 \(O(n*k)\)

查询 \(O(k^2\log m)\)

代码
C++
const i64 P = 998244353;
namespace BerlekampMassey {
    using i64 = long long;
    vector<int> v, a, d;
    vector<i64> res, tempc, tempd, base;
    i64 qpow(i64 a, i64 b = P - 2, i64 res = 1) {
        while (b) {
            if (b & 1) res = res * a % P;
            a = a * a % P;
            b >>= 1;
        }
        return res;
    }
    vector<int> berlekamp_massey() {
        vector<int> v, last;
        int k = -1, delta = 0;
        for (int i = 0 ; i < (int) a.size() ; i++) {
            int temp = 0;
            for (int j = 0 ; j < (int) v.size() ; j++) {
                temp = (temp + 1LL * a[i - j - 1] * v[j]) % P;
            }
            if (a[i] == temp) continue;
            if (k < 0) {
                k = i;
                delta = (a[i] - temp + P) % P;
                v = vector<int>(i + 1);
                continue;
            }
            vector<int> u = v;
            i64 val = (a[i] - temp + P) % P * qpow(delta) % P;
            if ((int) v.size() < (int) last.size() + i - k) {
                v.resize((int) last.size() + i - k);
            }
            v[i - k - 1] = (v[i - k - 1] + val) % P;
            for (int j = 0 ; j < (int) last.size() ; j++) {
                v[i - k + j] = (v[i - k + j] - val * last[j] % P + P) % P;
            }
            if ((int) u.size() - i < (int) last.size() - k) {
                last = u;
                k = i;
                delta = (a[i] - temp + P) % P;
            }
        }
        return v;
    }
    void init(const vector<int> &a) {
        Linear_Solve::a = a;
        int n = a.size() + 1;
        res.assign(n, 0);
        base.assign(n, 0);
        tempc.assign(n, 0);
        tempd.assign(n, 0);
        v = berlekamp_massey();
        while (v.size() > 1 && v.back() == 0) v.pop_back();
    }
    void mul(vector<i64> &a, vector<i64> &b, int k) {
        fill(begin(tempc), begin(tempc) + k * 2, 0);
        for (int i = 0 ; i < k ; i++) {
            if (a[i]) {
                for (int j = 0 ; j < k ; j++) {
                    tempc[i + j] = (tempc[i + j] + a[i] * b[j]) % P;
                }
            }
        }
        for (int i = (k * 2) - 1 ; i >= k ; i--) {
            if (tempc[i]) {
                for (int j = 0 ; j < (int) d.size() ; j++) {
                    tempc[i - k + d[j]] = (tempc[i - k + d[j]] - tempc[i] * tempd[d[j]] % P + P) % P;
                }
            }
        }
        for (int i = 0 ; i < k ; i++) {
            a[i] = tempc[i];
        }
    }
    i64 get(i64 n) {
        i64 ans = 0, cnt = 0;
        int k = v.size();
        vector<int> b = a;
        b.resize(k);
        for (int i = 0 ; i < k ; i++) {
            tempd[k - 1 - i] = -v[i];
        }
        tempd[k] = 1;
        d.clear();
        for (int i = 0 ; i < k ; i++) {
            if (tempd[i]) {
                d.emplace_back(i);
            }
        }
        for (int i = 0 ; i < k ; i++) {
            res[i] = base[i] = 0;
        }
        res[0] = 1;
        while (1LL << cnt <= n) cnt++;
        for (int i = cnt ; i >= 0 ; i--) {
            mul(res, res, k);
            if (n >> i & 1) {
                for (int j = k - 1 ; j >= 0 ; j--) {
                    res[j + 1] = res[j];
                }
                res[0] = 0;
                for (auto &j : d) {
                    res[j] = (res[j] - res[k] * tempd[j]) % P;
                }
            }
        }
        for (int i = 0 ; i < k ; i++) {
            ans = (ans + res[i] * b[i]) % P;
        }
        return ans;
    }
} // BerlekampMassey

查询 \(O(n\log n\log m)\)

代码
C++
{ 多项式基础模版 }
Poly BerlekampMassey(const Poly &a) {
    int n = a.size();
    Poly temp(n + 1), f(n + 1), g(n + 1);
    int k = 0, last_k = 0, last_delta, last = -1;
    for (int i = 0 ; i < n ; i++) {
        int delta = P - a[i];
        for (int j = 1 ; j <= k ; j++) {
            mod_add(delta, mod_mul_t(f[j], a[i - j]));
        }
        if (delta == 0) continue;
        if (last == -1) {
            k = i + 1;
        } else {
            int t = mod_mul_t(delta, qpow(last_delta));
            mod_add(temp[i - last], t);
            for (int j = 1 ; j <= last_k ; j++) {
                mod_sub(temp[i - last + j], mod_mul_t(t, g[j]));
            }
            int p = last_k;
            last_k = k;
            k = max(k, i - last + p);
            for (int j = 1 ; j <= last_k ; j++) {
                g[j] = f[j];
            }
            for (int j = 1 ; j <= k ; j++) {
                f[j] = temp[j];
            }
        }
        last_delta = delta;
        last = i;
    }
    return f.resize(k + 1), f;
}
i64 LinearRec(const Poly &a, const Poly &r, i64 n) {
    Poly F(begin(a), begin(a) + r.size() - 1), R = -r;
    R[0] = 1;
    F *= R;
    F.resize(R.size() - 1);
    for ( ; n ; n /= 2) {
        Poly rR = R;
        for (int i = 0 ; i < (int) rR.size() ; i += 2) {
            rR[i] = mod_sub_t(P, rR[i]);
        }
        Poly FR = F * rR, RR = R * rR;
        F.resize(0); R.resize(0);
        for (int i = (n & 1) ; i < (int) FR.size() ; i+= 2) {
            F.push_back(FR[i]);
        }
        for (int i = 0 ; i < (int) RR.size() ; i += 2) {
            R.push_back(RR[i]);
        }
    }
    return F.empty() ? 0 : mod_mul_t(F[0],qpow(R[0]));
}
例题

P5487 【模板】Berlekamp–Massey 算法 - 洛谷