Stern-Brocot Tree
(Math/sternbrocot.hpp)

• View this file on GitHub
• Last update: 2024-04-26 03:18:17+09:00
• Include: #include "Math/sternbrocot.hpp"

Depends on

• Fraction (Math/fraction.hpp)

Code

#pragma once
#include "Math/fraction.hpp"

namespace SternBrocotTree {
// R start
static vector<int> encode(Frac x) {
    return get(x).first;
}
static Frac decode(vector<int> &v) {
    Frac L(0, 1), R(1, 0);
    rep(i, 0, v.size()) {
        if (i & 1)
            R = Frac(L.a * v[i] + R.a, L.b * v[i] + R.b);
        else
            L = Frac(R.a * v[i] + L.a, R.b * v[i] + L.b);
    }
    return Frac(L.a + R.a, L.b + R.b);
}
static Frac lca(Frac x, Frac y) {
    auto px = encode(x), py = encode(y);
    vector<int> q;
    rep(i, 0, min(px.size(), py.size())) {
        q.push_back(min(px[i], py[i]));
        if (q.back() != px[i] or q.back() != py[i])
            break;
    }
    return decode(q);
}
static pair<Frac, Frac> child(Frac x) {
    auto [L, R] = subtree(x);
    Frac lb(L.a + x.a, L.b + x.b), rb(R.a + x.a, R.b + x.b);
    return {lb, rb};
}
static Frac la(Frac x, ll k) {
    auto path = encode(x);
    for (;;) {
        if (path.empty())
            return Frac(-1, 1);
        if (k <= path.back()) {
            path.back() -= k;
            break;
        } else {
            k = path.back();
            path.pop_back();
        }
    }
    return decode(path);
}
static pair<Frac, Frac> subtree(Frac x) {
    return get(x).second;
}

private:
static ll ceil(ll a, ll b) {
    return (a + b - 1) / b;
}
static pair<vector<int>, pair<Frac, Frac>> get(Frac &x) {
    Frac L(0, 1), R(1, 0), mid(1, 1);
    vector<int> path;
    for (;;) {
        if (mid == x)
            break;
        ll k = ceil(x.a * L.b - x.b * L.a, x.b * R.a - x.a * R.b) - 1;
        L = Frac(L.a + k * R.a, L.b + k * R.b);
        mid = Frac(L.a + R.a, L.b + R.b);
        path.push_back(k);
        if (mid == x)
            break;
        k = ceil(x.b * R.a - x.a * R.b, x.a * L.b - x.b * L.a) - 1;
        R = Frac(R.a + k * L.a, R.b + k * L.b);
        mid = Frac(L.a + R.a, L.b + R.b);
        path.push_back(k);
    }
    return {path, {L, R}};
}
}; // namespace SternBrocotTree

/**
 * @brief Stern-Brocot Tree
 */

#line 2 "Math/fraction.hpp"

template <typename T> struct Frac {
    T a, b;
    Frac(T _a = 0) { init(_a, 1); }
    Frac(T _a, T _b) { init(_a, _b); }
    Frac &init(T _a, T _b) {
        T g = GCD(_a, _b);
        a = _a / g, b = _b / g;
        if (b < 0)
            a = -a, b = -b;
        return *this;
    }
    Frac inv() const { return Frac(b, a); }
    Frac operator-() const { return Frac(-a, b); }
    Frac &operator+=(const Frac &x) { return init(a * x.b + x.a * b, b * x.b); }
    Frac &operator-=(const Frac &x) { return init(a * x.b - x.a * b, b * x.b); }
    Frac &operator*=(const Frac &x) { return init(a * x.a, b * x.b); }
    Frac &operator/=(const Frac &x) { return init(a * x.b, b * x.a); }
    Frac operator+(const Frac &x) const { return Frac(*this) += x; }
    Frac operator-(const Frac &x) const { return Frac(*this) -= x; }
    Frac operator*(const Frac &x) const { return Frac(*this) *= x; }
    Frac operator/(const Frac &x) const { return Frac(*this) /= x; }
    bool operator<(const Frac &x) const { return a * x.b < b * x.a; }
    bool operator>(const Frac &x) const { return x < *this; }
    bool operator<=(const Frac &x) const { return !(x < *this); }
    bool operator>=(const Frac &x) const { return !(*this < x); }
    bool operator==(const Frac &x) const { return (*this <= x and x <= *this); }
    bool operator!=(const Frac &x) const { return !(*this == x); }
    T GCD(T a, T b) {
        if (b == 0)
            return a;
        else
            return GCD(b, a % b);
    }
};
template <typename T> Frac<T> between(const Frac<T> &x, const Frac<T> &y) {
    if (x.a < x.b and y.b < y.a)
        return Frac(1);
    else if (x.b <= x.a) {
        T add = floor(x.a / x.b);
        return between(x - add, y - add) + add;
    } else
        return between(y.inv(), x.inv()).inv();
}

/**
 * @brief Fraction
 * @docs docs/fraction.md
 */
#line 3 "Math/sternbrocot.hpp"

namespace SternBrocotTree {
// R start
static vector<int> encode(Frac x) {
    return get(x).first;
}
static Frac decode(vector<int> &v) {
    Frac L(0, 1), R(1, 0);
    rep(i, 0, v.size()) {
        if (i & 1)
            R = Frac(L.a * v[i] + R.a, L.b * v[i] + R.b);
        else
            L = Frac(R.a * v[i] + L.a, R.b * v[i] + L.b);
    }
    return Frac(L.a + R.a, L.b + R.b);
}
static Frac lca(Frac x, Frac y) {
    auto px = encode(x), py = encode(y);
    vector<int> q;
    rep(i, 0, min(px.size(), py.size())) {
        q.push_back(min(px[i], py[i]));
        if (q.back() != px[i] or q.back() != py[i])
            break;
    }
    return decode(q);
}
static pair<Frac, Frac> child(Frac x) {
    auto [L, R] = subtree(x);
    Frac lb(L.a + x.a, L.b + x.b), rb(R.a + x.a, R.b + x.b);
    return {lb, rb};
}
static Frac la(Frac x, ll k) {
    auto path = encode(x);
    for (;;) {
        if (path.empty())
            return Frac(-1, 1);
        if (k <= path.back()) {
            path.back() -= k;
            break;
        } else {
            k = path.back();
            path.pop_back();
        }
    }
    return decode(path);
}
static pair<Frac, Frac> subtree(Frac x) {
    return get(x).second;
}

private:
static ll ceil(ll a, ll b) {
    return (a + b - 1) / b;
}
static pair<vector<int>, pair<Frac, Frac>> get(Frac &x) {
    Frac L(0, 1), R(1, 0), mid(1, 1);
    vector<int> path;
    for (;;) {
        if (mid == x)
            break;
        ll k = ceil(x.a * L.b - x.b * L.a, x.b * R.a - x.a * R.b) - 1;
        L = Frac(L.a + k * R.a, L.b + k * R.b);
        mid = Frac(L.a + R.a, L.b + R.b);
        path.push_back(k);
        if (mid == x)
            break;
        k = ceil(x.b * R.a - x.a * R.b, x.a * L.b - x.b * L.a) - 1;
        R = Frac(R.a + k * L.a, R.b + k * L.b);
        mid = Frac(L.a + R.a, L.b + R.b);
        path.push_back(k);
    }
    return {path, {L, R}};
}
}; // namespace SternBrocotTree

/**
 * @brief Stern-Brocot Tree
 */

Back to top page