Bostan-Mori Algorithm
(fps/bostan_mori.hpp)

• View this file on GitHub
• Last update: 2023-11-14 20:46:04+09:00
• Include: #include "fps/bostan_mori.hpp"

説明

$P(x)$ を高々 $d-1$ 次の多項式、 $Q(x)$ を $d$ 次の多項式とし、 $[x^0]Q(x) = 1$ とする。

$[x^N] \frac{P(x)}{Q(x)}$

を求める。$d$ 次の多項式積を $M(d)$ でできるとすると、計算量は $O(M(d) \log N)$ となる。

bostan_mori_algorithm(std::int64_t n, std::vector p, std::vector q)

$[x^N] \frac{P(x)}{Q(x)}$ を求める。

kitamasa(std::int64_t n, std::vector a, std::vector c)

$a$ を 数列の先頭 $d$ 項として、線形漸化式 $a_i = \sum_{j = 1}^d c_j a_{i-j}$ を満たす数列の第 $n$ 項を求める。 $O(M(d) \log N)$

Depends on

• テンプレート (template/template.hpp)

Verified with

• test/polynomial/Kth_term_of_Linearly_Recurrent_Sequence.test.cpp

Code

#pragma once

#include "../template/template.hpp"

namespace lib {

template <class T, std::vector<T> (*convolution)(const std::vector<T> &,
                                                 const std::vector<T> &)>
T bostan_mori_algorithm(long long n, std::vector<T> p, std::vector<T> q) {
    while (n > 0) {
        auto q_neg = q;
        for (int i = 1; i < (int)q_neg.size(); i += 2) q_neg[i] = -q_neg[i];
        p = convolution(p, q_neg);
        q = convolution(q, q_neg);
        for (int i = (n & 1LL); i < (int)p.size(); i += 2) p[i >> 1] = p[i];
        p.resize(((int)p.size() + 1 - (n & 1LL)) / 2);
        for (int i = 0; i < (int)q.size(); i += 2) q[i >> 1] = q[i];
        q.resize(((int)q.size() + 1) / 2);
        n >>= 1;
    }
    return p[0] / q[0];
}

template <class T, std::vector<T> (*convolution)(const std::vector<T> &,
                                                 const std::vector<T> &)>
T kitamasa(std::int64_t n, std::vector<T> a, std::vector<T> c) {
    if (n < (int)a.size()) return a[n];
    const int d = c.size();
    for (auto &val : c) val = -val;
    c.insert(c.begin(), 1);
    auto p = convolution(a, c);
    p.resize(d);
    return bostan_mori_algorithm<T, convolution>(n, p, c);
}

}  // namespace lib

#line 2 "fps/bostan_mori.hpp"

#line 2 "template/template.hpp"

#include <bits/stdc++.h>

#define rep(i, s, n) for (int i = (int)(s); i < (int)(n); i++)
#define rrep(i, s, n) for (int i = (int)(n)-1; i >= (int)(s); i--)
#define all(v) v.begin(), v.end()

using ll = long long;
using ld = long double;
using ull = unsigned long long;

template <typename T> bool chmin(T &a, const T &b) {
    if (a <= b) return false;
    a = b;
    return true;
}
template <typename T> bool chmax(T &a, const T &b) {
    if (a >= b) return false;
    a = b;
    return true;
}

namespace lib {

using namespace std;

}  // namespace lib

// using namespace lib;
#line 4 "fps/bostan_mori.hpp"

namespace lib {

template <class T, std::vector<T> (*convolution)(const std::vector<T> &,
                                                 const std::vector<T> &)>
T bostan_mori_algorithm(long long n, std::vector<T> p, std::vector<T> q) {
    while (n > 0) {
        auto q_neg = q;
        for (int i = 1; i < (int)q_neg.size(); i += 2) q_neg[i] = -q_neg[i];
        p = convolution(p, q_neg);
        q = convolution(q, q_neg);
        for (int i = (n & 1LL); i < (int)p.size(); i += 2) p[i >> 1] = p[i];
        p.resize(((int)p.size() + 1 - (n & 1LL)) / 2);
        for (int i = 0; i < (int)q.size(); i += 2) q[i >> 1] = q[i];
        q.resize(((int)q.size() + 1) / 2);
        n >>= 1;
    }
    return p[0] / q[0];
}

template <class T, std::vector<T> (*convolution)(const std::vector<T> &,
                                                 const std::vector<T> &)>
T kitamasa(std::int64_t n, std::vector<T> a, std::vector<T> c) {
    if (n < (int)a.size()) return a[n];
    const int d = c.size();
    for (auto &val : c) val = -val;
    c.insert(c.begin(), 1);
    auto p = convolution(a, c);
    p.resize(d);
    return bostan_mori_algorithm<T, convolution>(n, p, c);
}

}  // namespace lib

Back to top page