Hello, codeforces!
Today, based on the article Triangulation of polygon from a book "Computational Geometry: Algorithms and applications", I implemented an algorithm of triangulation of any polygon without self intersections (non-convex polygons too) in C++.
$$$\text{ }$$$ $$$\text{ }$$$
Why it was needed for me? Probably, you know about a competition RuCode 7.0. In this year, qualification round contained next problem H: you are given an endless plane, painted in black and white stripes with height 1 and infinity width. These stripes alternate: black, white, black, white... You are given a polygon on this plane as list of $$$n\text{ }\left(3\le n \le 100000\right)$$$ 2D points with integer coordinates $$$x_i, y_i \left(-10^9 \le x_i, y_i \le 10^9\right)$$$ in some traversal order. Need to calculate a total square inside this polygon which is painted in black color.
Firstly, I solved for any triangle ($$$n=3$$$). Secondly, I thought about the triangulation of polygon and accumulating answers for each triangle.
So, the heaviest part of solution is the triangulation of polygon. I implemented this, stress-tested and fixed all of incomprehensible and dubious places of pseudo code that can be found in above book. For example, the author of this book states that vertex which has a type "merge" should be connected with another "merge" vertex, and this is used in pseudo-code, but there is an example, where it is not true. Stress-testing shows all of corner cases, which are not covered by the pseudo code.
I fixed all of mistakes which the stress-testing has found. You can just call this function:
std::vector<Triangle> triangulate(vpii points) {
makeCounterClockwise(points);
return triangulateFast(points);
}
The output of this function is a vector of triangles. Be careful, because I'm erasing non-meaningful vertices like if $$$P_{i-1}$$$, $$$P_{i}$$$ and $$$P_{i+1}$$$ lie on the same line than erase $$$P_i$$$.
Source code of solution for given problem
#pragma GCC optimize("Ofast")
#ifndef __TEMPLATE_HPP__
#define __TEMPLATE_HPP__
#include <climits>
#include <unordered_map>
#include <random>
#include <chrono>
#include <numeric>
#include <iostream>
#include <vector>
#include <algorithm>
#include <map>
#include <queue>
#include <deque>
#include <stack>
#include <functional>
#include <bitset>
#include <string>
#include <sstream>
#include <fstream>
#include <iomanip>
#include <cmath>
#include <cassert>
#include <list>
#include <forward_list>
#include <set>
#include <unordered_set>
#include <cstdint>
#define all(x) std::begin(x), std::end(x)
#define isz(x) (int)std::size(x)
// marco for random variable name:
#define CONCAT(a, b) CONCAT_INNER(a, b)
#define CONCAT_INNER(a, b) a ## b
#define SomeVar CONCAT(var, CONCAT(_, CONCAT(__LINE__, CONCAT(_, __COUNTER__))))
const auto ready = [](){
std::ios_base::sync_with_stdio(0); std::cin.tie(0);
std::cout << std::fixed << std::setprecision(12);
return 1;
}();
// Defines:
using ll = long long;
using ull = unsigned long long;
using ld = long double;
using vi = std::vector<int>;
using vl = std::vector<ll>;
using vvi = std::vector<vi>;
using vvl = std::vector<vl>;
using pii = std::pair<int,int>;
using pil = std::pair<int,ll>;
using pli = std::pair<ll,int>;
using pll = std::pair<ll,ll>;
using vpii = std::vector<pii>;
using vvpii = std::vector<vpii>;
using vpll = std::vector<pll>;
using vvpll = std::vector<vpll>;
using vs = std::vector<std::string>;
using tiii = std::tuple<int,int,int>;
using tiil = std::tuple<int,int,ll>;
using vtiii = std::vector<tiii>;
using vtiil = std::vector<tiil>;
// Comparators:
#define GEN_COMPARATORS(CLASS) \
inline bool operator >(const CLASS& a, const CLASS& b) { return b < a; } \
inline bool operator<=(const CLASS& a, const CLASS& b) { return !(a > b); } \
inline bool operator>=(const CLASS& a, const CLASS& b) { return !(a < b); } \
inline bool operator!=(const CLASS& a, const CLASS& b) { return a < b || b < a; } \
inline bool operator==(const CLASS& a, const CLASS& b) { return !(a != b); }
#define GEN_COMPARATORS_MEMBERS(CLASS) \
bool operator >(const CLASS &other) const { return other < (*this); } \
bool operator<=(const CLASS &other) const { return !((*this) > other); } \
bool operator>=(const CLASS &other) const { return !((*this) < other); } \
bool operator!=(const CLASS &other) const { return (*this) < other || other < (*this); } \
bool operator==(const CLASS &other) const { return !((*this) != other); }
namespace std {
#if __cplusplus < 201703L
// Containers:
template <typename T> constexpr auto size(const T& t) -> decltype(t.size()) { return t.size(); }
// 1D arrays:
template<typename T, int N> auto size(const T (&)[N]) { return N; }
#endif
// 2D arrays:
template<typename T, int N, int M> auto size(const T (&)[N][M]) { return N * M; }
template<typename T, int N, int M> auto begin(T (&a)[N][M]) { return &a[0][0]; }
template<typename T, int N, int M> auto end(T (&a)[N][M]) { return &a[0][0] + N * M; }
// 3D arrays:
template<typename T, int N, int M, int K> auto size(const T (&)[N][M][K]) { return N * M * K; }
template<typename T, int N, int M, int K> auto begin(T (&a)[N][M][K]) { return &a[0][0][0]; }
template<typename T, int N, int M, int K> auto end(T (&a)[N][M][K]) { return &a[0][0][0] + N * M * K; }
}
// Algorithms:
template<typename C> void reuniq(C& c) { c.erase(unique(all(c)), end(c)); }
template<typename C, typename X> int lowpos(const C& c, X x) {
return int(lower_bound(all(c), x) - begin(c));
}
template<typename C, typename X> int uppos(const C& c, X x) {
return int(upper_bound(all(c), x) - begin(c));
}
template<typename X, typename Y> X& remin(X& x, const Y& y) { return x = (y < x) ? y : x; }
template<typename X, typename Y> X& remax(X& x, const Y& y) { return x = (x < y) ? y : x; }
// Input:
template<typename T> std::istream& operator>>(std::istream& is, std::vector<T>& vec) {
for (auto &it : vec) is >> it;
return is;
}
// ---- ---- ---- ---- ---- ---- OPERATORS FOR STL CONTAINERS ---- ---- ---- ---- ---- ----
#define INSERT(cont, to_front, to_back) \
template<typename X, typename Y, typename... T> cont<X,T...>& operator<<(cont<X,T...>& c, const Y& x) { return to_back, c; } \
template<typename X, typename Y, typename... T> cont<X,T...>& operator>>(const Y& x, cont<X,T...>& c) { return to_front, c; }
INSERT(std::vector, (c.insert(c.begin(), x)), (c.push_back(x)))
INSERT(std::queue, (c.push(x)), (c.push(x)))
INSERT(std::stack, (c.push(x)), (c.push(x)))
INSERT(std::priority_queue, (c.push(x)), (c.push(x)))
INSERT(std::deque, (c.push_front(x)), (c.push_back(x)))
INSERT(std::list, (c.push_front(x)), (c.push_back(x)))
INSERT(std::set, (c.insert(x)), (c.insert(x)))
INSERT(std::unordered_set, (c.insert(x)), (c.insert(x)))
INSERT(std::multiset, (c.insert(x)), (c.insert(x)))
#undef INSERT
#define REMOVE(cont, from_front, from_back) \
template<typename X, typename... T> X operator--(cont<X,T...>& c,int) { X x; from_back; return x; } \
template<typename X, typename... T> X operator--(cont<X,T...>& c) { X x; from_front; return x; }
REMOVE(std::vector,(x=c.front(),c.erase(c.begin())),(x=c.back(),c.pop_back()))
REMOVE(std::deque,(x=c.front(),c.pop_front()),(x=c.back(),c.pop_back()))
REMOVE(std::queue,(x=c.front(),c.pop()),(x=c.front(),c.pop()))
REMOVE(std::stack,(x=c.top(),c.pop()),(x=c.top(),c.pop()))
REMOVE(std::priority_queue,(x=c.top(),c.pop()),(x=c.top(),c.pop()))
REMOVE(std::list,(x=c.front(),c.pop_front()),(x=c.back(),c.pop_back()))
REMOVE(std::set,(x=*c.begin(),c.erase(c.begin())),(x=*std::prev(c.end()),c.erase(std::prev(c.end()))))
REMOVE(std::multiset,(x=*c.begin(),c.erase(c.begin())),(x=*std::prev(c.end()),c.erase(std::prev(c.end()))))
REMOVE(std::unordered_set,(x=*c.begin(),c.erase(c.begin())),(x=*std::prev(c.end()),c.erase(std::prev(c.end()))))
REMOVE(std::map,(x=*c.begin(),c.erase(c.begin())),(x=*std::prev(c.end()),c.erase(std::prev(c.end()))))
REMOVE(std::unordered_map,(x=*c.begin(),c.erase(c.begin())),(x=*std::prev(c.end()),c.erase(std::prev(c.end()))))
#undef REMOVE
#define EXTRACT(cont) \
template<typename X, typename... T> cont<X,T...>& operator>>(cont<X,T...>& c, X& x) { return x=c--,c; } \
template<typename X, typename... T> cont<X,T...>& operator<<(X& x, cont<X,T...>& c) { return x=--c,c; }
EXTRACT(std::vector) EXTRACT(std::list) EXTRACT(std::deque) EXTRACT(std::queue)
EXTRACT(std::stack) EXTRACT(std::priority_queue) EXTRACT(std::set)
EXTRACT(std::map) EXTRACT(std::multiset) EXTRACT(std::unordered_set)
#undef EXTRACT
#define MOVE(cont) \
template<typename X, typename... T> \
cont<X,T...>& operator>>(cont<X,T...>& l, cont<X,T...>& r) \
{{ while (isz(l)) l-- >> r;} return r; } \
template<typename X, typename... T> \
cont<X,T...>& operator<<(cont<X,T...>& l, cont<X,T...>& r) \
{{ while (isz(r)) l << --r;} return l; }
MOVE(std::list) MOVE(std::deque) MOVE(std::queue) MOVE(std::set) MOVE(std::multiset)
MOVE(std::stack) MOVE(std::priority_queue) MOVE(std::unordered_set)
#undef MOVE
template<typename X, typename... T>
std::vector<X,T...>& operator>>(std::vector<X,T...>& l, std::vector<X,T...>& r)
{ return r.insert(r.begin(),all(l)),l.clear(), r; }
template<typename X, typename... T>
std::vector<X,T...>& operator<<(std::vector<X,T...>& l, std::vector<X,T...>& r)
{ return l.insert(l.end(),all(r)),r.clear(), l; }
// Auto-revert:
template<typename F1, typename F2>
struct AutoRevert {
F1 f1; F2 f2;
AutoRevert(const F1 &f1_, const F2 &f2_) : f1(f1_), f2(f2_) { f1(); }
~AutoRevert() { f2(); }
};
template<typename T>
struct AutoRecover {
T &ref, val;
AutoRecover(T &x) : ref(x), val(x) { }
~AutoRecover() { ref = val; }
};
// Operations with bits:
template<typename T> void setbit(T &mask, int bit, bool x) { (mask &= ~(T(1) << bit)) |= (T(x) << bit); }
template<typename T> bool getbit(T &mask, int bit) { return (mask >> bit & 1); }
template<typename T> void flipbit(T &mask, int bit) { mask ^= (T(1) << bit); }
// Convert vector to tuple: auto [x,y,z] = to_tuple<3>(vec);
template <typename F, size_t... Is>
auto gen_tuple_impl(F func, std::index_sequence<Is...> ) { return std::make_tuple(func(Is)...); }
template <size_t N, typename F>
auto gen_tuple(F func) { return gen_tuple_impl(func, std::make_index_sequence<N>{} ); }
template<size_t N, typename ... A>
auto to_tuple(const std::vector<A...> &v) { return gen_tuple<N>([&v](size_t i){return v[i];});}
// pack / unpack vector: unpack(vec, x, y, z)
void unpack(const auto &vec, auto &...b) { int i = -1; ((b = vec[++i]),...); }
void pack(auto &vec, const auto &...b) { int i = -1; ((vec[++i] = b),...); }
// -----------------------------------------------------------------------------
#endif // __TEMPLATE_HPP__
#ifndef __DEBUG_HPP__
#define __DEBUG_HPP__
// ---- ---- ---- ---- ---- ---- DEBUG LIBRARY ---- ---- ---- ---- ---- ----
#define watch(...) debug && std::cerr << "{" << #__VA_ARGS__ << "} = " \
<< std::make_tuple(__VA_ARGS__) << std::endl
template<typename... X>
std::ostream& operator<<(std::ostream& os, const std::pair<X...>& p);
template<std::size_t I = 0, typename FuncT, typename... Tp>
inline typename std::enable_if<I == sizeof...(Tp), void>::type
for_each_const(const std::tuple<Tp...> &, FuncT) { }
template<std::size_t I = 0, typename FuncT, typename... Tp>
inline typename std::enable_if<I < sizeof...(Tp), void>::type
for_each_const(const std::tuple<Tp...>& t, FuncT f)
{ f(std::get<I>(t)),for_each_const<I + 1, FuncT, Tp...>(t, f); }
template<std::size_t I = 0, typename FuncT, typename... Tp>
inline typename std::enable_if<I == sizeof...(Tp), void>::type
for_each(std::tuple<Tp...> &, FuncT) { }
template<std::size_t I = 0, typename FuncT, typename... Tp>
inline typename std::enable_if<I < sizeof...(Tp), void>::type
for_each(std::tuple<Tp...>& t, FuncT f)
{ f(std::get<I>(t)); for_each<I + 1, FuncT, Tp...>(t, f); }
struct Printer {
std::ostream& os; bool was{0};
Printer(std::ostream& os_) : os(os_) { }
template<typename X> void operator()(X x);
};
template<typename Iterator>
std::ostream& print(std::ostream& os, Iterator begin, Iterator end)
{ return os << "{", std::for_each(begin,end,Printer(os)), os << "}"; }
#define OUTPUT(container) template<typename X, typename... T> \
std::ostream& operator<<(std::ostream& os, const container<X,T...>& c) \
{ return print(os, all(c)); }
OUTPUT(std::vector) OUTPUT(std::list) OUTPUT(std::deque)
OUTPUT(std::set) OUTPUT(std::unordered_set)
OUTPUT(std::multiset) OUTPUT(std::unordered_multiset)
OUTPUT(std::map) OUTPUT(std::multimap) OUTPUT(std::unordered_map)
#undef OUTPUT
#define OUTPUT2(container, get, pop) template<typename X, typename... T> \
std::ostream& operator<<(std::ostream& os, container<X,T...> c) { \
std::vector<X> v(c.size()); \
for (unsigned i = 0; i != v.size(); v[i++] = c.get(),c.pop()); \
return os << v; }
OUTPUT2(std::queue,front,pop)
OUTPUT2(std::stack,top,pop)
OUTPUT2(std::priority_queue,top,pop)
#undef OUTPUT2
template<typename... X>
std::ostream& operator<<(std::ostream& os, const std::tuple<X...>& t)
{ return os << "{", for_each_const(t, Printer(os)), os << "}"; }
template<typename X>
void Printer::operator()(X x)
{ os << (was?", ":(was=1,"")) << x; }
template<typename... X>
std::ostream& operator<<(std::ostream& os, const std::pair<X...>& p)
{ return os << std::make_tuple(std::get<0>(p), std::get<1>(p)); }
int debug = 0;
#endif // __DEBUG_HPP__
#ifndef __GEOMA_HPP__
#define __GEOMA_HPP__
namespace algos {
namespace geoma {
template<typename A, typename B>
inline std::pair<A,B> operator-(const std::pair<A,B> &a, const std::pair<A,B> &b) {
return std::pair<A,B>(a.first - b.first, a.second - b.second);
}
template<typename A, typename B>
inline std::pair<A,B> operator+(const std::pair<A,B> &a, const std::pair<A,B> &b) {
return std::pair<A,B>(a.first + b.first, a.second + b.second);
}
template<typename T, typename A, typename B>
inline T cross(const std::pair<A,B> &a, const std::pair<A,B> &b) {
return T(a.first) * b.second - T(a.second) * b.first;
}
template<typename T, typename A, typename B>
inline T doubleSquare(const std::pair<A,B> &a, const std::pair<A,B> &b, const std::pair<A,B> &c) {
T res = cross<T>(a-b, c-b);
if (res < 0)
res = -res;
return res;
}
template<typename T, typename A, typename B>
inline T dist2(const std::pair<A,B> &a) {
return T(a.first)*a.first + T(a.second)*a.second;
}
template<typename T, typename A, typename B>
inline T dist2(const std::pair<A, B> &a, const std::pair<A, B> &b) {
return dist2<T>(a-b);
}
using Triangle = std::tuple<pii,pii,pii>;
inline bool isPointInsideTriangle(Triangle tr, pii p) {
auto [A, B, C] = tr;
auto S = doubleSquare<ll>(A,B,C);
auto S1 = doubleSquare<ll>(p,B,C);
auto S2 = doubleSquare<ll>(A,p,C);
auto S3 = doubleSquare<ll>(A,B,p);
return S == S1 + S2 + S3;
}
inline ll orientedSquare(const std::vector<pii> &points) {
ll res = cross<ll>(points.back(), points.front());
for (int i = 0; i + 1 < isz(points); i++)
res += cross<ll>(points[i], points[i+1]);
return res;
}
inline vpii readPolygon() {
int n;
vpii points;
std::cin >> n;
points.resize(n);
for (auto &[x,y] : points)
std::cin >> x >> y;
return points;
}
inline void moveToPositive(vpii &points) {
int minX = INT_MAX, minY = INT_MAX;
for (auto &[x,y] : points) {
remin(minX, x);
remin(minY, y);
}
if (minX % 2 != 0) minX--;
if (minY % 2 != 0) minY--;
for (auto &[x,y] : points)
x -= minX, y -= minY;
}
inline void transpose(vpii &points) {
for (auto &[x,y] : points)
std::swap(x, y);
}
inline void makeCounterClockwise(vpii &points) {
ll sq = orientedSquare(points);
if (sq < 0) std::reverse(all(points));
}
inline int sign(auto x) { return x < 0 ? -1 : (x > 0 ? +1 : 0); }
const int LEFT = -1, RIGHT = +1;
inline int rotation(pii first, pii mid, pii last) {
// > 0 --> right
// < 0 --> left
return sign(cross<ll>(first-mid,last-mid));
}
inline bool isPointOnSegment(pii P, pii A, pii B) {
auto [x, y] = P;
auto [xmin,xmax] = std::minmax({A.first, B.first});
auto [ymin,ymax] = std::minmax({A.second, B.second});
if (x < xmin || x > xmax || y < ymin || y > ymax)
return false;
return cross<ll>(A-P,B-P) == 0;
}
const int NONE = 0, POINT = 1, INSIDE = 2;
inline int intersects(pii A, pii B, pii C, pii D) {
if (A == C || A == D || B == C || B == D)
return POINT;
if (isPointOnSegment(A,C,D) || isPointOnSegment(B,C,D) ||
isPointOnSegment(C,A,B) || isPointOnSegment(D,A,B))
return POINT;
auto [xmin1,xmax1] = std::minmax({A.first, B.first});
auto [ymin1,ymax1] = std::minmax({A.second, B.second});
auto [xmin2,xmax2] = std::minmax({C.first, D.first});
auto [ymin2,ymax2] = std::minmax({C.second, D.second});
if (xmax1 < xmin2 || xmin1 > xmax2)
return NONE;
if (ymax1 < ymin2 || ymin1 > ymax2)
return NONE;
auto rot1 = rotation(A,B,C);
auto rot2 = rotation(A,B,D);
if (sign(rot1) * sign(rot2) < 0)
return INSIDE;
return NONE;
}
inline bool isAlmostEqual(ld x, ld y) {
const ld EPS = 1e-10;
return std::abs(x-y) / (1 + std::abs(x)) < EPS;
}
inline bool angleLessThanPi(pii L, pii M, pii R) {
return rotation(L,M,R) == LEFT;
}
inline bool angleGreaterThanPi(pii L, pii M, pii R) {
return rotation(L,M,R) == RIGHT;
}
inline auto getPrev(const auto &vec, int i) {
i--;
if (i < 0) i = isz(vec)-1;
return vec[i];
}
inline auto getNext(const auto &vec, int i) {
i++;
if (i == isz(vec))
i = 0;
return vec[i];
}
inline bool isConvex(const vpii &P) {
for (int i = 0; i < isz(P); i++)
if (rotation(getPrev(P, i), P[i], getNext(P, i)) == RIGHT)
return false;
return true;
}
inline void removeSameLine(vpii &P) {
static std::vector<bool> keep;
keep.assign(isz(P), false);
for (int i = 0; i < isz(P); i++)
keep[i] = (rotation(getPrev(P,i), P[i], getNext(P,i)) != 0);
int p = 0;
for (int i = 0; i < isz(P); i++)
if (keep[i]) P[p++] = P[i];
P.resize(p);
}
} // namespace geoma
} // namespace algos
#endif // __GEOMA_HPP__
#ifndef __TRIANGULATE_HPP__
#define __TRIANGULATE_HPP__
namespace algos {
namespace triangulate {
using namespace geoma;
inline bool liesBelow(pii P1, pii P2) {
auto [x1,y1] = P1;
auto [x2,y2] = P2;
return (y1 < y2 || (y1 == y2 && x1 > x2));
}
inline bool liesAbove(pii P1, pii P2) {
auto [x1,y1] = P1;
auto [x2,y2] = P2;
return (y1 > y2 || (y1 == y2 && x1 < x2));
}
const int START = 0, SPLIT = 1, END = 2, MERGE = 3, REGULAR = 4;
inline auto getTypeOfVertex(const auto &P, int i) {
pii L = getPrev(P, i), M = P[i], R = getNext(P, i);
if (liesBelow(L,M) && liesBelow(R,M) && angleLessThanPi(L,M,R))
return START;
if (liesBelow(L,M) && liesBelow(R,M) && angleGreaterThanPi(L,M,R))
return SPLIT;
if (liesAbove(L,M) && liesAbove(R,M) && angleLessThanPi(L,M,R))
return END;
if (liesAbove(L,M) && liesAbove(R,M) && angleGreaterThanPi(L,M,R))
return MERGE;
return REGULAR;
}
inline std::string type2str(int t) {
switch(t) {
case START: return "start";
case SPLIT: return "split";
case END: return "end";
case MERGE: return "merge";
case REGULAR: return "regular";
default: return "none";
};
}
inline bool isMonotone(const vpii &P) {
if (isz(P) <= 3)
return true;
// we should start from max vertex
int imax = -1, ymax = INT_MIN;
for (int i = 0; i < isz(P); i++)
if (P[i].second >= ymax) {
imax = i;
ymax = P[i].second;
}
int len = 0;
int i = imax;
while (len <= isz(P) && P[i].second >= getNext(P, i).second) {
len++;
i++;
if (i == isz(P)) i = 0;
}
while (len <= isz(P) && P[i].second <= getNext(P, i).second) {
len++;
i++;
if (i == isz(P)) i = 0;
}
return len >= isz(P);
}
inline vpii findDiagonals(const vpii &P) {
static vtiii pq;
pq.clear();
for (int i = 0; i < isz(P); i++) {
auto [x,y] = P[i];
pq << tiii(y,-x,i);
}
std::sort(all(pq), std::greater<>());
int currY, currX;
auto setCurrY = [&](int newCurrY) { currY = newCurrY; };
auto setCurrX = [&](int newCurrX) { currX = newCurrX; };
static vi helper;
helper.assign(isz(P),-1);
auto getHelper = [&](int i) {
if (i < 0) i += isz(P);
assert(0 <= i && i < isz(P));
assert(helper[i] != -1);
return helper[i];
};
auto setHelper = [&](int i, int j) {
if (i < 0) i += isz(P);
assert(0 <= i && i < isz(P));
if (helper[i] == i)
helper[i] = j;
else if (helper[i] != -1) {
int k = helper[i];
if (P[k].second == P[j].second) {
//helper[i];
} else if (P[k].second > P[j].second) {
helper[i] = j;
}
} else
helper[i] = j;
};
auto getEdge = [&](int i) {
if (i < 0) i += isz(P);
if (i >= isz(P)) i -= isz(P);
auto L = P[i], R = getNext(P, i);
auto [ymin, ymax] = std::minmax({L.second, R.second});
return tiii(ymin, ymax, i);
};
auto getEdgeId = [&](int i) {
auto [y1,y2,j] = getEdge(i);
if (y1 == y2) return -1;
return j;
};
int cachedId = -1;
ld cachedX = -1;
auto isCached = [&](int id) { return id == cachedId; };
auto setCached = [&](int id, ld x) { cachedId = id; cachedX = x; };
auto clearCached = [&](){cachedId = -1; cachedX = -1;};
auto coordOfIntersection = [&](int i, int y) -> ld {
if (i < 0) i += isz(P);
if (isCached(i))
return cachedX;
assert(0 <= i && i < isz(P));
auto [x1,y1] = P[i];
auto [x2,y2] = getNext(P, i);
//ld B = ld(x2-x1)/(y2-y1);
//ld A = x1 - y1 * B;
// x1 - y1 * (x2-x1)/(y2-y1)
// x1 - y1 * (x2-x1)/(y2-y1) + y * (x2-x1)/(y2-y1)
// (x1 * (y2 - y1) + (y - y1) * (x2 - x1)) / (y2 - y1)
// (x1 * y2 + y * x2 - y * x1 - y1 * x2) / (y2 - y1)
assert(y2 != y1);
return ld((ld)x1 * y2 + (ld)y * x2 - (ld)y * x1 - (ld)y1 * x2) / (y2 - y1);
};
auto LessByIntersection = [&](const int &j1, const int &j2) {
ld x1 = coordOfIntersection(j1, currY);
ld x2 = coordOfIntersection(j2, currY);
if (x1 > x2) return false;
if (x1 < x2) return true;
return j1 < j2;
};
// y1 + t * (y2-y1) = y => t = (y - y1)/(y2 - y1)
// x1 + t * (x2-x1) = x => x1 + (y - y1) * (x2-x1) / (y2 - y1)
// coeff = (x2 - x1) / (y2 - y1)
// x1 + (y - y1) * coeff <= x
// x1 - y1 * coeff + y * coeff <= x
// A + y * B <= x -- need to sort by these values
std::set<int, decltype(LessByIntersection)> T(LessByIntersection);
auto insertEdgeInT = [&](int i) {
int ei = getEdgeId(i);
if (ei == -1) return;
auto y = P[i].second;
setCurrY(y);
setCached(ei, coordOfIntersection(ei, currY));
T.insert(ei);
clearCached();
};
auto deleteEdgeFromT = [&](int i) {
int j = getEdgeId(i);
if (j == -1) return;
auto it = T.find(j);
assert(it != T.end());
T.erase(it);
};
auto searchEdgeInT = [&](int i) {
assert(0 <= i && i < isz(P));
const auto [x, y] = P[i];
setCurrY(y); setCurrX(x);
int cand;
setCached(isz(P), x);
auto it = T.upper_bound(isz(P));
assert(it != T.begin());
cand = *std::prev(it);
clearCached();
return cand;
};
vpii diagonals;
auto addEdgeInTriangles = [&](int u, int v) {
if (getNext(P, u) != P[v] && getPrev(P, u) != P[v])
diagonals.emplace_back(u,v);
};
for (auto [y,x,i] : pq) {
x *= -1;
setCurrY(y);
switch (getTypeOfVertex(P,i)) {
case START: {
insertEdgeInT(i);
setHelper(i, i);
break;
}
case END: {
addEdgeInTriangles(i, getHelper(i-1));
deleteEdgeFromT(i-1);
break;
}
case SPLIT: {
int j = searchEdgeInT(i);
addEdgeInTriangles(i, getHelper(j));
setHelper(j, i);
insertEdgeInT(i);
setHelper(i, i);
break;
}
case MERGE: {
addEdgeInTriangles(i, getHelper(i-1));
deleteEdgeFromT(i-1);
int j = searchEdgeInT(i);
addEdgeInTriangles(i, getHelper(j));
setHelper(j, i);
break;
}
case REGULAR: {
if (liesAbove(getPrev(P,i), P[i]) && liesBelow(getNext(P,i), P[i])) {
addEdgeInTriangles(i, getHelper(i-1));
deleteEdgeFromT(i-1);
insertEdgeInT(i);
setHelper(i, i);
} else {
int j = searchEdgeInT(i);
addEdgeInTriangles(i, getHelper(j));
setHelper(j, i);
}
break;
}
};
}
return diagonals;
}
struct Node {
const vpii * parent = nullptr;
Node * next = nullptr;
Node * prev = nullptr;
int where{}, id{};
std::list<int> diagonals{};
void setNext(Node *n) {
next = n;
n->prev = this;
}
void setPrev(Node *n) {
prev = n;
n->next = this;
}
vpii to_vector() const {
vpii result(1, (*parent)[id]);
for (Node *it = next; it != this; it = it->next)
result.push_back((*parent)[it->id]);
return result;
}
pii to_point() const {
return (*parent)[id];
}
};
inline std::ostream &operator<<(std::ostream &os, const Node * node) {
return os << "(node #" << node->id+1 << " in list " << node->where << ")";
}
template<typename Func>
inline void splitByDiag(const vpii &P, const vpii &D, const Func &f) {
if (isz(P) <= 2)
return;
if (isz(P) == 3) {
assert(D.empty());
f(P);
return;
}
static std::vector<Node *> nodes;
nodes.clear();
for (int i = 0; i < isz(P); i++)
nodes.emplace_back(new Node{&P, nullptr, nullptr, 0, i, {}});
for (int i = 0; i+1 < isz(nodes); i++)
nodes[i]->setNext(nodes[i+1]);
nodes.front()->setPrev(nodes.back());
static std::vector<std::pair<Node*, Node*>> diagonals;
diagonals.clear();
for (const auto &[u,v] : D) {
assert(0 <= u && u < isz(nodes));
assert(0 <= v && v < isz(nodes));
int currDiag = isz(diagonals);
diagonals.push_back({nodes[u], nodes[v]});
nodes[u]->diagonals.push_back(currDiag);
nodes[v]->diagonals.push_back(currDiag);
}
int cntPoly = 1;
for (int id = 0; id < isz(diagonals); id++) {
// need to split by this diagonal
auto [startFi, startSe] = diagonals[id];
auto currFi = startFi, currSe = startSe;
while (currFi != startSe && currSe != startFi) {
currFi = currFi->next;
currSe = currSe->next;
}
// need to copy ends of diagonal:
nodes.push_back(new Node{&P, nullptr, nullptr, cntPoly, startFi->id, {}});
Node *copyFi = nodes.back();
nodes.push_back(new Node{&P, nullptr, nullptr, cntPoly, startSe->id, {}});
Node *copySe = nodes.back();
// choose a list with min len:
if (currFi != startSe) {
std::swap(copyFi, copySe);
std::swap(startFi, startSe);
std::swap(currFi, currSe);
}
assert(currFi == startSe);
// iterate over all vertices and updates diagonals
for (currFi = currFi->prev; currFi != startFi; currFi = currFi->prev) {
auto &dd = currFi->diagonals;
for (auto it = dd.begin(); it != dd.end(); ) {
// this diagonal still exists?
if (*it <= id) {
it = dd.erase(it);
continue;
}
auto &[lhs, rhs] = diagonals[*it];
if (lhs == startFi) (lhs = copyFi)->diagonals.push_back(*it);
if (rhs == startFi) (rhs = copyFi)->diagonals.push_back(*it);
if (lhs == startSe) (lhs = copySe)->diagonals.push_back(*it);
if (rhs == startSe) (rhs = copySe)->diagonals.push_back(*it);
it = std::next(it);
}
// set where this vertex lies
currFi->where = cntPoly;
}
// disconnect less half:
copyFi->setNext(currFi->next);
copySe->next = copyFi;
copySe->setNext(copyFi);
startSe->prev->setNext(copySe);
startFi->setNext(startSe);
cntPoly++;
}
static std::vector<Node *> anyNode;
anyNode.resize(cntPoly);
for (auto it : nodes)
anyNode[it->where] = it;
for (int i = 0; i < cntPoly; i++) {
vpii temp = anyNode[i]->to_vector();
removeSameLine(temp);
//assert(isMonotone(temp));
f(temp);
}
for (auto node : nodes)
delete node;
}
inline vvpii splitByDiag(const vpii &P, const vpii &D) {
vvpii result;
splitByDiag(P, D, [&](const vpii &PP){ result.push_back(PP); });
return result;
}
inline vpii TriangulateMonotonePolygon(const vpii &P) {
if (isz(P) <= 3)
return {};
// we need to create chains
int ymax = INT_MIN, ymin = INT_MAX, imax = -1;
int leftChainBegin = -1, leftChainEnd = -1;
for (int i = 0; i < isz(P); i++) {
if (P[i].second >= ymax) {
imax = i;
ymax = P[i].second;
}
if (P[i].second == ymax && getNext(P, i).second != ymax) {
leftChainBegin = i;
}
remin(ymin, P[i].second);
}
const int LEFT_CHAIN = 0, RIGHT_CHAIN = 1;
assert(imax != -1);
static vi chainType; chainType.assign(isz(P), -1);
static vi chainPos; chainPos.assign(isz(P), -1);
assert(leftChainBegin != -1);
leftChainEnd = leftChainBegin;
int pos = 0;
while (!(P[leftChainEnd].second == ymin && getPrev(P, leftChainEnd).second != ymin)) {
chainType[leftChainEnd] = LEFT_CHAIN;
chainPos[leftChainEnd] = pos++;
leftChainEnd++;
if (leftChainEnd == isz(P))
leftChainEnd = 0;
}
chainType[leftChainEnd] = LEFT_CHAIN;
if (chainPos[leftChainEnd] == -1)
chainPos[leftChainEnd] = pos++;
pos = 0;
for (int i = (leftChainBegin-1+isz(P)) % isz(P); chainType[i] == -1; ) {
chainType[i] = RIGHT_CHAIN;
chainPos[i] = pos++;
i--;
if (i < 0) i = isz(P)-1;
}
// {y-coord, chainId, posInChain}
auto cmp = [&](pii a, pii b) {
auto [ya, ida] = a;
auto [yb, idb] = b;
if (ya < yb) return true;
if (ya > yb) return false;
if (chainType[ida] != chainType[idb]) {
return chainType[ida] < chainType[idb];
}
return chainPos[ida] > chainPos[idb];
};
static vpii pq;
pq.clear();
for (int i = 0; i < isz(P); i++) {
auto [x,y] = P[i];
pq.push_back({y,i});
}
std::sort(pq.rbegin(), pq.rend(), cmp);
static vi S;
S.clear();
S.push_back(std::get<1>(pq[0]));
S.push_back(std::get<1>(pq[1]));
vpii diagonals;
auto addDiagonal = [&](int i, int j) {
if (getNext(P, i) != P[j] && getPrev(P, i) != P[j]) {
diagonals.emplace_back(i, j);
return true;
}
return false;
};
auto isSameChain = [&](int i, int j) {
return chainType[i] == chainType[j];
};
for (const auto &[y, id] : pq) {
if (!isSameChain(id,S.back())) {
int last = S.back();
while (isz(S) > 1) {
addDiagonal(id, S.back());
S.pop_back();
}
S = {last, id};
} else {
int want = LEFT, skip = RIGHT;
if (chainType[id] == LEFT_CHAIN)
std::swap(want, skip);
pii first = P[id];
while (isz(S) >= 2) {
pii mid = P[S[isz(S)-1]];
pii last = P[S[isz(S)-2]];
ll rot = rotation(first, mid, last);
if (rot == skip)
break;
if (rot != want) {
S.pop_back();
continue;
}
assert(rot == want);
addDiagonal(id, S[isz(S)-2]);
S.pop_back();
}
S.push_back(id);
}
}
return diagonals;
}
template<typename Func>
inline void triangulateFast(vpii inP, const Func &f) {
ll fullSquare = orientedSquare(inP);
removeSameLine(inP);
vpii D = findDiagonals(inP);
splitByDiag(inP, D, [&](vpii PP){
removeSameLine(PP);
vpii DD = TriangulateMonotonePolygon(PP);
auto ff = [&](const vpii &t){
fullSquare -= orientedSquare(t);
f(t);
};
splitByDiag(PP, DD, ff);
});
assert(fullSquare == 0);
}
inline std::vector<Triangle> triangulateFast(const vpii &inP) {
std::vector<Triangle> answ;
triangulateFast(inP, [&](const vpii &t){
assert(isz(t) <= 3);
if (isz(t) == 3) answ.emplace_back(t[0], t[1], t[2]);
});
return answ;
}
inline std::vector<Triangle> triangulateSlow(const vpii &points)
{
// slow triangulation in time O(n^2)
std::vector<Triangle> triangles;
std::list<pii> list(all(points));
while (isz(list) > 3) {
// search a vertex which we can cut off:
for (auto it = list.begin(); it != list.end(); ) {
auto next = (std::next(it) == list.end() ? list.begin() : std::next(it));
auto prev = (it == list.begin() ? std::prev(list.end()) : std::prev(it));
ll rot = rotation(*prev, *it, *next);
if (rot == 0) {
it = list.erase(it);
} else if (rot < 0) {
// check this vertex:
bool ok = true;
for (auto p : list) {
bool bad = isPointInsideTriangle({*prev, *it, *next}, p);
bad = bad && !isPointOnSegment(p, *prev, *it);
bad = bad && !isPointOnSegment(p, *it, *next);
bad = bad && !isPointOnSegment(p, *prev, *next);
if (bad) {
ok = false;
break;
}
}
if (ok) {
triangles.emplace_back(*prev, *it, *next);
it = list.erase(it);
} else {
it++;
}
} else {
it++;
}
}
}
if (isz(list) == 3)
triangles.emplace_back(*list.begin(), *(++list.begin()), *(++++list.begin()));
else assert(isz(list) < 3);
return triangles;
}
inline std::vector<Triangle> triangulateConvex(const vpii &points) {
// points is convex polygon
// take any vertex and make triangles
std::vector<Triangle> triangles;
for (int i = 1; i + 1 < isz(points); i++)
triangles.emplace_back(points[0], points[i], points[i+1]);
return triangles;
}
inline std::vector<Triangle> triangulate(vpii points) {
makeCounterClockwise(points);
return triangulateFast(points);
}
} // namespace triangulate
} // namespace algos
#endif // __TRIANGULATE_HPP__
using namespace algos::triangulate;
using namespace algos::geoma;
pii pickPointOutside(pii A, pii B) {
pii C1(A.first, B.second);
pii C2(B.first, A.second);
pii res = cross<ll>(A-C1, B-C1) < 0 ? C1 : C2;
watch("pickPointOutside", A, B, res);
return res;
}
std::pair<ld,ld> rectangle(int xmin, int ymin, int xmax, int ymax) {
ll dx = xmax - xmin, dy = ymax - ymin;
ld fi = (dx + 1) / 2 * dy;
ld se = dx / 2 * dy;
if (xmin % 2 != 0) std::swap(fi,se);
watch("rectangle", xmin, ymin, xmax, ymax, fi, se);
return {fi, se};
}
ll sum(ll n) {
// n^2 - (n-1)^2 + (n-2)^2 - (n-3)^2 + (n-4)^2 + ...
// wolfram: n*(n+1)/2
// check: 4^2-3^2+2^1-1=16-9+4-1=10
// check: 4 * 5 / 2 = 10
return n * (n+1) / 2;
}
std::pair<ld, ld> solveBaseCaseLD(ll x, ll y) {
// triangle
// *
// **
// ***
// ****
watch("solveBaseCaseLD",x,y);
ld c = y / (2.0L * x);
std::pair<ld, ld> res{c * sum(x), c * sum(x-1)};
watch(res);
return res;
}
std::pair<ld, ld> solveBaseCaseRU(ll x, ll y) {
// triangle
// ****
// ***
// **
// *
auto res = solveBaseCaseLD(x, y);
if (x % 2 == 0)
std::swap(res.first, res.second); // проверено
return res;
}
std::pair<ld, ld> solveBaseCaseRU(ll xmin, ll ymin, ll xmax, ll ymax) {
// triangle
// ****
// ***
// **
// *
auto res = solveBaseCaseRU(xmax-xmin, ymax-ymin);
if (xmin % 2 != 0) std::swap(res.first, res.second);
return res;
}
std::pair<ld, ld> solveBaseCaseLD(ll xmin, ll ymin, ll xmax, ll ymax) {
// triangle
// *
// **
// ***
// ****
auto res = solveBaseCaseLD(xmax-xmin, ymax-ymin);
if (xmin % 2 != 0) std::swap(res.first, res.second);
return res;
}
std::pair<ld, ld> solveBaseCaseRD(ll x, ll y) {
// triangle
// *
// **
// ***
// ****
auto res = solveBaseCaseLD(x, y);
if (x % 2 == 0) // проверено
std::swap(res.first, res.second);
return res;
}
std::pair<ld, ld> solveBaseCaseLU(ll x, ll y) {
// triangle
// ****
// ***
// **
// *
auto res = solveBaseCaseLD(x, y);
return res;
}
std::pair<ld, ld> solveBaseCaseLU(ll xmin, ll ymin, ll xmax, ll ymax) {
// triangle
// ****
// ***
// **
// *
auto res = solveBaseCaseLU(xmax-xmin, ymax-ymin);
if (xmin % 2 != 0) std::swap(res.first, res.second);
return res;
}
std::pair<ld, ld> solveBaseCaseRD(ll xmin, ll ymin, ll xmax, ll ymax) {
// triangle
// *
// **
// ***
// ****
auto res = solveBaseCaseRD(xmax-xmin,ymax-ymin);
if (xmin % 2 != 0)
std::swap(res.first, res.second);
return res;
}
std::pair<ld, ld> triangle(pii A, pii B, pii C) {
// B - corner
// AC - diagonal
if (doubleSquare<ll>(A,B,C) == 0)
return {0, 0};
if (A.first > C.first)
std::swap(A, C);
ll xmin = std::min({A.first, B.first, C.first});
ll xmax = std::max({A.first, B.first, C.first});
ll ymin = std::min({A.second, B.second, C.second});
ll ymax = std::max({A.second, B.second, C.second});
pii LD(xmin,ymin), LU(xmin,ymax), RD(xmax,ymin), RU(xmax,ymax);
std::pair<ld,ld> res;
if (A.second > C.second) {
if (B == LD) res = solveBaseCaseLD(xmin,ymin,xmax,ymax);
else res = solveBaseCaseRU(xmin,ymin,xmax,ymax);
} else {
if (B == RD) res = solveBaseCaseRD(xmin,ymin,xmax,ymax);
else res = solveBaseCaseLU(xmin,ymin,xmax,ymax);
}
return res;
}
bool inter(pii A, pii B, int x, ld &y) {
if (A.first > B.first) std::swap(A, B);
if (x < A.first || x > B.first)
return false;
if (x == A.first) {
y = A.second;
return true;
}
if (x == B.first) {
y = B.second;
return true;
}
pii v = B - A;
// A + t * (B - A) == {x, y}
// A.first + t * v.first == x
ld t = (x - A.first) / ld(v.first);
y = A.second + t * v.second;
return true;
}
ld solveVerticalCase(pii A, pii B, pii C, int x) {
std::vector<ld> y1, y2;
ld y;
if (inter(A,B,x,y)) y1 << y;
if (inter(B,C,x,y)) y1 << y;
if (inter(A,C,x,y)) y1 << y;
if (inter(A,B,x+1,y)) y2 << y;
if (inter(B,C,x+1,y)) y2 << y;
if (inter(A,C,x+1,y)) y2 << y;
std::sort(all(y1));
std::sort(all(y2));
ld res = (y1.back() - y1.front() + y2.back() - y2.front()) / 2;
return res;
}
std::pair<ld,ld> twoRectangles(pii A, pii B, pii C) {
if (A.first > B.first) std::swap(A, B);
if (A.first > C.first) std::swap(A, C);
if (B.first > C.first) std::swap(B, C);
assert(A.first <= B.first && B.first <= C.first);
auto [ymin, ymax] = std::minmax({A.second,C.second});
if (B.second >= ymax || B.second <= ymin || B.first == A.first || B.first == C.first) {
remin(ymin, B.second);
remax(ymax, B.second);
return rectangle(A.first, ymin, C.first, ymax);
}
std::pair<ld,ld> rect = rectangle(A.first, ymin, C.first, ymax);
// above or below than diagonal?
if (rotation(A,C,B) <= 0) { // above
// is A below than C?
if (A.second <= C.second)
rect = rect - rectangle(A.first, B.second, B.first, ymax);
else
rect = rect - rectangle(B.first, B.second, C.first, ymax);
} else { // below
// is A below than C?
if (A.second <= C.second)
rect = rect - rectangle(B.first, ymin, C.first, B.second);
else
rect = rect - rectangle(A.first, ymin, B.first, B.second);
}
return rect;
}
std::pair<ld,ld> calcTriangle(pii A, pii B, pii C) {
if (cross<ll>(A-B,C-B) > 0)
std::swap(A, C);
auto temp = twoRectangles(A,B,C);
auto P1 = pickPointOutside(A,B);
auto P2 = pickPointOutside(B,C);
auto P3 = pickPointOutside(C,A);
if (!(A.first == B.first || A.second == B.second))
temp = temp - triangle(A,P1,B);
if (!(B.first == C.first || B.second == C.second))
temp = temp - triangle(B,P2,C);
if (!(C.first == A.first || C.second == A.second))
temp = temp - triangle(C,P3,A);
return temp;
}
ld calcTriangleSlow(pii A, pii B, pii C) {
auto [xMin, xMax] = std::minmax({A.first, B.first, C.first});
ld ans{};
int x = xMin;
if (x % 2 != 0)
x++;
while (x+1 <= xMax) {
assert(x % 2 == 0);
ans += solveVerticalCase(A,B,C,x);
x += 2;
}
return ans;
}
void stress(int xs, int xf) {
//const int X = 10;
for (int x1 = xs; x1 <= xf; x1++)
for (int y1 = xs; y1 <= xf; y1++)
for (int x2 = xs; x2 <= xf; x2++)
for (int y2 = xs; y2 <= xf; y2++)
for (int x3 = xs; x3 <= xf; x3++)
for (int y3 = xs; y3 <= xf; y3++) {
pii A(x1,y1), B(x2,y2), C(x3,y3);
ll s = doubleSquare<ll>(A,B,C);
if (s == 0) continue;
ld fast = calcTriangle(A,B,C).first;
ld slow = calcTriangleSlow(A,B,C);
if (!isAlmostEqual(fast,slow)) {
debug = 1;
fast = calcTriangle(A,B,C).first;
watch(A,B,C);
watch(fast);
watch(slow);
std::exit(0);
}
}
}
std::mt19937 gen;
int randInt(int a, int b) {
assert(a <= b);
return std::uniform_int_distribution<int>(a,b)(gen);
}
vpii generatePoly(int xmin, int xmax, int ymin, int ymax, int n, int limit = 1e6) {
//bool debug = 1;
std::set<pii> used;
auto randPoint = [&](){
while (true) {
int x = randInt(xmin, xmax);
int y = randInt(ymin, ymax);
if (!used.count(pii(x,y)))
return pii(x, y);
}
};
vpii kek;
auto pushPoint = [&](pii p) {
kek << p;
used << p;
};
pushPoint(randPoint());
pushPoint(randPoint());
pushPoint(randPoint());
while (isz(kek) < n) {
int attempts = 0;
while (attempts < limit) {
attempts++;
auto p = randPoint();
int nP{}, nI{};
for (int i = 0; i + 1 < isz(kek); i++) {
int t = intersects(kek.back(), p, kek[i], kek[i+1]);
nP += (t == POINT);
nI += (t == INSIDE);
}
if (!(nI == 0 && nP == 1))
continue;
nP = nI = 0;
for (int i = 0; i + 1 < isz(kek); i++) {
int t = intersects(kek.front(), p, kek[i], kek[i+1]);
nP += (t == POINT);
nI += (t == INSIDE);
}
if (nI == 0 && nP == 1) {
pushPoint(p);
break;
}
}
if (attempts >= 1000)
break;
}
if (orientedSquare(kek) < 0)
std::reverse(all(kek));
return kek;
}
void testMonotone() {
//bool debug = 1;
//vpii P{{14,5},{12,6},{11,11},{8,10},{7,11},{3,9},{6,8},{5,5},{2,7},{1,3},{4,1},{7,2},{10,1},{9,4},{13,3}};
//vpii P{{-3, 5}, {-1, 5}, {4, 5}, {-4, 8}, {-7, 9}, {-10, 5}, {6, -3}, {5, -2}, {-4, 4
//}, {9, 4}};
//vpii P{{-7, 8}, {-4, 2}, {-9, -9}, {2, 9}, {1, 10}, {-1, 10}, {-5, 9}, {-4, 10}, {-6,
//10}, {-9, 5}};
//vpii P{{-10, 6}, {-9, 1}, {-6, -9}, {-6, -3}, {-4, -5}, {-3, -4}, {9, -3}, {-2, 4}, {
//-1, 0}, {-7, -1}};
//vpii P{{7, -6}, {6, -6}, {6, -4}, {3, 9}, {4, -7}, {1, -10}, {8, -9}, {8, 0}, {7, -1}
//, {6, 7}};
//vpii P{{3, 5}, {-2, 9}, {3, -5}, {3, -4}, {4, -7}, {9, -6}, {9, -5}, {8, -2}, {2, 7},
//{7, -1}};
vpii P{{9, 1}, {6, 1}, {10, 2}, {6, 2}, {10, 4}, {7, 5}, {5, 5}, {4, 10}, {-2, -1}, {
10, 1}};
//vpii P{{5, -6}, {10, 9}, {7, 0}};
std::cout << triangulateFast(P).size() << std::endl;
for (int nTests = 0; nTests <= 10000; nTests++) {
watch("start generate P");
while (true) {
//P = generatePoly(-10,10,-10,10,10);
//P = generatePoly(-100,100,-100,100,100);
P = generatePoly(-1e9,1e9,-1e9,1e9,2000);
ll s = orientedSquare(P);
if (s != 0)
break;
}
watch("P: generated");
//print("poly2.txt", {P},P,{});
//watch(P);
//break;
std::cout << triangulateFast(P).size() << " ";
//print("poly.txt", PP, P, D);
}
std::cout << "OK" << std::endl;
std::exit(0);
}
vpii genWorst1() {
const int XMIN = -(int)1e9;
const int XMAX = +(int)1e9;
const int YMAX = +(int)1e9;
vpii left;
for (int i = 0; i < 25000; i++) {
left.emplace_back(XMIN, YMAX-2*i);
left.emplace_back(XMAX-1, YMAX-(2*i+1));
}
vpii right = left;
for (auto &[x,y] : right) x++;
std::reverse(all(right));
return left << right;
}
vpii genWorst2() {
const int XMAX = +(int)1e9;
const int YMIN = -(int)1e9;
const int YMAX = +(int)1e9;
vpii top;
for (int i = 0; i < 25000; i++) {
top.emplace_back(XMAX-2*i, YMIN);
top.emplace_back(XMAX-(2*i+1), YMAX-1);
}
vpii bottom = top;
for (auto &[x,y] : bottom) y++;
std::reverse(all(bottom));
return top << bottom;
}
ld accumulateTriangles(const std::vector<Triangle> &triangles) {
ld ans{};
for (auto &[A,B,C] : triangles)
ans += calcTriangle(A,B,C).first;
return ans;
}
ld solveAllSubtasks(vpii points) {
moveToPositive(points);
transpose(points);
makeCounterClockwise(points);
std::vector<Triangle> triangles = isz(points) <= 1000 ? triangulateSlow(points)
: isConvex(points) ? triangulateConvex(points)
: triangulateFast(points);
return accumulateTriangles(triangles);
}
ld solveFast(vpii points) {
moveToPositive(points);
transpose(points);
makeCounterClockwise(points);
ld answ = 0;
triangulateFast(points, [&](const vpii &tr){
assert(isz(tr) <= 3);
if (isz(tr) == 3)
answ += calcTriangle(tr[0],tr[1],tr[2]).first;
});
return answ;
}
#ifndef __TIMER_HPP__
#define __TIMER_HPP__
class Timer {
std::chrono::time_point<std::chrono::steady_clock> timePoint;
size_t value;
public:
void start() { timePoint = std::chrono::steady_clock::now(); }
void finish() {
auto curr = std::chrono::steady_clock::now();
auto elapsed = std::chrono::duration_cast<std::chrono::milliseconds>(curr - timePoint);
value = elapsed.count();
}
size_t operator()() const { return value; }
};
#endif // __TIMER_HPP__
void testWorst1() {
std::cout << __FUNCTION__ << ": ";
Timer timer;
timer.start();
vpii P = genWorst1();
ld answ = solveFast(P);
timer.finish();
std::cout << "answ = " << answ << ", runtime = " << timer() << " ms" << std::endl;
}
void testWorst2() {
std::cout << __FUNCTION__ << ": ";
Timer timer;
timer.start();
vpii P = genWorst2();
ld answ = solveFast(P);
timer.finish();
std::cout << "answ = " << answ << ", runtime = " << timer() << " ms" << std::endl;
}
void testAll() {
testWorst1();
testWorst2();
std::exit(0);
}
int main() {
//testMonotone();
//testAll();
//testMonotone();
auto points = readPolygon();
std::cout << solveFast(points) << std::endl;
}
This solution works in $$$296$$$ ms in Yandex.Contest when $$$n=10^5$$$.
