cp-library

競技プログラミング用 C++ ライブラリ

View the Project on GitHub qandori/cp-library

:heavy_check_mark: monoid floor sum

math/monoid_floor_sum.hpp · GitHub で見る

共通ページ

#include "lib/math/monoid_floor_sum.hpp"
#include <string>

std::string op(std::string a, std::string b) { return a + b; }
std::string e() { return ""; }

// main 内
auto result = monoid_floor_sum<std::string, op, e>(3, 5, 2, 1, "X", "Y");
// result = "XXYX"
操作 内容・返り値 計算量
monoid_floor_sum<S,op,e>(n,m,a,b,x,y) 下記の floor sum に対応する列を op で集約して返す \(O(\log(n+a+b+m))\)

op は結合則を満たす演算、e() は単位元です。\(n,m,a,b\) は uint64_t、\(m > 0\) とします。

最初の 4 引数は項数、法、係数、定数項、最後の 2 引数は列の要素 x,y です。例は文字列の連結で要素の順番を確認しています。実際には行列の積など、順序に意味のある演算にも使えます。

\(\displaystyle h(i) = \left\lfloor\frac{ai+b}{m}\right\rfloor\) として、次の列を左から op で集約します。累乗は同じ要素の繰り返しです。

\[y^{h(0)},\ x,\ y^{h(1)-h(0)},\ x,\ y^{h(2)-h(1)},\ \ldots,\ x,\ y^{h(n)-h(n-1)}\]

時間は \(O(\log(n+a+b+m))\)、追加空間は \(O(1)\) です。op・e と S のコピーは \(O(1)\) とします。

検証コード

コード

#ifndef MONOID_FLOOR_SUM_HPP
#define MONOID_FLOOR_SUM_HPP

#include <cassert>
#include <concepts>
#include <cstdint>
#include <utility>

template <class S, auto op, auto e>
S monoid_floor_sum(std::uint64_t n, std::uint64_t m, std::uint64_t a, std::uint64_t b, S x, S y) {
	assert(m > 0);
	auto append = [](S value, const S &base, std::uint64_t k) {
		if (!k) return value;
		S p = base;
		for (;;) {
			if (k & 1) value = op(std::move(value), p);
			k >>= 1;
			if (!k) return value;
			p = op(p, p);
		}
	};
	S prefix = e();
	if (b >= m) {
		prefix = append(std::move(prefix), y, b / m);
		b %= m;
	}
	if (!n) return prefix;
	S suffix = e();
	for (;;) {
		if (a >= m) {
			x = append(std::move(x), y, a / m);
			a %= m;
		}
		__uint128_t z = __uint128_t(a) * n + b;
		if (z < m) return op(append(std::move(prefix), x, n), std::move(suffix));
		std::uint64_t k, r;
		if (z >> 64) k = z / m, r = z % m;
		else {
			std::uint64_t v = z;
			k = v / m, r = v % m;
		}
		suffix = op(append(y, x, r / a), std::move(suffix));
		n = k - 1;
		std::uint64_t t = m - b - 1;
		prefix = append(std::move(prefix), x, t / a + 1);
		b = t % a;
		if (!n) return op(std::move(prefix), std::move(suffix));
		std::swap(a, m);
		std::ranges::swap(x, y);
	}
}

#endif
Back to top page