ICPC Competitive Programming
← กลับหน้าหลัก

พื้นฐาน Geometry

พื้นฐาน Computational Geometry สำหรับ Competitive Programming

พื้นฐาน Geometry

Basic Types

Point / Vector

struct Point {
    double x, y;
    
    Point() : x(0), y(0) {}
    Point(double x, double y) : x(x), y(y) {}
    
    Point operator+(const Point& p) const { return Point(x + p.x, y + p.y); }
    Point operator-(const Point& p) const { return Point(x - p.x, y - p.y); }
    Point operator*(double t) const { return Point(x * t, y * t); }
    Point operator/(double t) const { return Point(x / t, y / t); }
    
    double dot(const Point& p) const { return x * p.x + y * p.y; }
    double cross(const Point& p) const { return x * p.y - y * p.x; }
    
    double norm() const { return sqrt(x * x + y * y); }
    double norm2() const { return x * x + y * y; }
    
    Point rotate(double angle) const {
        return Point(x * cos(angle) - y * sin(angle),
                     x * sin(angle) + y * cos(angle));
    }
    
    bool operator<(const Point& p) const {
        if (x != p.x) return x < p.x;
        return y < p.y;
    }
    
    bool operator==(const Point& p) const {
        return abs(x - p.x) < EPS && abs(y - p.y) < EPS;
    }
};

const double EPS = 1e-9;

Line

struct Line {
    Point a, b;  // line through points a and b
    
    // หรือใช้ ax + by + c = 0
    double A, B, C;
    
    Line(Point a, Point b) : a(a), b(b) {
        A = b.y - a.y;
        B = a.x - b.x;
        C = -(A * a.x + B * a.y);
    }
    
    // Distance from point to line
    double dist(Point p) const {
        return abs(A * p.x + B * p.y + C) / sqrt(A * A + B * B);
    }
    
    // Projection of point onto line
    Point project(Point p) const {
        Point v = b - a;
        double t = (p - a).dot(v) / v.norm2();
        return a + v * t;
    }
    
    // Reflection of point across line
    Point reflect(Point p) const {
        return project(p) * 2 - p;
    }
};

Basic Operations

Cross Product และ Orientation

// cross(b-a, c-a) > 0: counter-clockwise
// cross(b-a, c-a) < 0: clockwise
// cross(b-a, c-a) = 0: collinear
int sign(double x) {
    if (x > EPS) return 1;
    if (x < -EPS) return -1;
    return 0;
}

int orientation(Point a, Point b, Point c) {
    return sign((b - a).cross(c - a));
}

Distance Functions

// Distance between two points
double dist(Point a, Point b) {
    return (a - b).norm();
}

// Distance from point to line segment
double distToSegment(Point p, Point a, Point b) {
    if ((p - a).dot(b - a) <= 0) return dist(p, a);
    if ((p - b).dot(a - b) <= 0) return dist(p, b);
    return abs((b - a).cross(p - a)) / dist(a, b);
}

Line Intersection

// Check if two lines intersect (not parallel)
bool lineIntersect(Line l1, Line l2, Point& result) {
    double det = l1.A * l2.B - l2.A * l1.B;
    if (abs(det) < EPS) return false;  // parallel
    
    result.x = (l1.B * l2.C - l2.B * l1.C) / det;
    result.y = (l2.A * l1.C - l1.A * l2.C) / det;
    return true;
}

// Alternative using parametric form
bool lineIntersect2(Point a, Point b, Point c, Point d, Point& result) {
    Point ab = b - a, cd = d - c, ac = c - a;
    double denom = ab.cross(cd);
    if (abs(denom) < EPS) return false;
    
    double t = ac.cross(cd) / denom;
    result = a + ab * t;
    return true;
}

Segment Intersection

bool onSegment(Point p, Point a, Point b) {
    return min(a.x, b.x) <= p.x + EPS && p.x <= max(a.x, b.x) + EPS &&
           min(a.y, b.y) <= p.y + EPS && p.y <= max(a.y, b.y) + EPS;
}

