rec_add_point_get.hpp¶
#include "noya/rec_add_point_get.hpp"
#ifndef NOYA_REC_ADD_POINT_GET_HPP
#define NOYA_REC_ADD_POINT_GET_HPP 1
#include "noya/persistent_segtree.hpp"
#include "noya/rectangle_sum.hpp"
namespace noya {
/// @brief Offline rectangle-add point-get: add weight w to all points in [l, r) x [d, u).
template <class T, class C = noya::fenwick<T>>
std::vector<T>
static_rectangle_add_point_get(std::vector<std::array<int, 5>> areas,
std::vector<std::array<int, 2>> points) {
std::vector<std::array<T, 3>> weighted_points;
for (auto &[l, r, d, u, w] : areas) {
weighted_points.push_back({l, d, w});
weighted_points.push_back({l, u, -w});
weighted_points.push_back({r, d, -w});
weighted_points.push_back({r, u, w});
}
std::vector<std::array<T, 4>> queries;
for (auto &[x, y] : points) {
queries.push_back({0, x + 1, 0, y + 1});
}
return rectangle_sum<T, C>(weighted_points, queries);
}
/// @brief Dynamic rectangle-add point-get with semi-offline rebuilds.
template <class T> struct dynamic_rectangle_add_point_get {
persistent_segtree<T> seg;
std::vector<std::array<T, 3>> weighted_points;
std::vector<std::array<T, 5>> buckets;
dynamic_rectangle_add_point_get() {}
dynamic_rectangle_add_point_get(std::vector<std::array<T, 5>> &areas) {
buckets = areas;
build();
};
const int B = 6000;
/// @brief Add weight w to all points in rectangle [l, r) x [d, u).
void add_rectangle(T l, T r, T d, T u, T w) {
buckets.push_back({l, r, d, u, w});
if (int(buckets.size()) >= B) {
build();
}
}
void build() {
if (!buckets.empty()) {
for (auto &[l, r, d, u, w] : buckets) {
weighted_points.push_back({l, d, w});
weighted_points.push_back({l, u, -w});
weighted_points.push_back({r, d, -w});
weighted_points.push_back({r, u, w});
}
buckets.clear();
}
seg.build(weighted_points);
}
/// @brief Query the accumulated weight at point (x, y).
T query(int x, int y) {
T ans = 0;
ans = seg.prod(x + 1, y + 1);
for (auto &[l, r, d, u, w] : buckets) {
if (l <= x && x < r && d <= y && y < u) {
ans += w;
}
}
return ans;
}
};
} // namespace noya
#endif // NOYA_REC_ADD_POINT_GET_HPP
#include <algorithm>
#include <array>
#include <cassert>
#include <cmath>
#include <numeric>
#include <type_traits>
#include <vector>
namespace noya {
/// @brief Persistent segment tree supporting 2D rectangle sum and k-th queries.
template <class T> struct persistent_segtree {
std::vector<int> rt;
std::vector<int> ls, rs;
int N, M;
std::vector<T> Xs, Ys;
std::vector<T> sum;
int cnt, C;
persistent_segtree() {}
persistent_segtree(std::vector<std::array<T, 3>> &points) {
if (!points.empty())
build(points);
}
void update(int pre, int &p, int l, int r, int x, T v) {
p = cnt++;
ls[p] = ls[pre];
rs[p] = rs[pre];
sum[p] = sum[pre] + v;
if (l + 1 == r)
return;
int m = (l + r) / 2;
if (x < m)
update(ls[pre], ls[p], l, m, x, v);
else
update(rs[pre], rs[p], m, r, x, v);
}
int query_kth(int L, int R, int l, int r, T k) {
if (l + 1 == r) {
return l;
}
int m = (l + r) / 2;
T val = sum[ls[R]] - sum[ls[L]];
if (val >= k) {
return query_kth(ls[L], ls[R], l, m, k);
} else {
return query_kth(rs[L], rs[R], m, r, k - val);
}
}
T query_prefix(int p, int l, int r, int x) {
if (!p)
return 0;
if (l + 1 == r)
return sum[p];
int m = (l + r) / 2;
if (x < m)
return query_prefix(ls[p], l, m, x);
else
return sum[ls[p]] + query_prefix(rs[p], m, r, x);
}
T query_single(int p, int l, int r, int x) {
if (!p)
return 0;
if (l + 1 == r)
return sum[p];
int m = (l + r) / 2;
if (x < m)
return query_single(ls[p], l, m, x);
else
return query_single(rs[p], m, r, x);
}
T query(int L, int R, int ql, int qr, int l, int r) {
if (ql <= l && r <= qr)
return sum[R] - sum[L];
int m = (l + r) / 2;
T res = 0;
if (ql < m)
res += query(ls[L], ls[R], ql, qr, l, m);
if (qr > m)
res += query(rs[L], rs[R], ql, qr, m, r);
return res;
}
/// @brief Query the prefix sum up to coordinate (r, u) in the original space.
T prod(T r, T u) {
r = std::lower_bound(Xs.begin(), Xs.end(), r) - Xs.begin();
u = std::lower_bound(Ys.begin(), Ys.end(), u) - Ys.begin();
assert(0 <= r && r <= N);
assert(0 <= u && u <= M);
if (u > 0)
return query_prefix(rt[r], 0, M, u - 1);
else
return 0;
}
/// @brief Query the sum of weights in rectangle [l, r) x [d, u).
T prod(T l, T r, T d, T u) {
l = std::lower_bound(Xs.begin(), Xs.end(), l) - Xs.begin();
r = std::lower_bound(Xs.begin(), Xs.end(), r) - Xs.begin();
d = std::lower_bound(Ys.begin(), Ys.end(), d) - Ys.begin();
u = std::lower_bound(Ys.begin(), Ys.end(), u) - Ys.begin();
assert(0 <= l && l <= r && r <= N);
assert(0 <= d && d <= u && u <= M);
if (l == r)
return 0;
else if (d == u)
return 0;
else
return query(rt[l], rt[r], d, u, 0, M);
}
/// @brief Return the k-th smallest Y-coordinate among points with X in [l, r).
T kth(T l, T r, T k) {
l = std::lower_bound(Xs.begin(), Xs.end(), l) - Xs.begin();
r = std::lower_bound(Xs.begin(), Xs.end(), r) - Xs.begin();
assert(0 <= l && l <= r && r <= N);
if (sum[rt[r]] - sum[rt[l]] < k) {
return std::numeric_limits<T>::max();
}
int res = query_kth(rt[l], rt[r], 0, M, k);
return Ys[res];
}
/// @brief Build the persistent segment tree from weighted 2D points {x, y, w}.
void build(std::vector<std::array<T, 3>> points) {
C = int(points.size()) * 30;
cnt = 1;
ls.assign(C, 0);
rs.assign(C, 0);
sum.assign(C, 0);
Xs.clear();
Ys.clear();
for (auto &[x, y, w] : points) {
Xs.push_back(x);
Ys.push_back(y);
}
std::sort(Xs.begin(), Xs.end());
std::sort(Ys.begin(), Ys.end());
Xs.erase(std::unique(Xs.begin(), Xs.end()), Xs.end());
Ys.erase(std::unique(Ys.begin(), Ys.end()), Ys.end());
N = int(Xs.size());
std::vector<std::vector<std::array<T, 2>>> add(N);
for (auto &[x, y, w] : points) {
x = std::lower_bound(Xs.begin(), Xs.end(), x) - Xs.begin();
y = std::lower_bound(Ys.begin(), Ys.end(), y) - Ys.begin();
add[x].push_back({y, w});
}
M = int(Ys.size());
rt.assign(N + 1, 0);
rt[0] = cnt++;
for (int i = 0; i < N; i++) {
rt[i + 1] = rt[i];
for (auto &[y, w] : add[i]) {
int new_rt = 0;
update(rt[i + 1], new_rt, 0, M, y, w);
rt[i + 1] = new_rt;
}
}
std::vector<std::vector<std::array<T, 2>>>().swap(add);
}
/// @brief @return vector of (ls, rs, len, value) for all nodes.
std::vector<std::array<int, 4>> all_nodes() {
std::vector<std::array<int, 4>> nodes(cnt, std::array<int, 4>{});
std::vector<bool> vis(cnt);
vis[0] = true;
for (int i = 0; i <= N; i++) {
auto dfs = [&](auto &dfs, int p, int l, int r) -> void {
if (vis[p])
return;
vis[p] = true;
nodes[p] = {ls[p], rs[p], r - l, sum[p]};
int m = (l + r) / 2;
dfs(dfs, ls[p], l, m);
dfs(dfs, rs[p], m, r);
};
dfs(dfs, rt[i], 0, M);
}
return nodes;
}
};
} // namespace noya
namespace atcoder {
namespace internal {
#ifndef _MSC_VER
template <class T>
using is_signed_int128 =
typename std::conditional<std::is_same<T, __int128_t>::value ||
std::is_same<T, __int128>::value,
std::true_type,
std::false_type>::type;
template <class T>
using is_unsigned_int128 =
typename std::conditional<std::is_same<T, __uint128_t>::value ||
std::is_same<T, unsigned __int128>::value,
std::true_type,
std::false_type>::type;
template <class T>
using make_unsigned_int128 =
typename std::conditional<std::is_same<T, __int128_t>::value,
__uint128_t,
unsigned __int128>;
template <class T>
using is_integral = typename std::conditional<std::is_integral<T>::value ||
is_signed_int128<T>::value ||
is_unsigned_int128<T>::value,
std::true_type,
std::false_type>::type;
template <class T>
using is_signed_int = typename std::conditional<(is_integral<T>::value &&
std::is_signed<T>::value) ||
is_signed_int128<T>::value,
std::true_type,
std::false_type>::type;
template <class T>
using is_unsigned_int =
typename std::conditional<(is_integral<T>::value &&
std::is_unsigned<T>::value) ||
is_unsigned_int128<T>::value,
std::true_type,
std::false_type>::type;
template <class T>
using to_unsigned = typename std::conditional<
is_signed_int128<T>::value,
make_unsigned_int128<T>,
typename std::conditional<std::is_signed<T>::value,
std::make_unsigned<T>,
std::common_type<T>>::type>::type;
#else
template <class T> using is_integral = typename std::is_integral<T>;
template <class T>
using is_signed_int =
typename std::conditional<is_integral<T>::value && std::is_signed<T>::value,
std::true_type,
std::false_type>::type;
template <class T>
using is_unsigned_int =
typename std::conditional<is_integral<T>::value &&
std::is_unsigned<T>::value,
std::true_type,
std::false_type>::type;
template <class T>
using to_unsigned = typename std::conditional<is_signed_int<T>::value,
std::make_unsigned<T>,
std::common_type<T>>::type;
#endif
template <class T>
using is_signed_int_t = std::enable_if_t<is_signed_int<T>::value>;
template <class T>
using is_unsigned_int_t = std::enable_if_t<is_unsigned_int<T>::value>;
template <class T> using to_unsigned_t = typename to_unsigned<T>::type;
} // namespace internal
} // namespace atcoder
namespace atcoder {
// Reference: https://en.wikipedia.org/wiki/Fenwick_tree
template <class T> struct fenwick_tree {
using U = internal::to_unsigned_t<T>;
public:
fenwick_tree() : _n(0) {}
explicit fenwick_tree(int n) : _n(n), data(n) {}
void add(int p, T x) {
assert(0 <= p && p < _n);
p++;
while (p <= _n) {
data[p - 1] += U(x);
p += p & -p;
}
}
T sum(int l, int r) {
assert(0 <= l && l <= r && r <= _n);
return sum(r) - sum(l);
}
private:
int _n;
std::vector<U> data;
U sum(int r) {
U s = 0;
while (r > 0) {
s += data[r - 1];
r -= r & -r;
}
return s;
}
};
} // namespace atcoder
namespace noya {
template <class T> struct block {
int V, sqrtV;
block() {}
block(const int &_V) {
if (_V > 0) {
build(_V);
}
}
std::vector<T> point, blo;
void build(const int &_V) {
V = _V;
sqrtV = sqrt(V);
point.assign(V, 0);
blo.assign(V / sqrtV + 1, 0);
}
void add(int x, T v) {
assert(0 <= x && x < V);
int bel = x / sqrtV;
blo[bel] += v;
point[x] += v;
}
T query(int x) const {
assert(0 <= x && x <= V);
T res = 0;
int bel = x / sqrtV;
for (int i = 0; i < bel; i++)
res += blo[i];
int start = bel * sqrtV;
int end = x;
for (int i = start; i < end; i++)
res += point[i];
return res;
}
/// @brief Sum of [l, r).
T prod(int l, int r) const {
assert(0 <= l && l <= r && r <= V);
return query(r) - query(l);
}
};
template <class T> struct fenwick : atcoder::fenwick_tree<T> {
using atcoder::fenwick_tree<T>::fenwick_tree;
using atcoder::fenwick_tree<T>::add;
T query(int x) { return this->sum(0, x); }
T prod(int l, int r) { return this->sum(l, r); }
};
} // namespace noya
namespace noya {
/// @brief Offline rectangle sum. Points (x, y, weight), queries [l, r) x [d, u).
template <class T, class C = noya::fenwick<T>>
std::vector<T> rectangle_sum(std::vector<std::array<int, 3>> points,
std::vector<std::array<int, 4>> queries) {
std::vector<int> Xs;
std::vector<int> Ys;
for (auto &[x, y, z] : points) {
Xs.push_back(x);
Ys.push_back(y);
}
for (auto &[l, r, d, u] : queries) {
Xs.push_back(l);
Xs.push_back(r);
Ys.push_back(d);
Ys.push_back(u);
}
std::sort(Xs.begin(), Xs.end());
Xs.erase(std::unique(Xs.begin(), Xs.end()), Xs.end());
std::sort(Ys.begin(), Ys.end());
Ys.erase(std::unique(Ys.begin(), Ys.end()), Ys.end());
for (auto &[x, y, z] : points) {
x = std::lower_bound(Xs.begin(), Xs.end(), x) - Xs.begin();
y = std::lower_bound(Ys.begin(), Ys.end(), y) - Ys.begin();
}
for (auto &[l, r, d, u] : queries) {
l = std::lower_bound(Xs.begin(), Xs.end(), l) - Xs.begin();
r = std::lower_bound(Xs.begin(), Xs.end(), r) - Xs.begin();
d = std::lower_bound(Ys.begin(), Ys.end(), d) - Ys.begin();
u = std::lower_bound(Ys.begin(), Ys.end(), u) - Ys.begin();
}
int X = int(Xs.size());
int Y = int(Ys.size());
C F(Y);
std::vector<std::vector<std::array<int, 3>>> Q(X + 1);
for (int i = 0; i < int(queries.size()); i++) {
auto &[l, r, d, u] = queries[i];
Q[r].push_back({i, d, u});
Q[l].push_back({~i, d, u});
}
std::vector<std::vector<std::array<int, 2>>> A(X + 1);
for (auto &[x, y, z] : points) {
A[x].push_back({y, z});
}
std::vector<T> ans(queries.size());
for (int x = 0; x < X; x++) {
for (auto &[y, z] : A[x]) {
F.add(y, z);
}
for (auto &[i, d, u] : Q[x + 1]) {
if (i >= 0)
ans[i] += F.prod(d, u);
else
ans[~i] -= F.prod(d, u);
}
}
return ans;
}
} // namespace noya
namespace noya {
/// @brief Offline rectangle-add point-get: add weight w to all points in [l, r) x [d, u).
template <class T, class C = noya::fenwick<T>>
std::vector<T>
static_rectangle_add_point_get(std::vector<std::array<int, 5>> areas,
std::vector<std::array<int, 2>> points) {
std::vector<std::array<T, 3>> weighted_points;
for (auto &[l, r, d, u, w] : areas) {
weighted_points.push_back({l, d, w});
weighted_points.push_back({l, u, -w});
weighted_points.push_back({r, d, -w});
weighted_points.push_back({r, u, w});
}
std::vector<std::array<T, 4>> queries;
for (auto &[x, y] : points) {
queries.push_back({0, x + 1, 0, y + 1});
}
return rectangle_sum<T, C>(weighted_points, queries);
}
/// @brief Dynamic rectangle-add point-get with semi-offline rebuilds.
template <class T> struct dynamic_rectangle_add_point_get {
persistent_segtree<T> seg;
std::vector<std::array<T, 3>> weighted_points;
std::vector<std::array<T, 5>> buckets;
dynamic_rectangle_add_point_get() {}
dynamic_rectangle_add_point_get(std::vector<std::array<T, 5>> &areas) {
buckets = areas;
build();
};
const int B = 6000;
/// @brief Add weight w to all points in rectangle [l, r) x [d, u).
void add_rectangle(T l, T r, T d, T u, T w) {
buckets.push_back({l, r, d, u, w});
if (int(buckets.size()) >= B) {
build();
}
}
void build() {
if (!buckets.empty()) {
for (auto &[l, r, d, u, w] : buckets) {
weighted_points.push_back({l, d, w});
weighted_points.push_back({l, u, -w});
weighted_points.push_back({r, d, -w});
weighted_points.push_back({r, u, w});
}
buckets.clear();
}
seg.build(weighted_points);
}
/// @brief Query the accumulated weight at point (x, y).
T query(int x, int y) {
T ans = 0;
ans = seg.prod(x + 1, y + 1);
for (auto &[l, r, d, u, w] : buckets) {
if (l <= x && x < r && d <= y && y < u) {
ans += w;
}
}
return ans;
}
};
} // namespace noya