Geometry: Polygon
(geometry/polygon.hpp)

View this file on GitHub
Last update: 2026-10-22 23:24:07+09:00
Include: #include "geometry/polygon.hpp"

概要

多角形と点集合のアルゴリズム。

使い方

geo::polygon<T>: std::vector<geo::point<T>>。頂点を順に並べる。
area2(poly): 符号付き面積の 2 倍を返す。反時計回りで正。
area(poly): 符号付き面積を返す。浮動小数点型のみ。
is_convex(poly): 反時計回りの凸多角形かを返す。$180^\circ$ の角を許す。
contains(poly, p): 点 p が多角形の内部なら $2$、辺上なら $1$、外部なら $0$ を返す。
convex_hull(ps, strict = true): 点集合 ps の凸包を反時計回りで返す。strict なら辺上の点を含まない。
farthest_pair(ps): 最遠点対の添字を返す。
closest_pair(ps): 最近点対の添字を返す。
convex_cut(poly, l): 凸多角形 poly を直線 l で切った左側を返す。浮動小数点型のみ。

制約

farthest_pair, closest_pair では $\lvert \mathrm{ps} \rvert \ge 2$
convex_cut では poly は反時計回りの凸多角形

計算量

area2, area, is_convex, contains, convex_cut: $\mathcal{O}(N)$
convex_hull, farthest_pair, closest_pair: $\mathcal{O}(N \log N)$

問題例

Library Checker - Static Convex Hull（問題）

std::vector<geo::point<long long>> ps(n);
for (auto& p : ps) std::cin >> p;
geo::polygon<long long> ch = geo::convex_hull(ps);

Depends on

Required by

Geometry: Circle (geometry/circle.hpp)

Verified with

Code

#pragma once
#include <algorithm>
#include <cassert>
#include <cmath>
#include <concepts>
#include <utility>
#include <vector>

