library
This documentation is automatically generated by online-judge-tools/verification-helper
test/1_mytest/polygon_triangulation.test.cpp
- View this file on GitHub
- Last update: 2025-05-18 17:51:29+09:00
- Problem: https://judge.yosupo.jp/problem/aplusb
Depends on
- alg/monoid/add.hpp
- ds/fenwicktree/fenwicktree.hpp
- ds/fenwicktree/fenwicktree_01.hpp
- ds/hashmap.hpp
- ds/splaytree/splaytree.hpp
- ds/splaytree/splaytree_basic.hpp
- geo/angle_sort.hpp
- geo/base.hpp
- geo/convex_hull.hpp
- geo/count_points_in_triangles.hpp
- geo/cross_point.hpp
- geo/polygon_triangulation.hpp
- graph/base.hpp
- graph/planar_graph.hpp
- my_template.hpp
- random/base.hpp
- random/random_polygon.hpp
Code
#define PROBLEM "https://judge.yosupo.jp/problem/aplusb"
#include "my_template.hpp"
#include "random/random_polygon.hpp"
#include "geo/polygon_triangulation.hpp"
void test() {
auto check = [&](vc<Point<ll>> point) -> void {
int N = len(point);
auto dat = polygon_triangulation(point);
assert(len(dat) == N - 2);
// 簡易テスト. 面積和がいい感じの N-2 個になってればいいことにする.
ll AREA = 0, AREA3 = 0;
FOR(i, N) AREA += point[i].det(point[(i + 1) % N]);
for (auto& [a, b, c]: dat) {
ll S = (point[b] - point[a]).det(point[c] - point[a]);
assert(S > 0);
AREA3 += S;
}
assert(AREA == AREA3);
};
FOR(10000) {
int N = RNG(3, 20);
int K = RNG(3, 10);
vc<Point<ll>> point = random_polygon(N, K);
check(point);
}
}
void solve() {
int a, b;
cin >> a >> b;
cout << a + b << "\n";
}
signed main() {
test();
solve();
return 0;
}#line 1 "test/1_mytest/polygon_triangulation.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/aplusb"
#line 1 "my_template.hpp"
#if defined(LOCAL)
#include <my_template_compiled.hpp>
#else
// https://codeforces.com/blog/entry/96344
// https://codeforces.com/blog/entry/126772?#comment-1154880
#include <bits/allocator.h>
#pragma GCC optimize("Ofast,unroll-loops")
#pragma GCC target("avx2,popcnt")
#include <bits/stdc++.h>
using namespace std;
using ll = long long;
using u8 = uint8_t;
using u16 = uint16_t;
using u32 = uint32_t;
using u64 = uint64_t;
using i128 = __int128;
using u128 = unsigned __int128;
using f128 = __float128;
template <class T>
constexpr T infty = 0;
template <>
constexpr int infty<int> = 1'010'000'000;
template <>
constexpr ll infty<ll> = 2'020'000'000'000'000'000;
template <>
constexpr u32 infty<u32> = infty<int>;
template <>
constexpr u64 infty<u64> = infty<ll>;
template <>
constexpr i128 infty<i128> = i128(infty<ll>) * 2'000'000'000'000'000'000;
template <>
constexpr double infty<double> = infty<ll>;
template <>
constexpr long double infty<long double> = infty<ll>;
using pi = pair<ll, ll>;
using vi = vector<ll>;
template <class T>
using vc = vector<T>;
template <class T>
using vvc = vector<vc<T>>;
template <class T>
using vvvc = vector<vvc<T>>;
template <class T>
using vvvvc = vector<vvvc<T>>;
template <class T>
using vvvvvc = vector<vvvvc<T>>;
template <class T>
using pq = priority_queue<T>;
template <class T>
using pqg = priority_queue<T, vector<T>, greater<T>>;
#define vv(type, name, h, ...) vector<vector<type>> name(h, vector<type>(__VA_ARGS__))
#define vvv(type, name, h, w, ...) vector<vector<vector<type>>> name(h, vector<vector<type>>(w, vector<type>(__VA_ARGS__)))
#define vvvv(type, name, a, b, c, ...) \
vector<vector<vector<vector<type>>>> name(a, vector<vector<vector<type>>>(b, vector<vector<type>>(c, vector<type>(__VA_ARGS__))))
// https://trap.jp/post/1224/
#define FOR1(a) for (ll _ = 0; _ < ll(a); ++_)
#define FOR2(i, a) for (ll i = 0; i < ll(a); ++i)
#define FOR3(i, a, b) for (ll i = a; i < ll(b); ++i)
#define FOR4(i, a, b, c) for (ll i = a; i < ll(b); i += (c))
#define FOR1_R(a) for (ll i = (a)-1; i >= ll(0); --i)
#define FOR2_R(i, a) for (ll i = (a)-1; i >= ll(0); --i)
#define FOR3_R(i, a, b) for (ll i = (b)-1; i >= ll(a); --i)
#define overload4(a, b, c, d, e, ...) e
#define overload3(a, b, c, d, ...) d
#define FOR(...) overload4(__VA_ARGS__, FOR4, FOR3, FOR2, FOR1)(__VA_ARGS__)
#define FOR_R(...) overload3(__VA_ARGS__, FOR3_R, FOR2_R, FOR1_R)(__VA_ARGS__)
#define all(x) x.begin(), x.end()
#define len(x) ll(x.size())
#define elif else if
#define eb emplace_back
#define mp make_pair
#define mt make_tuple
#define fi first
#define se second
#define stoi stoll
int popcnt(int x) { return __builtin_popcount(x); }
int popcnt(u32 x) { return __builtin_popcount(x); }
int popcnt(ll x) { return __builtin_popcountll(x); }
int popcnt(u64 x) { return __builtin_popcountll(x); }
int popcnt_sgn(int x) { return (__builtin_parity(unsigned(x)) & 1 ? -1 : 1); }
int popcnt_sgn(u32 x) { return (__builtin_parity(x) & 1 ? -1 : 1); }
int popcnt_sgn(ll x) { return (__builtin_parityll(x) & 1 ? -1 : 1); }
int popcnt_sgn(u64 x) { return (__builtin_parityll(x) & 1 ? -1 : 1); }
// (0, 1, 2, 3, 4) -> (-1, 0, 1, 1, 2)
int topbit(int x) { return (x == 0 ? -1 : 31 - __builtin_clz(x)); }
int topbit(u32 x) { return (x == 0 ? -1 : 31 - __builtin_clz(x)); }
int topbit(ll x) { return (x == 0 ? -1 : 63 - __builtin_clzll(x)); }
int topbit(u64 x) { return (x == 0 ? -1 : 63 - __builtin_clzll(x)); }
// (0, 1, 2, 3, 4) -> (-1, 0, 1, 0, 2)
int lowbit(int x) { return (x == 0 ? -1 : __builtin_ctz(x)); }
int lowbit(u32 x) { return (x == 0 ? -1 : __builtin_ctz(x)); }
int lowbit(ll x) { return (x == 0 ? -1 : __builtin_ctzll(x)); }
int lowbit(u64 x) { return (x == 0 ? -1 : __builtin_ctzll(x)); }
template <typename T>
T kth_bit(int k) {
return T(1) << k;
}
template <typename T>
bool has_kth_bit(T x, int k) {
return x >> k & 1;
}
template <typename UINT>
struct all_bit {
struct iter {
UINT s;
iter(UINT s) : s(s) {}
int operator*() const { return lowbit(s); }
iter &operator++() {
s &= s - 1;
return *this;
}
bool operator!=(const iter) const { return s != 0; }
};
UINT s;
all_bit(UINT s) : s(s) {}
iter begin() const { return iter(s); }
iter end() const { return iter(0); }
};
template <typename UINT>
struct all_subset {
static_assert(is_unsigned<UINT>::value);
struct iter {
UINT s, t;
bool ed;
iter(UINT s) : s(s), t(s), ed(0) {}
int operator*() const { return s ^ t; }
iter &operator++() {
(t == 0 ? ed = 1 : t = (t - 1) & s);
return *this;
}
bool operator!=(const iter) const { return !ed; }
};
UINT s;
all_subset(UINT s) : s(s) {}
iter begin() const { return iter(s); }
iter end() const { return iter(0); }
};
template <typename T>
T floor(T a, T b) {
return a / b - (a % b && (a ^ b) < 0);
}
template <typename T>
T ceil(T x, T y) {
return floor(x + y - 1, y);
}
template <typename T>
T bmod(T x, T y) {
return x - y * floor(x, y);
}
template <typename T>
pair<T, T> divmod(T x, T y) {
T q = floor(x, y);
return {q, x - q * y};
}
template <typename T, typename U>
T SUM(const vector<U> &A) {
T sm = 0;
for (auto &&a: A) sm += a;
return sm;
}
#define MIN(v) *min_element(all(v))
#define MAX(v) *max_element(all(v))
#define LB(c, x) distance((c).begin(), lower_bound(all(c), (x)))
#define UB(c, x) distance((c).begin(), upper_bound(all(c), (x)))
#define UNIQUE(x) sort(all(x)), x.erase(unique(all(x)), x.end()), x.shrink_to_fit()
template <typename T>
T POP(deque<T> &que) {
T a = que.front();
que.pop_front();
return a;
}
template <typename T>
T POP(pq<T> &que) {
T a = que.top();
que.pop();
return a;
}
template <typename T>
T POP(pqg<T> &que) {
T a = que.top();
que.pop();
return a;
}
template <typename T>
T POP(vc<T> &que) {
T a = que.back();
que.pop_back();
return a;
}
template <typename F>
ll binary_search(F check, ll ok, ll ng, bool check_ok = true) {
if (check_ok) assert(check(ok));
while (abs(ok - ng) > 1) {
auto x = (ng + ok) / 2;
(check(x) ? ok : ng) = x;
}
return ok;
}
template <typename F>
double binary_search_real(F check, double ok, double ng, int iter = 100) {
FOR(iter) {
double x = (ok + ng) / 2;
(check(x) ? ok : ng) = x;
}
return (ok + ng) / 2;
}
template <class T, class S>
inline bool chmax(T &a, const S &b) {
return (a < b ? a = b, 1 : 0);
}
template <class T, class S>
inline bool chmin(T &a, const S &b) {
return (a > b ? a = b, 1 : 0);
}
// ? は -1
vc<int> s_to_vi(const string &S, char first_char) {
vc<int> A(S.size());
FOR(i, S.size()) { A[i] = (S[i] != '?' ? S[i] - first_char : -1); }
return A;
}
template <typename T, typename U>
vector<T> cumsum(vector<U> &A, int off = 1) {
int N = A.size();
vector<T> B(N + 1);
FOR(i, N) { B[i + 1] = B[i] + A[i]; }
if (off == 0) B.erase(B.begin());
return B;
}
// stable sort
template <typename T>
vector<int> argsort(const vector<T> &A) {
vector<int> ids(len(A));
iota(all(ids), 0);
sort(all(ids), [&](int i, int j) { return (A[i] == A[j] ? i < j : A[i] < A[j]); });
return ids;
}
// A[I[0]], A[I[1]], ...
template <typename T>
vc<T> rearrange(const vc<T> &A, const vc<int> &I) {
vc<T> B(len(I));
FOR(i, len(I)) B[i] = A[I[i]];
return B;
}
template <typename T, typename... Vectors>
void concat(vc<T> &first, const Vectors &... others) {
vc<T> &res = first;
(res.insert(res.end(), others.begin(), others.end()), ...);
}
#endif
#line 3 "test/1_mytest/polygon_triangulation.test.cpp"
#line 1 "random/random_polygon.hpp"
#line 2 "random/base.hpp"
u64 RNG_64() {
static u64 x_ = u64(chrono::duration_cast<chrono::nanoseconds>(chrono::high_resolution_clock::now().time_since_epoch()).count()) * 10150724397891781847ULL;
x_ ^= x_ << 7;
return x_ ^= x_ >> 9;
}
u64 RNG(u64 lim) { return RNG_64() % lim; }
ll RNG(ll l, ll r) { return l + RNG_64() % (r - l); }
#line 2 "geo/base.hpp"
template <typename T>
struct Point {
T x, y;
Point() : x(0), y(0) {}
template <typename A, typename B>
Point(A x, B y) : x(x), y(y) {}
template <typename A, typename B>
Point(pair<A, B> p) : x(p.fi), y(p.se) {}
Point operator+=(const Point p) {
x += p.x, y += p.y;
return *this;
}
Point operator-=(const Point p) {
x -= p.x, y -= p.y;
return *this;
}
Point operator+(Point p) const { return {x + p.x, y + p.y}; }
Point operator-(Point p) const { return {x - p.x, y - p.y}; }
bool operator==(Point p) const { return x == p.x && y == p.y; }
bool operator!=(Point p) const { return x != p.x || y != p.y; }
Point operator-() const { return {-x, -y}; }
Point operator*(T t) const { return {x * t, y * t}; }
Point operator/(T t) const { return {x / t, y / t}; }
bool operator<(Point p) const {
if (x != p.x) return x < p.x;
return y < p.y;
}
T dot(const Point& other) const { return x * other.x + y * other.y; }
T det(const Point& other) const { return x * other.y - y * other.x; }
double norm() { return sqrtl(x * x + y * y); }
double angle() { return atan2(y, x); }
Point rotate(double theta) {
static_assert(!is_integral<T>::value);
double c = cos(theta), s = sin(theta);
return Point{c * x - s * y, s * x + c * y};
}
Point rot90(bool ccw) { return (ccw ? Point{-y, x} : Point{y, -x}); }
};
#ifdef FASTIO
template <typename T>
void rd(Point<T>& p) {
fastio::rd(p.x), fastio::rd(p.y);
}
template <typename T>
void wt(Point<T>& p) {
fastio::wt(p.x);
fastio::wt(' ');
fastio::wt(p.y);
}
#endif
// A -> B -> C と進むときに、左に曲がるならば +1、右に曲がるならば -1
template <typename T>
int ccw(Point<T> A, Point<T> B, Point<T> C) {
T x = (B - A).det(C - A);
if (x > 0) return 1;
if (x < 0) return -1;
return 0;
}
template <typename REAL, typename T, typename U>
REAL dist(Point<T> A, Point<U> B) {
REAL dx = REAL(A.x) - REAL(B.x);
REAL dy = REAL(A.y) - REAL(B.y);
return sqrt(dx * dx + dy * dy);
}
// ax+by+c
template <typename T>
struct Line {
T a, b, c;
Line(T a, T b, T c) : a(a), b(b), c(c) {}
Line(Point<T> A, Point<T> B) { a = A.y - B.y, b = B.x - A.x, c = A.x * B.y - A.y * B.x; }
Line(T x1, T y1, T x2, T y2) : Line(Point<T>(x1, y1), Point<T>(x2, y2)) {}
template <typename U>
U eval(Point<U> P) {
return U(a) * P.x + U(b) * P.y + U(c);
}
template <typename U>
T eval(U x, U y) {
return a * x + b * y + c;
}
// 同じ直線が同じ a,b,c で表現されるようにする
void normalize() {
static_assert(is_same_v<T, int> || is_same_v<T, long long>);
T g = gcd(gcd(abs(a), abs(b)), abs(c));
a /= g, b /= g, c /= g;
if (b < 0) { a = -a, b = -b, c = -c; }
if (b == 0 && a < 0) { a = -a, b = -b, c = -c; }
}
bool is_parallel(Line other) { return a * other.b - b * other.a == 0; }
bool is_orthogonal(Line other) { return a * other.a + b * other.b == 0; }
};
template <typename T>
struct Segment {
Point<T> A, B;
Segment(Point<T> A, Point<T> B) : A(A), B(B) {}
Segment(T x1, T y1, T x2, T y2) : Segment(Point<T>(x1, y1), Point<T>(x2, y2)) {}
bool contain(Point<T> C) {
T det = (C - A).det(B - A);
if (det != 0) return 0;
return (C - A).dot(B - A) >= 0 && (C - B).dot(A - B) >= 0;
}
Line<T> to_Line() { return Line(A, B); }
};
template <typename REAL>
struct Circle {
Point<REAL> O;
REAL r;
Circle() {}
Circle(Point<REAL> O, REAL r) : O(O), r(r) {}
Circle(REAL x, REAL y, REAL r) : O(x, y), r(r) {}
template <typename T>
bool contain(Point<T> p) {
REAL dx = p.x - O.x, dy = p.y - O.y;
return dx * dx + dy * dy <= r * r;
}
};
#line 2 "geo/convex_hull.hpp"
#line 4 "geo/convex_hull.hpp"
// allow_180=true で同一座標点があるとこわれる
// full なら I[0] が sorted で min になる
template <typename T, bool allow_180 = false>
vector<int> ConvexHull(vector<Point<T>>& XY, string mode = "full", bool sorted = false) {
assert(mode == "full" || mode == "lower" || mode == "upper");
ll N = XY.size();
if (N == 1) return {0};
if (N == 2) {
if (XY[0] < XY[1]) return {0, 1};
if (XY[1] < XY[0]) return {1, 0};
return {0};
}
vc<int> I(N);
if (sorted) {
FOR(i, N) I[i] = i;
} else {
I = argsort(XY);
}
if constexpr (allow_180) { FOR(i, N - 1) assert(XY[i] != XY[i + 1]); }
auto check = [&](ll i, ll j, ll k) -> bool {
T det = (XY[j] - XY[i]).det(XY[k] - XY[i]);
if constexpr (allow_180) return det >= 0;
return det > T(0);
};
auto calc = [&]() {
vector<int> P;
for (auto&& k: I) {
while (P.size() > 1) {
auto i = P[P.size() - 2];
auto j = P[P.size() - 1];
if (check(i, j, k)) break;
P.pop_back();
}
P.eb(k);
}
return P;
};
vc<int> P;
if (mode == "full" || mode == "lower") {
vc<int> Q = calc();
P.insert(P.end(), all(Q));
}
if (mode == "full" || mode == "upper") {
if (!P.empty()) P.pop_back();
reverse(all(I));
vc<int> Q = calc();
P.insert(P.end(), all(Q));
}
if (mode == "upper") reverse(all(P));
while (len(P) >= 2 && XY[P[0]] == XY[P.back()]) P.pop_back();
return P;
}
#line 2 "geo/cross_point.hpp"
#line 4 "geo/cross_point.hpp"
// 平行でないことを仮定
template <typename REAL, typename T>
Point<REAL> cross_point(const Line<T> L1, const Line<T> L2) {
T det = L1.a * L2.b - L1.b * L2.a;
assert(det != 0);
REAL x = -REAL(L1.c) * L2.b + REAL(L1.b) * L2.c;
REAL y = -REAL(L1.a) * L2.c + REAL(L1.c) * L2.a;
return Point<REAL>(x / det, y / det);
}
// 浮動小数点数はエラー
// 0: 交点なし
// 1: 一意な交点
// 2：2 つ以上の交点（整数型を利用して厳密にやる）
template <typename T>
int count_cross(Segment<T> S1, Segment<T> S2, bool include_ends) {
static_assert(!std::is_floating_point<T>::value);
Line<T> L1 = S1.to_Line();
Line<T> L2 = S2.to_Line();
if (L1.is_parallel(L2)) {
if (L1.eval(S2.A) != 0) return 0;
// 4 点とも同一直線上にある
T a1 = S1.A.x, b1 = S1.B.x;
T a2 = S2.A.x, b2 = S2.B.x;
if (a1 == b1) {
a1 = S1.A.y, b1 = S1.B.y;
a2 = S2.A.y, b2 = S2.B.y;
}
if (a1 > b1) swap(a1, b1);
if (a2 > b2) swap(a2, b2);
T a = max(a1, a2);
T b = min(b1, b2);
if (a < b) return 2;
if (a > b) return 0;
return (include_ends ? 1 : 0);
}
// 平行でない場合
T a1 = L2.eval(S1.A), b1 = L2.eval(S1.B);
T a2 = L1.eval(S2.A), b2 = L1.eval(S2.B);
if (a1 > b1) swap(a1, b1);
if (a2 > b2) swap(a2, b2);
bool ok1 = 0, ok2 = 0;
if (include_ends) {
ok1 = (a1 <= T(0)) && (T(0) <= b1);
ok2 = (a2 <= T(0)) && (T(0) <= b2);
} else {
ok1 = (a1 < T(0)) && (T(0) < b1);
ok2 = (a2 < T(0)) && (T(0) < b2);
}
return (ok1 && ok2 ? 1 : 0);
}
// 4 次式まで登場している、オーバーフロー注意！
// https://codeforces.com/contest/607/problem/E
template <typename REAL, typename T>
vc<Point<REAL>> cross_point(const Circle<T> C, const Line<T> L) {
T a = L.a, b = L.b, c = L.a * (C.O.x) + L.b * (C.O.y) + L.c;
T r = C.r;
bool SW = 0;
T abs_a = (a < 0 ? -a : a);
T abs_b = (b < 0 ? -b : b);
if (abs_a < abs_b) {
swap(a, b);
SW = 1;
}
// ax+by+c=0, x^2+y^2=r^2
T D = 4 * c * c * b * b - 4 * (a * a + b * b) * (c * c - a * a * r * r);
if (D < 0) return {};
REAL sqD = sqrtl(D);
REAL y1 = (-2 * b * c + sqD) / (2 * (a * a + b * b));
REAL y2 = (-2 * b * c - sqD) / (2 * (a * a + b * b));
REAL x1 = (-b * y1 - c) / a;
REAL x2 = (-b * y2 - c) / a;
if (SW) swap(x1, y1), swap(x2, y2);
x1 += C.O.x, x2 += C.O.x;
y1 += C.O.y, y2 += C.O.y;
if (D == 0) return {Point<REAL>(x1, y1)};
return {Point<REAL>(x1, y1), Point<REAL>(x2, y2)};
}
// https://codeforces.com/contest/2/problem/C
template <typename REAL, typename T>
tuple<bool, Point<T>, Point<T>> cross_point_circle(Circle<T> C1, Circle<T> C2) {
using P = Point<T>;
P O{0, 0};
P A = C1.O, B = C2.O;
if (A == B) return {false, O, O};
T d = (B - A).norm();
REAL cos_val = (C1.r * C1.r + d * d - C2.r * C2.r) / (2 * C1.r * d);
if (cos_val < -1 || 1 < cos_val) return {false, O, O};
REAL t = acos(cos_val);
REAL u = (B - A).angle();
P X = A + P{C1.r * cos(u + t), C1.r * sin(u + t)};
P Y = A + P{C1.r * cos(u - t), C1.r * sin(u - t)};
return {true, X, Y};
}
#line 1 "geo/count_points_in_triangles.hpp"
#line 2 "geo/angle_sort.hpp"
#line 4 "geo/angle_sort.hpp"
// lower: -1, origin: 0, upper: 1, (-pi,pi]
template <typename T> int lower_or_upper(const Point<T> &p) {
if (p.y != 0)
return (p.y > 0 ? 1 : -1);
if (p.x > 0)
return -1;
if (p.x < 0)
return 1;
return 0;
}
// L<R:-1, L==R:0, L>R:1, (-pi,pi]
template <typename T> int angle_comp_3(const Point<T> &L, const Point<T> &R) {
int a = lower_or_upper(L), b = lower_or_upper(R);
if (a != b)
return (a < b ? -1 : +1);
T det = L.det(R);
if (det > 0)
return -1;
if (det < 0)
return 1;
return 0;
}
// 偏角ソートに対する argsort, (-pi,pi]
template <typename T> vector<int> angle_sort(vector<Point<T>> &P) {
vc<int> I(len(P));
FOR(i, len(P)) I[i] = i;
sort(all(I), [&](auto &L, auto &R) -> bool {
return angle_comp_3(P[L], P[R]) == -1;
});
return I;
}
// 偏角ソートに対する argsort, (-pi,pi]
template <typename T> vector<int> angle_sort(vector<pair<T, T>> &P) {
vc<Point<T>> tmp(len(P));
FOR(i, len(P)) tmp[i] = Point<T>(P[i]);
return angle_sort<T>(tmp);
}
#line 2 "ds/fenwicktree/fenwicktree_01.hpp"
#line 2 "alg/monoid/add.hpp"
template <typename E>
struct Monoid_Add {
using X = E;
using value_type = X;
static constexpr X op(const X &x, const X &y) noexcept { return x + y; }
static constexpr X inverse(const X &x) noexcept { return -x; }
static constexpr X power(const X &x, ll n) noexcept { return X(n) * x; }
static constexpr X unit() { return X(0); }
static constexpr bool commute = true;
};
#line 3 "ds/fenwicktree/fenwicktree.hpp"
template <typename Monoid>
struct FenwickTree {
using G = Monoid;
using MX = Monoid;
using E = typename G::value_type;
int n;
vector<E> dat;
E total;
FenwickTree() {}
FenwickTree(int n) { build(n); }
template <typename F>
FenwickTree(int n, F f) {
build(n, f);
}
FenwickTree(const vc<E>& v) { build(v); }
void build(int m) {
n = m;
dat.assign(m, G::unit());
total = G::unit();
}
void build(const vc<E>& v) {
build(len(v), [&](int i) -> E { return v[i]; });
}
template <typename F>
void build(int m, F f) {
n = m;
dat.clear();
dat.reserve(n);
total = G::unit();
FOR(i, n) { dat.eb(f(i)); }
for (int i = 1; i <= n; ++i) {
int j = i + (i & -i);
if (j <= n) dat[j - 1] = G::op(dat[i - 1], dat[j - 1]);
}
total = prefix_sum(m);
}
E prod_all() { return total; }
E sum_all() { return total; }
E sum(int k) { return prefix_sum(k); }
E prod(int k) { return prefix_prod(k); }
E prefix_sum(int k) { return prefix_prod(k); }
E prefix_prod(int k) {
chmin(k, n);
E ret = G::unit();
for (; k > 0; k -= k & -k) ret = G::op(ret, dat[k - 1]);
return ret;
}
E sum(int L, int R) { return prod(L, R); }
E prod(int L, int R) {
chmax(L, 0), chmin(R, n);
if (L == 0) return prefix_prod(R);
assert(0 <= L && L <= R && R <= n);
E pos = G::unit(), neg = G::unit();
while (L < R) { pos = G::op(pos, dat[R - 1]), R -= R & -R; }
while (R < L) { neg = G::op(neg, dat[L - 1]), L -= L & -L; }
return G::op(pos, G::inverse(neg));
}
vc<E> get_all() {
vc<E> res(n);
FOR(i, n) res[i] = prod(i, i + 1);
return res;
}
void add(int k, E x) { multiply(k, x); }
void multiply(int k, E x) {
static_assert(G::commute);
total = G::op(total, x);
for (++k; k <= n; k += k & -k) dat[k - 1] = G::op(dat[k - 1], x);
}
void set(int k, E x) { add(k, G::op(G::inverse(prod(k, k + 1)), x)); }
template <class F>
int max_right(const F check, int L = 0) {
assert(check(G::unit()));
E s = G::unit();
int i = L;
// 2^k 進むとダメ
int k = [&]() {
while (1) {
if (i % 2 == 1) { s = G::op(s, G::inverse(dat[i - 1])), i -= 1; }
if (i == 0) { return topbit(n) + 1; }
int k = lowbit(i) - 1;
if (i + (1 << k) > n) return k;
E t = G::op(s, dat[i + (1 << k) - 1]);
if (!check(t)) { return k; }
s = G::op(s, G::inverse(dat[i - 1])), i -= i & -i;
}
}();
while (k) {
--k;
if (i + (1 << k) - 1 < len(dat)) {
E t = G::op(s, dat[i + (1 << k) - 1]);
if (check(t)) { i += (1 << k), s = t; }
}
}
return i;
}
// check(i, x)
template <class F>
int max_right_with_index(const F check, int L = 0) {
assert(check(L, G::unit()));
E s = G::unit();
int i = L;
// 2^k 進むとダメ
int k = [&]() {
while (1) {
if (i % 2 == 1) { s = G::op(s, G::inverse(dat[i - 1])), i -= 1; }
if (i == 0) { return topbit(n) + 1; }
int k = lowbit(i) - 1;
if (i + (1 << k) > n) return k;
E t = G::op(s, dat[i + (1 << k) - 1]);
if (!check(i + (1 << k), t)) { return k; }
s = G::op(s, G::inverse(dat[i - 1])), i -= i & -i;
}
}();
while (k) {
--k;
if (i + (1 << k) - 1 < len(dat)) {
E t = G::op(s, dat[i + (1 << k) - 1]);
if (check(i + (1 << k), t)) { i += (1 << k), s = t; }
}
}
return i;
}
template <class F>
int min_left(const F check, int R) {
assert(check(G::unit()));
E s = G::unit();
int i = R;
// false になるところまで戻る
int k = 0;
while (i > 0 && check(s)) {
s = G::op(s, dat[i - 1]);
k = lowbit(i);
i -= i & -i;
}
if (check(s)) {
assert(i == 0);
return 0;
}
// 2^k 進むと ok になる
// false を維持して進む
while (k) {
--k;
E t = G::op(s, G::inverse(dat[i + (1 << k) - 1]));
if (!check(t)) { i += (1 << k), s = t; }
}
return i + 1;
}
int kth(E k, int L = 0) {
return max_right([&k](E x) -> bool { return x <= k; }, L);
}
};
#line 4 "ds/fenwicktree/fenwicktree_01.hpp"
struct FenwickTree_01 {
using MX = Monoid_Add<int>;
int N, n;
vc<u64> dat;
FenwickTree<Monoid_Add<int>> bit;
FenwickTree_01() {}
FenwickTree_01(int n) { build(n); }
template <typename F>
FenwickTree_01(int n, F f) {
build(n, f);
}
void build(int m) {
N = m;
n = ceil<int>(N + 1, 64);
dat.assign(n, u64(0));
bit.build(n);
}
void build(vc<int> dat) {
build(len(dat), [&](int i) -> int { return dat[i]; });
}
template <typename F>
void build(int m, F f) {
N = m;
n = ceil<int>(N + 1, 64);
dat.assign(n, u64(0));
FOR(i, N) { dat[i / 64] |= u64(f(i)) << (i % 64); }
bit.build(n, [&](int i) -> int { return popcnt(dat[i]); });
}
int sum_all() { return bit.sum_all(); }
int sum(int k) { return prefix_sum(k); }
int prefix_sum(int k) {
int ans = bit.sum(k / 64);
ans += popcnt(dat[k / 64] & ((u64(1) << (k % 64)) - 1));
return ans;
}
int sum(int L, int R) {
if (L == 0) return prefix_sum(R);
int ans = 0;
ans -= popcnt(dat[L / 64] & ((u64(1) << (L % 64)) - 1));
ans += popcnt(dat[R / 64] & ((u64(1) << (R % 64)) - 1));
ans += bit.sum(L / 64, R / 64);
return ans;
}
int prod(int L, int R) { return sum(L, R); }
void add(int k, int x) {
if (x == 1) add(k);
elif (x == -1) remove(k);
else assert(0);
}
void multiply(int k, int x) { add(k, x); }
void add(int k) {
dat[k / 64] |= u64(1) << (k % 64);
bit.add(k / 64, 1);
}
void remove(int k) {
dat[k / 64] &= ~(u64(1) << (k % 64));
bit.add(k / 64, -1);
}
int kth(int k, int L = 0) {
if (k >= sum_all()) return N;
k += popcnt(dat[L / 64] & ((u64(1) << (L % 64)) - 1));
L /= 64;
int mid = 0;
auto check = [&](auto e) -> bool {
if (e <= k) chmax(mid, e);
return e <= k;
};
int idx = bit.max_right(check, L);
if (idx == n) return N;
k -= mid;
u64 x = dat[idx];
int p = popcnt(x);
if (p <= k) return N;
k = binary_search([&](int n) -> bool { return (p - popcnt(x >> n)) <= k; }, 0, 64, 0);
return 64 * idx + k;
}
int next(int k) {
int idx = k / 64;
k %= 64;
u64 x = dat[idx] & ~((u64(1) << k) - 1);
if (x) return 64 * idx + lowbit(x);
idx = bit.kth(0, idx + 1);
if (idx == n || !dat[idx]) return N;
return 64 * idx + lowbit(dat[idx]);
}
int prev(int k) {
if (k == N) --k;
int idx = k / 64;
k %= 64;
u64 x = dat[idx];
if (k < 63) x &= (u64(1) << (k + 1)) - 1;
if (x) return 64 * idx + topbit(x);
idx = bit.min_left([&](auto e) -> bool { return e <= 0; }, idx) - 1;
if (idx == -1) return -1;
return 64 * idx + topbit(dat[idx]);
}
};
#line 6 "geo/count_points_in_triangles.hpp"
// 点群 A, B を入力 （Point<ll>）
// query(i,j,k)：三角形 AiAjAk 内部の Bl の個数（非負）を返す
// 前計算 O(NMlogM)、クエリ O(1)
// https://codeforces.com/contest/13/problem/D
// https://codeforces.com/contest/852/problem/H
struct Count_Points_In_Triangles {
using P = Point<ll>;
const int LIM = 1'000'000'000 + 10;
vc<P> A, B;
vc<int> new_idx; // O から見た偏角ソート順を管理
vc<int> point; // A[i] と一致する B[j] の数え上げ
vvc<int> seg; // 線分 A[i]A[j] 内にある B[k] の数え上げ
vvc<int> tri; // OA[i]A[j] 内部にある B[k] の数え上げ
Count_Points_In_Triangles(const vc<P>& A, const vc<P>& B) : A(A), B(B) {
for (auto&& p: A) assert(max(abs(p.x), abs(p.y)) < LIM);
for (auto&& p: B) assert(max(abs(p.x), abs(p.y)) < LIM);
build();
}
int count3(int i, int j, int k) {
i = new_idx[i], j = new_idx[j], k = new_idx[k];
if (i > j) swap(i, j);
if (j > k) swap(j, k);
if (i > j) swap(i, j);
assert(i <= j && j <= k);
ll d = (A[j] - A[i]).det(A[k] - A[i]);
if (d == 0) return 0;
if (d > 0) { return tri[i][j] + tri[j][k] - tri[i][k] - seg[i][k]; }
int x = tri[i][k] - tri[i][j] - tri[j][k];
return x - seg[i][j] - seg[j][k] - point[j];
}
// segment
int count2(int i, int j) {
i = new_idx[i], j = new_idx[j];
if (i > j) swap(i, j);
return seg[i][j];
}
private:
P take_origin() {
// OAiAj, OAiBj が同一直線上にならないようにする
// fail prob: at most N(N+M)/LIM
return P{-LIM, RNG(-LIM, LIM)};
}
void build() {
P O = take_origin();
for (auto&& p: A) p = p - O;
for (auto&& p: B) p = p - O;
int N = len(A), M = len(B);
vc<int> I = angle_sort(A);
A = rearrange(A, I);
new_idx.resize(N);
FOR(i, N) new_idx[I[i]] = i;
I = angle_sort(B);
B = rearrange(B, I);
point.assign(N, 0);
seg.assign(N, vc<int>(N));
tri.assign(N, vc<int>(N));
// point
FOR(i, N) FOR(j, M) if (A[i] == B[j])++ point[i];
int m = 0;
FOR(j, N) {
// OA[i]A[j], B[k]
while (m < M && A[j].det(B[m]) < 0) ++m;
vc<P> C(m);
FOR(k, m) C[k] = B[k] - A[j];
vc<int> I(m);
FOR(i, m) I[i] = i;
sort(all(I), [&](auto& a, auto& b) -> bool { return C[a].det(C[b]) > 0; });
C = rearrange(C, I);
vc<int> rk(m);
FOR(k, m) rk[I[k]] = k;
FenwickTree_01 bit(m);
int k = m;
FOR_R(i, j) {
while (k > 0 && A[i].det(B[k - 1]) > 0) { bit.add(rk[--k], 1); }
P p = A[i] - A[j];
int lb = binary_search([&](int n) -> bool { return (n == 0 ? true : C[n - 1].det(p) > 0); }, 0, m + 1);
int ub = binary_search([&](int n) -> bool { return (n == 0 ? true : C[n - 1].det(p) >= 0); }, 0, m + 1);
seg[i][j] += bit.sum(lb, ub), tri[i][j] += bit.sum(lb);
}
}
}
};
#line 7 "random/random_polygon.hpp"
vc<Point<ll>> random_polygon(int N, int XY_ABS_MAX = 10) {
assert(N >= 3);
using P = Point<ll>;
auto trial = [&]() -> vc<P> {
set<Point<ll>> S;
while (len(S) < N) {
int x = RNG(-XY_ABS_MAX, XY_ABS_MAX + 1);
int y = RNG(-XY_ABS_MAX, XY_ABS_MAX + 1);
S.insert(Point<ll>(x, y));
}
vc<P> point(all(S));
auto I = ConvexHull<ll, true>(point);
Count_Points_In_Triangles CT(point, point);
vc<int> other;
vc<int> done(N);
for (auto& i: I) done[i]++;
if (MAX(done) >= 2) return {};
FOR(i, N) if (!done[i]) other.eb(i);
int fail = 0;
while (len(other)) {
if (fail > 1000) return {};
++fail;
int i = RNG(0, len(I)), j = RNG(0, len(other));
swap(other[j], other.back());
int a = I[i], b = I[(i + 1) % len(I)], c = other.back();
if ((point[b] - point[a]).det(point[c] - point[a]) < 0) continue;
if (CT.count3(a, b, c)) continue;
if (CT.count2(a, c) + CT.count2(b, c)) continue;
bool ok = 1;
for (auto& v: {a, b}) {
FOR(i, len(I)) {
Segment<ll> S1(point[v], point[c]);
Segment<ll> S2(point[I[i]], point[I[(i + 1) % len(I)]]);
if (count_cross(S1, S2, false)) ok = 0;
}
}
if (!ok) continue;
fail = 0;
I.insert(I.begin() + i + 1, POP(other));
}
point = rearrange(point, I);
FOR(i, N) {
if ((point[(i + 2) % N] - point[i]).det(point[(i + 1) % N] - point[i]) == 0) return {};
}
return point;
};
while (1) {
vc<P> ANS = trial();
if (ANS.empty()) continue;
int k = RNG(0, len(ANS));
rotate(ANS.begin(), ANS.begin() + k, ANS.end());
return ANS;
}
}
#line 2 "ds/splaytree/splaytree.hpp"
/*
update でちゃんと prod が計算されてくれれば prod は op(lprod,x,rprod) でなくてもよい.
*/
// Node 型を別に定義して使う
template <typename Node>
struct SplayTree {
Node *pool;
const int NODES;
int pid;
using np = Node *;
using X = typename Node::value_type;
using A = typename Node::operator_type;
vc<np> FREE;
SplayTree(int NODES) : NODES(NODES), pid(0) { pool = new Node[NODES]; }
~SplayTree() { delete[] pool; }
void free_subtree(np c) {
if (!c) return;
auto dfs = [&](auto &dfs, np c) -> void {
if (c->l) dfs(dfs, c->l);
if (c->r) dfs(dfs, c->r);
FREE.eb(c);
};
dfs(dfs, c);
}
void reset() {
pid = 0;
FREE.clear();
}
np new_root() { return nullptr; }
np new_node(const X &x) {
assert(!FREE.empty() || pid < NODES);
np n = (FREE.empty() ? &(pool[pid++]) : POP(FREE));
Node::new_node(n, x);
return n;
}
np new_node(const vc<X> &dat) {
auto dfs = [&](auto &dfs, int l, int r) -> np {
if (l == r) return nullptr;
if (r == l + 1) return new_node(dat[l]);
int m = (l + r) / 2;
np l_root = dfs(dfs, l, m);
np r_root = dfs(dfs, m + 1, r);
np root = new_node(dat[m]);
root->l = l_root, root->r = r_root;
if (l_root) l_root->p = root;
if (r_root) r_root->p = root;
root->update();
return root;
};
return dfs(dfs, 0, len(dat));
}
u32 get_size(np root) { return (root ? root->size : 0); }
np merge(np l_root, np r_root) {
if (!l_root) return r_root;
if (!r_root) return l_root;
assert((!l_root->p) && (!r_root->p));
splay_kth(r_root, 0); // splay したので prop 済
r_root->l = l_root;
l_root->p = r_root;
r_root->update();
return r_root;
}
np merge3(np a, np b, np c) { return merge(merge(a, b), c); }
np merge4(np a, np b, np c, np d) { return merge(merge(merge(a, b), c), d); }
pair<np, np> split(np root, u32 k) {
assert(!root || !root->p);
if (k == 0) return {nullptr, root};
if (k == (root->size)) return {root, nullptr};
splay_kth(root, k - 1);
np right = root->r;
root->r = nullptr, right->p = nullptr;
root->update();
return {root, right};
}
tuple<np, np, np> split3(np root, u32 l, u32 r) {
np nm, nr;
tie(root, nr) = split(root, r);
tie(root, nm) = split(root, l);
return {root, nm, nr};
}
tuple<np, np, np, np> split4(np root, u32 i, u32 j, u32 k) {
np d;
tie(root, d) = split(root, k);
auto [a, b, c] = split3(root, i, j);
return {a, b, c, d};
}
tuple<np, np, np> split_L_root_R(np root) {
u32 s = (root->l ? root->l->size : 0);
return split3(root, s, s + 1);
}
// 部分木が区間 [l,r) に対応するようなノードを作って返す
// そのノードが root になるわけではないので、
// このノードを参照した後にすぐに splay して根に持ち上げること
void goto_between(np &root, u32 l, u32 r) {
if (l == 0 && r == root->size) return;
if (l == 0) {
splay_kth(root, r);
root = root->l;
return;
}
if (r == root->size) {
splay_kth(root, l - 1);
root = root->r;
return;
}
splay_kth(root, r);
np rp = root;
root = rp->l;
root->p = nullptr;
splay_kth(root, l - 1);
root->p = rp;
rp->l = root;
rp->update();
root = root->r;
}
vc<X> get_all(const np &root) {
vc<X> res;
auto dfs = [&](auto &dfs, np root) -> void {
if (!root) return;
root->prop();
dfs(dfs, root->l);
res.eb(root->get());
dfs(dfs, root->r);
};
dfs(dfs, root);
return res;
}
X get(np &root, u32 k) {
assert(root == nullptr || !root->p);
splay_kth(root, k);
return root->get();
}
void set(np &root, u32 k, const X &x) {
assert(root != nullptr && !root->p);
splay_kth(root, k);
root->set(x);
}
void multiply(np &root, u32 k, const X &x) {
assert(root != nullptr && !root->p);
splay_kth(root, k);
root->multiply(x);
}
X prod(np &root, u32 l, u32 r) {
assert(root == nullptr || !root->p);
using Mono = typename Node::Monoid_X;
if (l == r) return Mono::unit();
assert(0 <= l && l < r && r <= root->size);
goto_between(root, l, r);
X res = root->prod;
splay(root, true);
return res;
}
X prod(np &root) {
assert(root == nullptr || !root->p);
using Mono = typename Node::Monoid_X;
return (root ? root->prod : Mono::unit());
}
void apply(np &root, u32 l, u32 r, const A &a) {
if (l == r) return;
assert(0 <= l && l < r && r <= root->size);
goto_between(root, l, r);
root->apply(a);
splay(root, true);
}
void apply(np &root, const A &a) {
if (!root) return;
root->apply(a);
}
void reverse(np &root, u32 l, u32 r) {
assert(root == nullptr || !root->p);
if (l == r) return;
assert(0 <= l && l < r && r <= root->size);
goto_between(root, l, r);
root->reverse();
splay(root, true);
}
void reverse(np root) {
if (!root) return;
root->reverse();
}
void rotate(Node *n) {
// n を根に近づける。prop, update は rotate の外で行う。
Node *pp, *p, *c;
p = n->p;
pp = p->p;
if (p->l == n) {
c = n->r;
n->r = p;
p->l = c;
} else {
c = n->l;
n->l = p;
p->r = c;
}
if (pp && pp->l == p) pp->l = n;
if (pp && pp->r == p) pp->r = n;
n->p = pp;
p->p = n;
if (c) c->p = p;
}
void prop_from_root(np c) {
if (!c->p) {
c->prop();
return;
}
prop_from_root(c->p);
c->prop();
}
void splay(Node *me, bool prop_from_root_done) {
// これを呼ぶ時点で、me の祖先（me を除く）は既に prop 済であることを仮定
// 特に、splay 終了時点で me は upd / prop 済である
if (!prop_from_root_done) prop_from_root(me);
me->prop();
while (me->p) {
np p = me->p;
np pp = p->p;
if (!pp) {
rotate(me);
p->update();
break;
}
bool same = (p->l == me && pp->l == p) || (p->r == me && pp->r == p);
if (same) rotate(p), rotate(me);
if (!same) rotate(me), rotate(me);
pp->update(), p->update();
}
// me の update は最後だけでよい
me->update();
}
void splay_kth(np &root, u32 k) {
assert(0 <= k && k < (root->size));
while (1) {
root->prop();
u32 sl = (root->l ? root->l->size : 0);
if (k == sl) break;
if (k < sl)
root = root->l;
else {
k -= sl + 1;
root = root->r;
}
}
splay(root, true);
}
// check(x), 左側のノード全体が check を満たすように切る
template <typename F>
pair<np, np> split_max_right(np root, F check) {
if (!root) return {nullptr, nullptr};
assert(!root->p);
np c = find_max_right(root, check);
if (!c) {
splay(root, true);
return {nullptr, root};
}
splay(c, true);
np right = c->r;
if (!right) return {c, nullptr};
right->p = nullptr;
c->r = nullptr;
c->update();
return {c, right};
}
// check(x, cnt), 左側のノード全体が check を満たすように切る
template <typename F>
pair<np, np> split_max_right_cnt(np root, F check) {
if (!root) return {nullptr, nullptr};
assert(!root->p);
np c = find_max_right_cnt(root, check);
if (!c) {
splay(root, true);
return {nullptr, root};
}
splay(c, true);
np right = c->r;
if (!right) return {c, nullptr};
right->p = nullptr;
c->r = nullptr;
c->update();
return {c, right};
}
// 左側のノード全体の prod が check を満たすように切る
template <typename F>
pair<np, np> split_max_right_prod(np root, F check) {
if (!root) return {nullptr, nullptr};
assert(!root->p);
np c = find_max_right_prod(root, check);
if (!c) {
splay(root, true);
return {nullptr, root};
}
splay(c, true);
np right = c->r;
if (!right) return {c, nullptr};
right->p = nullptr;
c->r = nullptr;
c->update();
return {c, right};
}
template <typename F>
np find_max_right(np root, const F &check) {
// 最後に見つけた ok の点、最後に探索した点
np last_ok = nullptr, last = nullptr;
while (root) {
last = root;
root->prop();
if (check(root->x)) {
last_ok = root;
root = root->r;
} else {
root = root->l;
}
}
splay(last, true);
return last_ok;
}
template <typename F>
np find_max_right_cnt(np root, const F &check) {
// 最後に見つけた ok の点、最後に探索した点
np last_ok = nullptr, last = nullptr;
ll n = 0;
while (root) {
last = root;
root->prop();
ll ns = (root->l ? root->l->size : 0);
if (check(root->x, n + ns + 1)) {
last_ok = root;
n += ns + 1;
root = root->r;
} else {
root = root->l;
}
}
splay(last, true);
return last_ok;
}
template <typename F>
np find_max_right_prod(np root, const F &check) {
using Mono = typename Node::Monoid_X;
X prod = Mono::unit();
// 最後に見つけた ok の点、最後に探索した点
np last_ok = nullptr, last = nullptr;
while (root) {
last = root;
root->prop();
np tmp = root->r;
root->r = nullptr;
root->update();
X lprod = Mono::op(prod, root->prod);
root->r = tmp;
root->update();
if (check(lprod)) {
prod = lprod;
last_ok = root;
root = root->r;
} else {
root = root->l;
}
}
splay(last, true);
return last_ok;
}
};
#line 2 "ds/splaytree/splaytree_basic.hpp"
namespace SplayTreeNodes {
template <typename S>
struct Node_Basic {
using value_type = S;
using operator_type = int;
using np = Node_Basic *;
np p, l, r;
bool rev;
S x;
u32 size;
static void new_node(np n, const S &x) {
n->p = n->l = n->r = nullptr;
n->x = x, n->size = 1, n->rev = 0;
}
void update() {
size = 1;
if (l) { size += l->size; }
if (r) { size += r->size; }
}
void prop() {
if (rev) {
if (l) { l->rev ^= 1, swap(l->l, l->r); }
if (r) { r->rev ^= 1, swap(r->l, r->r); }
rev = 0;
}
}
// update, prop 以外で呼ばれるものは、splay 後であることが想定されている。
// したがってその時点で update, prop 済であることを仮定してよい。
S get() { return x; }
void set(const S &xx) {
x = xx;
update();
}
void reverse() {
swap(l, r);
rev ^= 1;
}
};
template <typename S>
using SplayTree_Basic = SplayTree<Node_Basic<S>>;
} // namespace SplayTreeNodes
using SplayTreeNodes::SplayTree_Basic;
#line 2 "ds/hashmap.hpp"
// u64 -> Val
template <typename Val>
struct HashMap {
// n は入れたいものの個数で ok
HashMap(u32 n = 0) { build(n); }
void build(u32 n) {
u32 k = 8;
while (k < n * 2) k *= 2;
cap = k / 2, mask = k - 1;
key.resize(k), val.resize(k), used.assign(k, 0);
}
// size を保ったまま. size=0 にするときは build すること.
void clear() {
used.assign(len(used), 0);
cap = (mask + 1) / 2;
}
int size() { return len(used) / 2 - cap; }
int index(const u64& k) {
int i = 0;
for (i = hash(k); used[i] && key[i] != k; i = (i + 1) & mask) {}
return i;
}
Val& operator[](const u64& k) {
if (cap == 0) extend();
int i = index(k);
if (!used[i]) { used[i] = 1, key[i] = k, val[i] = Val{}, --cap; }
return val[i];
}
Val get(const u64& k, Val default_value) {
int i = index(k);
return (used[i] ? val[i] : default_value);
}
bool count(const u64& k) {
int i = index(k);
return used[i] && key[i] == k;
}
// f(key, val)
template <typename F>
void enumerate_all(F f) {
FOR(i, len(used)) if (used[i]) f(key[i], val[i]);
}
private:
u32 cap, mask;
vc<u64> key;
vc<Val> val;
vc<bool> used;
u64 hash(u64 x) {
static const u64 FIXED_RANDOM = std::chrono::steady_clock::now().time_since_epoch().count();
x += FIXED_RANDOM;
x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9;
x = (x ^ (x >> 27)) * 0x94d049bb133111eb;
return (x ^ (x >> 31)) & mask;
}
void extend() {
vc<pair<u64, Val>> dat;
dat.reserve(len(used) / 2 - cap);
FOR(i, len(used)) {
if (used[i]) dat.eb(key[i], val[i]);
}
build(2 * len(dat));
for (auto& [a, b]: dat) (*this)[a] = b;
}
};
#line 3 "graph/base.hpp"
template <typename T>
struct Edge {
int frm, to;
T cost;
int id;
};
template <typename T = int, bool directed = false>
struct Graph {
static constexpr bool is_directed = directed;
int N, M;
using cost_type = T;
using edge_type = Edge<T>;
vector<edge_type> edges;
vector<int> indptr;
vector<edge_type> csr_edges;
vc<int> vc_deg, vc_indeg, vc_outdeg;
bool prepared;
class OutgoingEdges {
public:
OutgoingEdges(const Graph* G, int l, int r) : G(G), l(l), r(r) {}
const edge_type* begin() const {
if (l == r) { return 0; }
return &G->csr_edges[l];
}
const edge_type* end() const {
if (l == r) { return 0; }
return &G->csr_edges[r];
}
private:
const Graph* G;
int l, r;
};
bool is_prepared() { return prepared; }
Graph() : N(0), M(0), prepared(0) {}
Graph(int N) : N(N), M(0), prepared(0) {}
void build(int n) {
N = n, M = 0;
prepared = 0;
edges.clear();
indptr.clear();
csr_edges.clear();
vc_deg.clear();
vc_indeg.clear();
vc_outdeg.clear();
}
void add(int frm, int to, T cost = 1, int i = -1) {
assert(!prepared);
assert(0 <= frm && 0 <= to && to < N);
if (i == -1) i = M;
auto e = edge_type({frm, to, cost, i});
edges.eb(e);
++M;
}
#ifdef FASTIO
// wt, off
void read_tree(bool wt = false, int off = 1) { read_graph(N - 1, wt, off); }
void read_graph(int M, bool wt = false, int off = 1) {
for (int m = 0; m < M; ++m) {
INT(a, b);
a -= off, b -= off;
if (!wt) {
add(a, b);
} else {
T c;
read(c);
add(a, b, c);
}
}
build();
}
#endif
void build() {
assert(!prepared);
prepared = true;
indptr.assign(N + 1, 0);
for (auto&& e: edges) {
indptr[e.frm + 1]++;
if (!directed) indptr[e.to + 1]++;
}
for (int v = 0; v < N; ++v) { indptr[v + 1] += indptr[v]; }
auto counter = indptr;
csr_edges.resize(indptr.back() + 1);
for (auto&& e: edges) {
csr_edges[counter[e.frm]++] = e;
if (!directed) csr_edges[counter[e.to]++] = edge_type({e.to, e.frm, e.cost, e.id});
}
}
OutgoingEdges operator[](int v) const {
assert(prepared);
return {this, indptr[v], indptr[v + 1]};
}
vc<int> deg_array() {
if (vc_deg.empty()) calc_deg();
return vc_deg;
}
pair<vc<int>, vc<int>> deg_array_inout() {
if (vc_indeg.empty()) calc_deg_inout();
return {vc_indeg, vc_outdeg};
}
int deg(int v) {
if (vc_deg.empty()) calc_deg();
return vc_deg[v];
}
int in_deg(int v) {
if (vc_indeg.empty()) calc_deg_inout();
return vc_indeg[v];
}
int out_deg(int v) {
if (vc_outdeg.empty()) calc_deg_inout();
return vc_outdeg[v];
}
#ifdef FASTIO
void debug() {
#ifdef LOCAL
print("Graph");
if (!prepared) {
print("frm to cost id");
for (auto&& e: edges) print(e.frm, e.to, e.cost, e.id);
} else {
print("indptr", indptr);
print("frm to cost id");
FOR(v, N) for (auto&& e: (*this)[v]) print(e.frm, e.to, e.cost, e.id);
}
#endif
}
#endif
vc<int> new_idx;
vc<bool> used_e;
// G における頂点 V[i] が、新しいグラフで i になるようにする
// {G, es}
// sum(deg(v)) の計算量になっていて、
// 新しいグラフの n+m より大きい可能性があるので注意
Graph<T, directed> rearrange(vc<int> V, bool keep_eid = 0) {
if (len(new_idx) != N) new_idx.assign(N, -1);
int n = len(V);
FOR(i, n) new_idx[V[i]] = i;
Graph<T, directed> G(n);
vc<int> history;
FOR(i, n) {
for (auto&& e: (*this)[V[i]]) {
if (len(used_e) <= e.id) used_e.resize(e.id + 1);
if (used_e[e.id]) continue;
int a = e.frm, b = e.to;
if (new_idx[a] != -1 && new_idx[b] != -1) {
history.eb(e.id);
used_e[e.id] = 1;
int eid = (keep_eid ? e.id : -1);
G.add(new_idx[a], new_idx[b], e.cost, eid);
}
}
}
FOR(i, n) new_idx[V[i]] = -1;
for (auto&& eid: history) used_e[eid] = 0;
G.build();
return G;
}
Graph<T, true> to_directed_tree(int root = -1) {
if (root == -1) root = 0;
assert(!is_directed && prepared && M == N - 1);
Graph<T, true> G1(N);
vc<int> par(N, -1);
auto dfs = [&](auto& dfs, int v) -> void {
for (auto& e: (*this)[v]) {
if (e.to == par[v]) continue;
par[e.to] = v, dfs(dfs, e.to);
}
};
dfs(dfs, root);
for (auto& e: edges) {
int a = e.frm, b = e.to;
if (par[a] == b) swap(a, b);
assert(par[b] == a);
G1.add(a, b, e.cost);
}
G1.build();
return G1;
}
HashMap<int> MP_FOR_EID;
int get_eid(u64 a, u64 b) {
if (len(MP_FOR_EID) == 0) {
MP_FOR_EID.build(N - 1);
for (auto& e: edges) {
u64 a = e.frm, b = e.to;
u64 k = to_eid_key(a, b);
MP_FOR_EID[k] = e.id;
}
}
return MP_FOR_EID.get(to_eid_key(a, b), -1);
}
u64 to_eid_key(u64 a, u64 b) {
if (!directed && a > b) swap(a, b);
return N * a + b;
}
private:
void calc_deg() {
assert(vc_deg.empty());
vc_deg.resize(N);
for (auto&& e: edges) vc_deg[e.frm]++, vc_deg[e.to]++;
}
void calc_deg_inout() {
assert(vc_indeg.empty());
vc_indeg.resize(N);
vc_outdeg.resize(N);
for (auto&& e: edges) { vc_indeg[e.to]++, vc_outdeg[e.frm]++; }
}
};
#line 4 "graph/planar_graph.hpp"
/*
・連結平面グラフになっていないときにどう動作するかは何も考えていない
・N=1 も扱わない
・0番目に外面が入る
・次数 1 の点とかはあっても大丈夫っぽい？
*/
template <typename XY>
struct Planar_Graph {
using P = Point<XY>;
int NV, NE, NF;
// 頂点, 辺からなるグラフ. 有向辺を 2 つ入れておく
Graph<int, 1> G;
// 頂点属性
vc<P> point; // 座標
// 辺属性
vc<int> left_face; // 有向辺の左にある面の番号
vc<int> nxt_edge; // 面を反時計回りにまわるときの次の辺
// 面属性
vc<int> first_edge;
Planar_Graph(int N, vc<P> point) : NV(N), G(N), point(point) { assert(N > 1); }
void add(int a, int b) { G.add(a, b), G.add(b, a); }
void build() {
G.build();
NE = G.M / 2;
nxt_edge.assign(G.M, -1);
left_face.assign(G.M, -1);
int v0 = 0;
int e0 = 0;
FOR(v, NV) {
if (point[v] < point[v0]) v0 = v;
vc<int> eid;
vc<P> dir;
for (auto& e: G[v]) {
eid.eb(e.id);
dir.eb(point[e.to] - point[e.frm]);
}
auto I = angle_sort(dir);
assert(len(I) > 0);
FOR(k, len(I)) {
int i = (k == 0 ? I.back() : I[k - 1]);
int j = I[k];
i = eid[i], j = eid[j];
nxt_edge[j ^ 1] = i;
}
if (v == v0) e0 = eid[I[0]] ^ 1;
}
for (auto& x: nxt_edge) assert(x != -1);
auto make_face = [&](int e) -> void {
int p = len(first_edge);
first_edge.eb(e);
while (left_face[e] == -1) {
left_face[e] = p;
e = nxt_edge[e];
}
};
make_face(e0);
FOR(e, 2 * NE) {
if (left_face[e] == -1) make_face(e);
}
NF = len(first_edge);
assert(NV - NE + NF == 2);
}
// return {vs, es}
// vs = [v0,v1,v2,v0], es = [e0,e1,e2]
pair<vc<int>, vc<int>> get_face_data(int fid) {
vc<int> eid = {first_edge[fid]};
while (1) {
int e = nxt_edge[eid.back()];
if (e == first_edge[fid]) break;
eid.eb(e);
}
vc<int> vid;
for (auto& e: eid) vid.eb(G.edges[e].frm);
vid.eb(vid[0]);
return {vid, eid};
}
};
#line 4 "geo/polygon_triangulation.hpp"
template <typename T>
vc<tuple<int, int, int>> monotone_polygon_triangulation(vc<Point<T>> point) {
int N = len(point);
int rot = min_element(all(point)) - point.begin();
rotate(point.begin(), point.begin() + rot, point.end());
int n = max_element(all(point)) - point.begin();
FOR(i, n) assert(point[i] < point[i + 1]);
FOR(i, n, N) assert(point[(i + 1) % N] < point[i]);
vc<tuple<int, int, int>> res;
auto side = [&](int i) -> int {
assert(i != 0 && i != n);
return (i < n ? 0 : 1);
};
vc<int> I = argsort(point);
vc<int> stack = {I[0], I[1]};
int s = side(I[1]);
FOR(i, 2, N - 1) {
int v = I[i], t = side(v);
if (s == 0 && t == 0) {
while (len(stack) >= 2 && ccw(point[stack[len(stack) - 2]], point[stack[len(stack) - 1]], point[v]) == 1) {
res.eb(stack[len(stack) - 2], stack[len(stack) - 1], v), POP(stack);
}
stack.eb(v);
}
elif (s == 1 && t == 1) {
while (len(stack) >= 2 && ccw(point[stack[len(stack) - 2]], point[stack[len(stack) - 1]], point[v]) == -1) {
res.eb(stack[len(stack) - 2], v, stack[len(stack) - 1]), POP(stack);
}
stack.eb(v);
}
elif (s == 0 && t == 1) {
FOR(j, len(stack) - 1) res.eb(stack[j], stack[j + 1], v);
stack = {stack.back(), v}, s = t;
}
elif (s == 1 && t == 0) {
FOR(j, len(stack) - 1) res.eb(stack[j], v, stack[j + 1]);
stack = {stack.back(), v}, s = t;
}
}
if (s == 0) { FOR(j, len(stack) - 1) res.eb(stack[j], stack[j + 1], n); }
elif (s == 1) { FOR(j, len(stack) - 1) res.eb(stack[j], n, stack[j + 1]); }
for (auto& [a, b, c]: res) a = (a + rot) % N, b = (b + rot) % N, c = (c + rot) % N;
return res;
}
// (i,j,k), ccw
template <typename T>
vc<tuple<int, int, int>> polygon_triangulation(vc<Point<T>> point) {
using P = Point<T>;
int N = len(point);
enum vtype { MERGE, SPLIT, START, END, UPPER, LOWER };
auto pre = [&](int i) -> int { return (i > 0 ? i - 1 : N - 1); };
auto nxt = [&](int i) -> int { return (i < N - 1 ? i + 1 : 0); };
auto get_vtype = [&](int i) -> vtype {
int l = pre(i), r = nxt(i);
if (point[i] < point[l] && point[i] < point[r]) { return (ccw(point[l], point[i], point[r]) == 1 ? START : SPLIT); }
if (point[l] < point[i] && point[r] < point[i]) { return (ccw(point[l], point[i], point[r]) == 1 ? END : MERGE); }
if (point[l] < point[i] && point[i] < point[r]) return LOWER;
if (point[r] < point[i] && point[i] < point[l]) return UPPER;
assert(0);
return END;
};
SplayTree_Basic<int> ST(N);
using np = decltype(ST)::np;
vc<np> nodes(N);
FOR(i, N) nodes[i] = ST.new_node(i);
np S = ST.new_root();
auto comp = [&](int i, P p) -> bool {
P A = point[i], B = point[nxt(i)];
return ccw(A, B, p) == -1;
};
vc<int> helper(N, -1);
vc<bool> merged(N);
Planar_Graph<T> G(N, point);
FOR(i, N) G.add(i, nxt(i));
auto add_edge = [&](int v, int w) -> void { merged[w] = 1, G.add(v, w); };
auto fix_up = [&](int v, int e) -> void {
int w = helper[e];
if (get_vtype(w) == vtype::MERGE && !merged[w]) { add_edge(v, w); }
};
auto I = argsort(point);
for (auto& i: I) {
vtype t = get_vtype(i);
if (t == vtype::MERGE) {
ST.splay(nodes[i], 1), S = nodes[i];
int n = (nodes[i]->l ? nodes[i]->l->size : 0);
auto [L, M, R] = ST.split3(S, n, n + 1);
int j = ST.get(R, 0);
S = ST.merge(L, R);
fix_up(i, i), fix_up(i, j);
helper[j] = i;
}
if (t == vtype::SPLIT) {
auto [L, R] = ST.split_max_right(S, [&](int k) -> bool { return comp(k, point[i]); });
int j = ST.get(R, 0);
add_edge(i, helper[j]);
helper[j] = i, helper[pre(i)] = i;
S = ST.merge3(L, nodes[pre(i)], R);
}
if (t == vtype::START) {
auto [L, R] = ST.split_max_right(S, [&](int k) -> bool { return comp(k, point[i]); });
S = ST.merge3(L, nodes[pre(i)], R), helper[pre(i)] = i;
}
if (t == vtype::END) {
ST.splay(nodes[i], 1), S = nodes[i];
int n = (nodes[i]->l ? nodes[i]->l->size : 0);
auto [L, M, R] = ST.split3(S, n, n + 1);
S = ST.merge(L, R);
fix_up(i, i);
}
if (t == vtype::UPPER) {
ST.splay(nodes[i], 1), S = nodes[i];
int n = (nodes[i]->l ? nodes[i]->l->size : 0);
auto [L, M, R] = ST.split3(S, n, n + 1);
S = ST.merge3(L, nodes[pre(i)], R);
fix_up(i, i);
helper[pre(i)] = i;
}
if (t == vtype::LOWER) {
auto [L, R] = ST.split_max_right(S, [&](int k) -> bool { return comp(k, point[i]); });
int j = ST.get(R, 0);
S = ST.merge(L, R);
fix_up(i, j);
helper[j] = i;
}
}
G.build();
vc<tuple<int, int, int>> ANS;
FOR(f, 1, G.NF) {
auto [vs, es] = G.get_face_data(f);
POP(vs);
vc<P> sub = rearrange(point, vs);
for (auto& [a, b, c]: monotone_polygon_triangulation(sub)) ANS.eb(vs[a], vs[b], vs[c]);
}
return ANS;
}
#line 6 "test/1_mytest/polygon_triangulation.test.cpp"
void test() {
auto check = [&](vc<Point<ll>> point) -> void {
int N = len(point);
auto dat = polygon_triangulation(point);
assert(len(dat) == N - 2);
// 簡易テスト. 面積和がいい感じの N-2 個になってればいいことにする.
ll AREA = 0, AREA3 = 0;
FOR(i, N) AREA += point[i].det(point[(i + 1) % N]);
for (auto& [a, b, c]: dat) {
ll S = (point[b] - point[a]).det(point[c] - point[a]);
assert(S > 0);
AREA3 += S;
}
assert(AREA == AREA3);
};
FOR(10000) {
int N = RNG(3, 20);
int K = RNG(3, 10);
vc<Point<ll>> point = random_polygon(N, K);
check(point);
}
}
void solve() {
int a, b;
cin >> a >> b;
cout << a + b << "\n";
}
signed main() {
test();
solve();
return 0;
}