bool segmentIntersect(Point a, Point b, Point c, Point d) {
    int d1 = orientation(c, d, a);
    int d2 = orientation(c, d, b);
    int d3 = orientation(a, b, c);
    int d4 = orientation(a, b, d);
    
    if (((d1 > 0 && d2 < 0) || (d1 < 0 && d2 > 0)) &&
        ((d3 > 0 && d4 < 0) || (d3 < 0 && d4 > 0)))
        return true;
    
    // Collinear cases
    if (d1 == 0 && onSegment(a, c, d)) return true;
    if (d2 == 0 && onSegment(b, c, d)) return true;
    if (d3 == 0 && onSegment(c, a, b)) return true;
    if (d4 == 0 && onSegment(d, a, b)) return true;
    
    return false;
}

Polygon Operations

Area of Polygon (Shoelace Formula)

double polygonArea(vector<Point>& poly) {
    double area = 0;
    int n = poly.size();
    
    for (int i = 0; i < n; i++) {
        area += poly[i].cross(poly[(i + 1) % n]);
    }
    
    return abs(area) / 2;
}

Point in Polygon

// Ray casting algorithm
bool pointInPolygon(Point p, vector<Point>& poly) {
    int n = poly.size();
    int count = 0;
    
    for (int i = 0; i < n; i++) {
        Point a = poly[i], b = poly[(i + 1) % n];
        
        if ((a.y <= p.y && p.y < b.y) || (b.y <= p.y && p.y < a.y)) {
            double x = a.x + (p.y - a.y) * (b.x - a.x) / (b.y - a.y);
            if (p.x < x) count++;
        }
    }
    
    return count % 2 == 1;
}

// Winding number (handles edges)
int windingNumber(Point p, vector<Point>& poly) {
    int n = poly.size();
    int wn = 0;
    
    for (int i = 0; i < n; i++) {
        Point a = poly[i], b = poly[(i + 1) % n];
        
        if (a.y <= p.y) {
            if (b.y > p.y && orientation(a, b, p) > 0) wn++;
        } else {
            if (b.y <= p.y && orientation(a, b, p) < 0) wn--;
        }
    }
    
    return wn;  // != 0 means inside
}

Polygon Convexity Check

bool isConvex(vector<Point>& poly) {
    int n = poly.size();
    bool hasPos = false, hasNeg = false;
    
    for (int i = 0; i < n; i++) {
        int o = orientation(poly[i], poly[(i+1)%n], poly[(i+2)%n]);
        if (o > 0) hasPos = true;
        if (o < 0) hasNeg = true;
    }
    
    return !(hasPos && hasNeg);
}

Circle Operations

Circle

struct Circle {
    Point c;
    double r;
    
    Circle() : r(0) {}
    Circle(Point c, double r) : c(c), r(r) {}
    
    // Check if point is inside
    bool contains(Point p) const {
        return dist(c, p) <= r + EPS;
    }
    
    // Circumference
    double circumference() const {
        return 2 * M_PI * r;
    }
    
    // Area
    double area() const {
        return M_PI * r * r;
    }
};

Circle-Line Intersection

vector<Point> circleLineIntersect(Circle c, Line l) {
    Point proj = l.project(c.c);
    double d = dist(c.c, proj);
    
    if (d > c.r + EPS) return {};
    if (abs(d - c.r) < EPS) return {proj};
    
    double h = sqrt(c.r * c.r - d * d);
    Point v = (l.b - l.a) / (l.b - l.a).norm();
    
    return {proj + v * h, proj - v * h};
}

Circle-Circle Intersection

vector<Point> circleCircleIntersect(Circle c1, Circle c2) {
    double d = dist(c1.c, c2.c);
    
    if (d > c1.r + c2.r + EPS) return {};  // too far
    if (d < abs(c1.r - c2.r) - EPS) return {};  // one inside other
    
    double a = (c1.r * c1.r - c2.r * c2.r + d * d) / (2 * d);
    double h = sqrt(max(0.0, c1.r * c1.r - a * a));
    
    Point p = c1.c + (c2.c - c1.c) * (a / d);
    Point v = (c2.c - c1.c).rotate(M_PI / 2) * (h / d);
    
    if (abs(h) < EPS) return {p};
    return {p + v, p - v};
}

Tips

  1. Epsilon comparison - floating point ต้องใช้ epsilon
  2. Integer geometry - ถ้าทำได้ ใช้ integer และหลีกเลี่ยง division
  3. Cross product sign - สำคัญมากสำหรับ orientation tests
  4. Degenerate cases - จัดการ cases พิเศษเช่น parallel lines, coincident points