#include "geometry/line.hpp"
#include "geometry/point.hpp"
namespace geo {
template <class T> using polygon = std::vector<point<T>>;
template <class T> T area2(const polygon<T>& poly) {
    int n = int(poly.size());
    if constexpr (std::is_integral_v<T>) {
        __int128 ret = 0;
        for (int i = 0; i < n; i++) ret += (__int128)cross(poly[i], poly[(i + 1) % n]);
        return T(ret);
    } else {
        T ret = 0;
        for (int i = 0; i < n; i++) ret += cross(poly[i], poly[(i + 1) % n]);
        return ret;
    }
}
template <std::floating_point T> T area(const polygon<T>& poly) {
    return area2(poly) / T(2);
}
template <class T> bool is_convex(const polygon<T>& poly) {
    int n = int(poly.size());
    for (int i = 0; i < n; i++) {
        if (ccw(poly[i], poly[(i + 1) % n], poly[(i + 2) % n]) < 0) return false;
    }
    return true;
}
template <class T> int contains(const polygon<T>& poly, const point<T>& p) {
    int n = int(poly.size());
    bool in = false;
    for (int i = 0; i < n; i++) {
        point<T> a = poly[i] - p, b = poly[(i + 1) % n] - p;
        if (ccw(a, b, point<T>()) == 0 && sgn(dot(a, b)) <= 0) return 1;
        if (sgn(a.y) > sgn(b.y)) std::swap(a, b);
        if (sgn(a.y) <= 0 && 0 < sgn(b.y) && sgn(cross(a, b)) > 0) in = !in;
    }
    return in ? 2 : 0;
}
template <class T> polygon<T> convex_hull(std::vector<point<T>> ps, bool strict = true) {
    int n = int(ps.size());
    std::sort(ps.begin(), ps.end());
    ps.erase(std::unique(ps.begin(), ps.end()), ps.end());
    n = int(ps.size());
    if (n <= 2) return ps;
    polygon<T> ret(2 * n);
    int k = 0;
    for (int i = 0; i < n; ret[k++] = ps[i++]) {
        while (k >= 2 && ccw(ret[k - 2], ret[k - 1], ps[i]) < (strict ? 1 : 0)) k--;
    }
    for (int i = n - 2, t = k + 1; i >= 0; ret[k++] = ps[i--]) {
        while (k >= t && ccw(ret[k - 2], ret[k - 1], ps[i]) < (strict ? 1 : 0)) k--;
    }
    ret.resize(k - 1);
    return ret;
}
template <class T> std::pair<int, int> farthest_pair(const std::vector<point<T>>& ps) {
    int n = int(ps.size());
    assert(n >= 2);
    polygon<T> ch = convex_hull(ps);
    int m = int(ch.size());
    auto index_of = [&](const point<T>& p) {
        return int(std::find(ps.begin(), ps.end(), p) - ps.begin());
    };
    if (m == 1) return {0, 1};
    if (m == 2) return {index_of(ch[0]), index_of(ch[1])};
    std::pair<int, int> ret = {0, 1};
    T best = norm(ch[0] - ch[1]);
    for (int i = 0, j = 1; i < m; i++) {
        while (norm(ch[i] - ch[(j + 1) % m]) > norm(ch[i] - ch[j])) j = (j + 1) % m;
        if (norm(ch[i] - ch[j]) > best) {
            best = norm(ch[i] - ch[j]);
            ret = {i, j};
        }
    }
    return {index_of(ch[ret.first]), index_of(ch[ret.second])};
}
template <class T> std::pair<int, int> closest_pair(const std::vector<point<T>>& ps) {
    int n = int(ps.size());
    assert(n >= 2);
    std::vector<int> id(n);
    for (int i = 0; i < n; i++) id[i] = i;
    std::sort(id.begin(), id.end(), [&](int i, int j) { return ps[i] < ps[j]; });
    std::pair<int, int> ret = {id[0], id[1]};
    T best = norm(ps[id[0]] - ps[id[1]]);
    std::vector<int> buf(n);
    auto rec = [&](auto rec, int l, int r) -> void {
        if (r - l <= 1) return;
        int m = (l + r) / 2;
        T x = ps[id[m]].x;
        rec(rec, l, m);
        rec(rec, m, r);
        std::inplace_merge(id.begin() + l, id.begin() + m, id.begin() + r,
                           [&](int i, int j) { return ps[i].y < ps[j].y; });
        int k = 0;
        for (int t = l; t < r; t++) {
            int i = id[t];
            T dx = ps[i].x - x;
            if (dx * dx >= best) continue;
            for (int u = k - 1; u >= 0; u--) {
                int j = buf[u];
                T dy = ps[i].y - ps[j].y;
                if (dy * dy >= best) break;
                T d = norm(ps[i] - ps[j]);
                if (d < best) {
                    best = d;
                    ret = {j, i};
                }
            }
            buf[k++] = i;
        }
    };
    rec(rec, 0, n);
    return ret;
}
template <std::floating_point T> polygon<T> convex_cut(const polygon<T>& poly, const line<T>& l) {
    int n = int(poly.size());
    polygon<T> ret;
    for (int i = 0; i < n; i++) {
        const point<T>& p = poly[i];
        const point<T>& q = poly[(i + 1) % n];
        int cp = ccw(l.a, l.b, p), cq = ccw(l.a, l.b, q);
        if (cp >= 0) ret.push_back(p);
        if (cp * cq < 0) ret.push_back(cross_point(l, line<T>(p, q)));
    }
    return ret;
}
}  // namespace geo

#line 2 "geometry/polygon.hpp"
#include <algorithm>
#include <cassert>
#include <cmath>
#include <concepts>
#include <utility>
#include <vector>

#line 6 "geometry/line.hpp"
#include <iostream>

#line 5 "geometry/point.hpp"
#include <type_traits>
namespace geo {
inline long double eps = 1e-9;
template <class T> int sgn(T a) {
    if constexpr (std::is_integral_v<T>) {
        return (a > 0) - (a < 0);
    } else {
        return a < -eps ? -1 : a > eps ? 1 : 0;
    }
}
template <class T> struct point {
    T x, y;
    point() : x(0), y(0) {
    }
    point(T x, T y) : x(x), y(y) {
    }
    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*=(T k) {
        x *= k, y *= k;
        return *this;
    }
    point& operator/=(T k) requires std::floating_point<T> {
        x /= k, y /= k;
        return *this;
    }
    point operator-() const {
        return point(-x, -y);
    }
    friend point operator+(point a, const point& b) {
        return a += b;
    }
    friend point operator-(point a, const point& b) {
        return a -= b;
    }
    friend point operator*(point a, T k) {
        return a *= k;
    }
    friend point operator*(T k, point a) {
        return a *= k;
    }
    friend point operator/(point a, T k) requires std::floating_point<T> {
        return a /= k;
    }
    friend bool operator==(const point& a, const point& b) {
        return sgn(a.x - b.x) == 0 && sgn(a.y - b.y) == 0;
    }
    friend bool operator!=(const point& a, const point& b) {
        return !(a == b);
    }
    friend bool operator<(const point& a, const point& b) {
        if (sgn(a.x - b.x) != 0) return a.x < b.x;
        return sgn(a.y - b.y) < 0;
    }
    friend bool operator>(const point& a, const point& b) {
        return b < a;
    }
    friend std::istream& operator>>(std::istream& is, point& p) {
        return is >> p.x >> p.y;
    }
    friend std::ostream& operator<<(std::ostream& os, const point& p) {
        return os << p.x << ' ' << p.y;
    }
};
template <class T> T dot(const point<T>& a, const point<T>& b) {
    return a.x * b.x + a.y * b.y;
}
template <class T> T cross(const point<T>& a, const point<T>& b) {
    return a.x * b.y - a.y * b.x;
}
template <class T> T norm(const point<T>& p) {
    return p.x * p.x + p.y * p.y;
}
template <std::floating_point T> T abs(const point<T>& p) {
    return std::sqrt(norm(p));
}
template <std::floating_point T> T arg(const point<T>& p) {
    return std::atan2(p.y, p.x);
}
template <class T> point<T> normal(const point<T>& p) {
    return point<T>(-p.y, p.x);
}
template <std::floating_point T> point<T> unit(const point<T>& p) {
    return p / abs(p);
}
template <std::floating_point T> point<T> rotate(const point<T>& p, T theta) {
    T c = std::cos(theta), s = std::sin(theta);
    return point<T>(p.x * c - p.y * s, p.x * s + p.y * c);
}
template <class T> int ccw(const point<T>& a, const point<T>& b, const point<T>& c) {
    return sgn(cross(b - a, c - a));
}
struct arg_less {
    template <class T> static int group(const point<T>& p) {
        if (sgn(p.y) < 0) return 0;
        if (sgn(p.y) == 0 && sgn(p.x) >= 0) return 1;
        return 2;
    }
    template <class T> bool operator()(const point<T>& a, const point<T>& b) const {
        if (group(a) != group(b)) return group(a) < group(b);
        return sgn(cross(a, b)) > 0;
    }
};
}  // namespace geo
#line 9 "geometry/line.hpp"
namespace geo {
template <class T> struct line {
    point<T> a, b;
    line() {
    }
    line(const point<T>& a, const point<T>& b) : a(a), b(b) {
    }
    friend std::istream& operator>>(std::istream& is, line& l) {
        return is >> l.a >> l.b;
    }
};
template <class T> struct segment {
    point<T> a, b;
    segment() {
    }
    segment(const point<T>& a, const point<T>& b) : a(a), b(b) {
    }
    friend std::istream& operator>>(std::istream& is, segment& s) {
        return is >> s.a >> s.b;
    }
};
template <class T> line<T> to_line(const segment<T>& s) {
    return line<T>(s.a, s.b);
}
template <class T> bool parallel(const line<T>& l, const line<T>& m) {
    return sgn(cross(l.b - l.a, m.b - m.a)) == 0;
}
template <class T> bool orthogonal(const line<T>& l, const line<T>& m) {
    return sgn(dot(l.b - l.a, m.b - m.a)) == 0;
}
template <class T> bool on_line(const line<T>& l, const point<T>& p) {
    return ccw(l.a, l.b, p) == 0;
}
template <class T> bool on_segment(const segment<T>& s, const point<T>& p) {
    return ccw(s.a, s.b, p) == 0 && sgn(dot(s.a - p, s.b - p)) <= 0;
}
template <class T> bool intersect(const line<T>& l, const line<T>& m) {
    return !parallel(l, m) || on_line(l, m.a);
}
template <class T> bool intersect(const line<T>& l, const segment<T>& s) {
    return ccw(l.a, l.b, s.a) * ccw(l.a, l.b, s.b) <= 0;
}
template <class T> bool intersect(const segment<T>& s, const line<T>& l) {
    return intersect(l, s);
}
template <class T> bool intersect(const segment<T>& s, const segment<T>& t) {
    if (ccw(s.a, s.b, t.a) * ccw(s.a, s.b, t.b) > 0) return false;
    if (ccw(t.a, t.b, s.a) * ccw(t.a, t.b, s.b) > 0) return false;
    if (ccw(s.a, s.b, t.a) != 0 || ccw(s.a, s.b, t.b) != 0) return true;
    return on_segment(s, t.a) || on_segment(s, t.b) || on_segment(t, s.a) || on_segment(t, s.b);
}
template <std::floating_point T> point<T> projection(const line<T>& l, const point<T>& p) {
    point<T> d = l.b - l.a;
    return l.a + d * (dot(d, p - l.a) / norm(d));
}
template <std::floating_point T> point<T> reflection(const line<T>& l, const point<T>& p) {
    return projection(l, p) * T(2) - p;
}
template <std::floating_point T> line<T> perpendicular_bisector(const point<T>& a, const point<T>& b) {
    point<T> m = (a + b) / T(2);
    return line<T>(m, m + normal(b - a));
}
template <std::floating_point T> T dist(const point<T>& a, const point<T>& b) {
    return abs(a - b);
}
template <std::floating_point T> T dist(const line<T>& l, const point<T>& p) {
    return std::abs(cross(l.b - l.a, p - l.a)) / abs(l.b - l.a);
}
template <std::floating_point T> T dist(const point<T>& p, const line<T>& l) {
    return dist(l, p);
}
template <std::floating_point T> T dist(const segment<T>& s, const point<T>& p) {
    if (sgn(dot(s.b - s.a, p - s.a)) < 0) return abs(p - s.a);
    if (sgn(dot(s.a - s.b, p - s.b)) < 0) return abs(p - s.b);
    return dist(to_line(s), p);
}
template <std::floating_point T> T dist(const point<T>& p, const segment<T>& s) {
    return dist(s, p);
}
template <std::floating_point T> T dist(const line<T>& l, const line<T>& m) {
    if (intersect(l, m)) return 0;
    return dist(l, m.a);
}
template <std::floating_point T> T dist(const line<T>& l, const segment<T>& s) {
    if (intersect(l, s)) return 0;
    return std::min(dist(l, s.a), dist(l, s.b));
}
template <std::floating_point T> T dist(const segment<T>& s, const line<T>& l) {
    return dist(l, s);
}
template <std::floating_point T> T dist(const segment<T>& s, const segment<T>& t) {
    if (intersect(s, t)) return 0;
    return std::min({dist(s, t.a), dist(s, t.b), dist(t, s.a), dist(t, s.b)});
}
template <std::floating_point T> point<T> cross_point(const line<T>& l, const line<T>& m) {
    assert(!parallel(l, m));
    point<T> d = l.b - l.a;
    return l.a + d * (cross(m.b - m.a, m.a - l.a) / cross(m.b - m.a, d));
}
template <std::floating_point T> point<T> cross_point(const line<T>& l, const segment<T>& s) {
    assert(intersect(l, s));
    return cross_point(l, to_line(s));
}
template <std::floating_point T> point<T> cross_point(const segment<T>& s, const line<T>& l) {
    return cross_point(l, s);
}
template <std::floating_point T> point<T> cross_point(const segment<T>& s, const segment<T>& t) {
    assert(intersect(s, t));
    return cross_point(to_line(s), to_line(t));
}
}  // namespace geo
#line 11 "geometry/polygon.hpp"
namespace geo {
template <class T> using polygon = std::vector<point<T>>;
template <class T> T area2(const polygon<T>& poly) {
    int n = int(poly.size());
    if constexpr (std::is_integral_v<T>) {
        __int128 ret = 0;
        for (int i = 0; i < n; i++) ret += (__int128)cross(poly[i], poly[(i + 1) % n]);
        return T(ret);
    } else {
        T ret = 0;
        for (int i = 0; i < n; i++) ret += cross(poly[i], poly[(i + 1) % n]);
        return ret;
    }
}
template <std::floating_point T> T area(const polygon<T>& poly) {
    return area2(poly) / T(2);
}
template <class T> bool is_convex(const polygon<T>& poly) {
    int n = int(poly.size());
    for (int i = 0; i < n; i++) {
        if (ccw(poly[i], poly[(i + 1) % n], poly[(i + 2) % n]) < 0) return false;
    }
    return true;
}
template <class T> int contains(const polygon<T>& poly, const point<T>& p) {
    int n = int(poly.size());
    bool in = false;
    for (int i = 0; i < n; i++) {
        point<T> a = poly[i] - p, b = poly[(i + 1) % n] - p;
        if (ccw(a, b, point<T>()) == 0 && sgn(dot(a, b)) <= 0) return 1;
        if (sgn(a.y) > sgn(b.y)) std::swap(a, b);
        if (sgn(a.y) <= 0 && 0 < sgn(b.y) && sgn(cross(a, b)) > 0) in = !in;
    }
    return in ? 2 : 0;
}
template <class T> polygon<T> convex_hull(std::vector<point<T>> ps, bool strict = true) {
    int n = int(ps.size());
    std::sort(ps.begin(), ps.end());
    ps.erase(std::unique(ps.begin(), ps.end()), ps.end());
    n = int(ps.size());
    if (n <= 2) return ps;
    polygon<T> ret(2 * n);
    int k = 0;
    for (int i = 0; i < n; ret[k++] = ps[i++]) {
        while (k >= 2 && ccw(ret[k - 2], ret[k - 1], ps[i]) < (strict ? 1 : 0)) k--;
    }
    for (int i = n - 2, t = k + 1; i >= 0; ret[k++] = ps[i--]) {
        while (k >= t && ccw(ret[k - 2], ret[k - 1], ps[i]) < (strict ? 1 : 0)) k--;
    }
    ret.resize(k - 1);
    return ret;
}
template <class T> std::pair<int, int> farthest_pair(const std::vector<point<T>>& ps) {
    int n = int(ps.size());
    assert(n >= 2);
    polygon<T> ch = convex_hull(ps);
    int m = int(ch.size());
    auto index_of = [&](const point<T>& p) {
        return int(std::find(ps.begin(), ps.end(), p) - ps.begin());
    };
    if (m == 1) return {0, 1};
    if (m == 2) return {index_of(ch[0]), index_of(ch[1])};
    std::pair<int, int> ret = {0, 1};
    T best = norm(ch[0] - ch[1]);
    for (int i = 0, j = 1; i < m; i++) {
        while (norm(ch[i] - ch[(j + 1) % m]) > norm(ch[i] - ch[j])) j = (j + 1) % m;
        if (norm(ch[i] - ch[j]) > best) {
            best = norm(ch[i] - ch[j]);
            ret = {i, j};
        }
    }
    return {index_of(ch[ret.first]), index_of(ch[ret.second])};
}
template <class T> std::pair<int, int> closest_pair(const std::vector<point<T>>& ps) {
    int n = int(ps.size());
    assert(n >= 2);
    std::vector<int> id(n);
    for (int i = 0; i < n; i++) id[i] = i;
    std::sort(id.begin(), id.end(), [&](int i, int j) { return ps[i] < ps[j]; });
    std::pair<int, int> ret = {id[0], id[1]};
    T best = norm(ps[id[0]] - ps[id[1]]);
    std::vector<int> buf(n);
    auto rec = [&](auto rec, int l, int r) -> void {
        if (r - l <= 1) return;
        int m = (l + r) / 2;
        T x = ps[id[m]].x;
        rec(rec, l, m);
        rec(rec, m, r);
        std::inplace_merge(id.begin() + l, id.begin() + m, id.begin() + r,
                           [&](int i, int j) { return ps[i].y < ps[j].y; });
        int k = 0;
        for (int t = l; t < r; t++) {
            int i = id[t];
            T dx = ps[i].x - x;
            if (dx * dx >= best) continue;
            for (int u = k - 1; u >= 0; u--) {
                int j = buf[u];
                T dy = ps[i].y - ps[j].y;
                if (dy * dy >= best) break;
                T d = norm(ps[i] - ps[j]);
                if (d < best) {
                    best = d;
                    ret = {j, i};
                }
            }
            buf[k++] = i;
        }
    };
    rec(rec, 0, n);
    return ret;
}
template <std::floating_point T> polygon<T> convex_cut(const polygon<T>& poly, const line<T>& l) {
    int n = int(poly.size());
    polygon<T> ret;
    for (int i = 0; i < n; i++) {
        const point<T>& p = poly[i];
        const point<T>& q = poly[(i + 1) % n];
        int cp = ccw(l.a, l.b, p), cq = ccw(l.a, l.b, q);
        if (cp >= 0) ret.push_back(p);
        if (cp * cq < 0) ret.push_back(cross_point(l, line<T>(p, q)));
    }
    return ret;
}
}  // namespace geo

Geometry: Polygon (geometry/polygon.hpp)

概要

使い方

制約

計算量

問題例

Depends on

Required by

Verified with

Code

Geometry: Polygon
(geometry/polygon.hpp)