1 幾何基礎

#include <cstdio>
#include <cmath>
using namespace std;
const double pi = acos(-1.0);
const double inf = 1e100;
const double eps = 1e-6;
int sgn(double d){
    if(fabs(d) < eps)
        return 0;
    if(d > 0)
        return 1;
    return -1;
}
int dcmp(double x, double y){
    if
 
(fabs(x - y) < eps)
        return 0;
    if(x > y)
        return 1;
    return -1;
}
int main() {
    double x = 1.49999;
    int fx = floor(x);//向下取整函式
    int cx = ceil(x);//向上取整函式
    int rx = round(x);//四捨五入函式
    printf("%f %d %d %d\n", x, fx, cx, rx);
    //輸出結果 1.499990 1 2 1
    return  0 ;
}

2 點與向量

2.1 手動實現

struct Point{
    double x, y;
    Point(double x = 0, double y = 0):x(x),y(y){}
};
typedef Point Vector;
Vector operator + (Vector A, Vector B){
    return Vector(A.x+B.x, A.y+B.y);
}
Vector operator - (Point A, Point B){
    return Vector(A.x-B.x, A.y-B.y);
}
Vector operator * (Vector A, double
 
 p){
    return Vector(A.x*p, A.y*p);
}
Vector operator / (Vector A, double p){
    return Vector(A.x/p, A.y/p);
}
bool operator < (const Point& a, const Point& b){
    if(a.x == b.x)
        return a.y < b.y;
    return a.x < b.x;
}
const double eps = 1e-6;
int sgn(double x){
    if(fabs(x) < eps)
        return 0;
    if(x < 0)
        return -1;
    return 1;
}
bool operator == (const Point& a, const Point& b){
    if(sgn(a.x-b.x) == 0 && sgn(a.y-b.y) == 0)
        return true;
    return false;
}
double Dot(Vector A, Vector B){
    return A.x*B.x + A.y*B.y;
}
double Length(Vector A){
    return sqrt(Dot(A, A));
}
double Angle(Vector A, Vector B){
    return acos(Dot(A, B)/Length(A)/Length(B));
}
double Cross(Vector A, Vector B){
    return A.x*B.y-A.y*B.x;
}
double Area2(Point A, Point B, Point C){
    return Cross(B-A, C-A);
}
Vector Rotate(Vector A, double rad){//rad為弧度 且為逆時針旋轉的角
    return Vector(A.x*cos(rad)-A.y*sin(rad), A.x*sin(rad)+A.y*cos(rad));
}
Vector Normal(Vector A){//向量A左轉90°的單位法向量
    double L = Length(A);
    return Vector(-A.y/L, A.x/L);
}
bool ToLeftTest(Point a, Point b, Point c){
    return Cross(b - a, c - b) > 0;
}

2.2 複數黑科技

#include <complex>
using namespace std;
typedef complex<double> Point;
typedef Point Vector;//複數定義向量後，自動擁有建構函式、加減法和數量積
const double eps = 1e-9;
int sgn(double x){
    if(fabs(x) < eps)
        return 0;
    if(x < 0)
        return -1;
    return 1;
}
double Length(Vector A){
    return abs(A);
}
double Dot(Vector A, Vector B){//conj(a+bi)返回共軛複數a-bi
    return real(conj(A)*B);
}
double Cross(Vector A, Vector B){
    return imag(conj(A)*B);
}
Vector Rotate(Vector A, double rad){
    return A*exp(Point(0, rad));//exp(p)返回以e為底複數的指數
}

3 點與線

3.1 直線定義

struct Line{//直線定義
    Point v, p;
    Line(Point v, Point p):v(v), p(p) {}
    Point point(double t){//返回點P = v + (p - v)*t
        return v + (p - v)*t;
    }
};

3.2 求兩直線交點

//呼叫前需保證 Cross(v, w) != 0
Point GetLineIntersection(Point P, Vector v, Point Q, Vector w){
    Vector u = P-Q;
    double t = Cross(w, u)/Cross(v, w);
    return P+v*t;
}

3.3 求點到直線距離

//點P到直線AB距離
double DistanceToLine(Point P, Point A, Point B){
    Vector v1 = B-A, v2 = P-A;
    return fabs(Cross(v1, v2)/Length(v1));
}//不去絕對值，得到的是有向距離

3.4 求點到線段距離

//點P到線段AB距離公式
double DistanceToSegment(Point P, Point A, Point B){
    if(A == B)
        return Length(P-A);
    Vector v1 = B-A, v2 = P-A, v3 = P-B;
    if(dcmp(Dot(v1, v2)) < 0)
        return Length(v2);
    if(dcmp(Dot(v1, v3)) > 0)
        return Length(v3);
    return DistanceToLine(P, A, B);
}

3.5 求點在直線上的投影點

//點P在直線AB上的投影點
Point GetLineProjection(Point P, Point A, Point B){
    Vector v = B-A;
    return A+v*(Dot(v, P-A)/Dot(v, v));
}

3.6 判斷點是否線上段上

//判斷p點是否線上段a1a2上
bool OnSegment(Point p, Point a1, Point a2){
    return dcmp(Cross(a1-p, a2-p)) == 0 && dcmp(Dot(a1-p, a2-p)) < 0;
}

3.7 判斷兩線段是否相交

//判斷兩線段是否相交
bool SegmentProperIntersection(Point a1, Point a2, Point b1, Point b2){
    double c1 = Cross(a2-a1, b1-a1), c2 = Cross(a2-a1, b2-a1);
    double c3 = Cross(b2-b1, a1-b1), c4 = Cross(b2-b1, a2-b1);
    //if判斷控制是否允許線段在端點處相交，根據需要新增
    if(!sgn(c1) || !sgn(c2) || !sgn(c3) || !sgn(c4)){
        bool f1 = OnSegment(b1, a1, a2);
        bool f2 = OnSegment(b2, a1, a2);
        bool f3 = OnSegment(a1, b1, b2);
        bool f4 = OnSegment(a2, b1, b2);
        bool f = (f1|f2|f3|f4);
        return f;
    }
    return (sgn(c1)*sgn(c2) < 0 && sgn(c3)*sgn(c4) < 0);
}

4 多邊形

4.1 求多邊形有向面積

//多邊形有向面積
double PolygonArea(Point* p, int n){//p為端點集合，n為端點個數
    double s = 0;
    for(int i = 1; i < n-1; ++i)
        s += Cross(p[i]-p[0], p[i+1]-p[0]);
    return s;
}

4.2 判斷點是否在多邊形內

//判斷點是否在多邊形內，若點在多邊形內返回1，在多邊形外部返回0，在多邊形上返回-1
int isPointInPolygon(Point p, vector<Point> poly){
    int wn = 0;
    int n = poly.size();
    for(int i = 0; i < n; ++i){
        if(OnSegment(p, poly[i], poly[(i+1)%n])) return -1;
        int k = sgn(Cross(poly[(i+1)%n] - poly[i], p - poly[i]));
        int d1 = sgn(poly[i].y - p.y);
        int d2 = sgn(poly[(i+1)%n].y - p.y);
        if(k > 0 && d1 <= 0 && d2 > 0) wn++;
        if(k < 0 && d2 <= 0 && d1 > 0) wn--;
    }
    if(wn != 0)
        return 1;
    return 0;
}

5 圓

5.1 求圓與直線交點

int getLineCircleIntersection(Line L, Circle C, double& t1, double& t2, vector<Point>& sol){
    double a = L.v.x, b = L.p.x - C.c.x, c = L.v.y, d = L.p.y - C.c.y;
    double e = a*a + c*c, f = 2*(a*b + c*d), g = b*b + d*d - C.r*C.r;
    double delta = f*f - 4*e*g;//判別式
    if(sgn(delta) < 0)//相離
        return 0;
    if(sgn(delta) == 0){//相切
        t1 = -f /(2*e);
        t2 = -f /(2*e);
        sol.push_back(L.point(t1));//sol存放交點本身
        return 1;
    }
    //相交
    t1 = (-f - sqrt(delta))/(2*e);
    sol.push_back(L.point(t1));
    t2 = (-f + sqrt(delta))/(2*e);
    sol.push_back(L.point(t2));
    return 2;
}

5.2 求兩圓交點

5.3 求點到圓的切線

5.4 求兩圓公切線

5.5 求兩圓相交面積

double AreaOfOverlap(Point c1, double r1, Point c2, double r2){
    double d = Length(c1 - c2);
    if(r1 + r2 < d + eps)
        return 0.0;
    if(d < fabs(r1 - r2) + eps){
        double r = min(r1, r2);
        return pi*r*r;
    }
    double x = (d*d + r1*r1 - r2*r2)/(2.0*d);
    double p = (r1 + r2 + d)/2.0;
    double t1 = acos(x/r1);
    double t2 = acos((d - x)/r2);
    double s1 = r1*r1*t1;
    double s2 = r2*r2*t2;
    double s3 = 2*sqrt(p*(p - r1)*(p - r2)*(p - d));
    return s1 + s2 - s3;
}

5.6 用給定半徑的圓覆蓋最多的點

6 凸包

6.1 GrahamScan O(nlogn)

const int maxn = 1e3 + 5;
const double eps = 1e-9;
struct Point {
    double x, y;
    Point(double x = 0, double y = 0):x(x),y(y){}
};
typedef Point Vector;
Point lst[maxn];
int stk[maxn], top;
Vector operator - (Point A, Point B){
    return Vector(A.x-B.x, A.y-B.y);
}
int sgn(double x){
    if(fabs(x) < eps)
        return 0;
    if(x < 0)
        return -1;
    return 1;
}
double Cross(Vector v0, Vector v1) {
    return v0.x*v1.y - v1.x*v0.y;
}
double Dis(Point p1, Point p2) { //計算 p1p2的 距離
    return sqrt((p2.x-p1.x)*(p2.x-p1.x)+(p2.y-p1.y)*(p2.y-p1.y));
}
bool cmp(Point p1, Point p2) { //極角排序函式 ，角度相同則距離小的在前面
    int tmp = sgn(Cross(p1 - lst[0], p2 - lst[0]));
    if(tmp > 0)
        return true;
    if(tmp == 0 && Dis(lst[0], p1) < Dis(lst[0], p2))
        return true;
    return false;
}
//點的編號0 ~ n - 1
//返回凸包結果stk[0 ~ top - 1]為凸包的編號
void Graham(int n) {
    int k = 0;
    Point p0;
    p0.x = lst[0].x;
    p0.y = lst[0].y;
    for(int i = 1; i < n; ++i) {
        if( (p0.y > lst[i].y) || ((p0.y == lst[i].y) && (p0.x > lst[i].x)) ) {
            p0.x = lst[i].x;
            p0.y = lst[i].y;
            k = i;
        }
    }
    lst[k] = lst[0];
    lst[0] = p0;
    sort(lst + 1, lst + n, cmp);
    if(n == 1) {
        top = 1;
        stk[0] = 0;
        return ;
    }
    if(n == 2) {
        top = 2;
        stk[0] = 0;
        stk[1] = 1;
        return ;
    }
    stk[0] = 0;
    stk[1] = 1;
    top = 2;
    for(int i = 2; i < n; ++i) {
        while(top > 1 && Cross(lst[stk[top - 1]] - lst[stk[top - 2]], lst[i] - lst[stk[top - 2]]) <= 0)
            --top;
        stk[top] = i;
        ++top;
    }
    return ;
}

6.2 Andrew O(nlogn)

struct Point {
    double x, y;
    Point(double x = 0, double y = 0):x(x),y(y){}
};
typedef Point Vector;
Vector operator - (Point A, Point B){
    return Vector(A.x-B.x, A.y-B.y);
}
bool operator < (const Point& a, const Point& b){
    if(a.x == b.x)
        return a.y < b.y;
    return a.x < b.x;
}
double Cross(Vector v0, Vector v1) {
    return v0.x*v1.y - v1.x*v0.y;
}
//計算凸包，輸入點陣列為 p，個數為 n， 輸出點陣列為 ch。函式返回凸包頂點數
//如果不希望凸包的邊上有輸入點，則把兩個 <= 改為 <
//在精度要求高時建議用dcmp比較
//輸入不能有重複點，函式執行完後輸入點的順序被破壞
int ConvexHull(Point* p, int n, Point* ch) {
    sort(p, p+n);
    int m = 0;
    for(int i = 0; i < n; ++i) {
        while(m > 1 && Cross(ch[m-1] - ch[m-2], p[i] - ch[m-2]) <= 0) {
            m--;
        }
        ch[m++] = p[i];
    }
    int k = m;
    for(int i = n-2; i>= 0; --i) {
        while(m > k && Cross(ch[m-1] - ch[m-2], p[i] - ch[m-2]) <= 0) {
            m--;
        }
        ch[m++] = p[i];
    }
    if(n > 1)
        --m;
    return m;
}

7 半平面交 O(nlogn)

const double eps = 1e-6;
struct Point{
    double x, y;
    Point(double x = 0, double y = 0):x(x),y(y){}
};
typedef Point Vector;
Vector operator + (Vector A, Vector B){
    return Vector(A.x+B.x, A.y+B.y);
}
Vector operator - (Point A, Point B){
    return Vector(A.x-B.x, A.y-B.y);
}
Vector operator * (Vector A, double p){
    return Vector(A.x*p, A.y*p);
}
int sgn(double x){
    if(fabs(x) < eps)
        return 0;
    if(x < 0)
        return -1;
    return 1;
}
double Dot(Vector A, Vector B){
    return A.x*B.x + A.y*B.y;
}
double Cross(Vector A, Vector B){
    return A.x*B.y-A.y*B.x;
}
double Length(Vector A){
    return sqrt(Dot(A, A));
}
Vector Normal(Vector A){//向量A左轉90°的單位法向量
    double L = Length(A);
    return Vector(-A.y/L, A.x/L);
}
struct Line{
    Point p;//直線上任意一點
    Vector v;//方向向量，它的左邊就是對應的半平面
    double ang;//極角，即從x軸正半軸旋轉到向量v所需要的角（弧度）
    Line(){}
    Line(Point p, Vector v) : p(p), v(v){
        ang = atan2(v.y, v.x);
    }
    bool operator < (const Line& L) const {//排序用的比較運算子
        return ang < L.ang;
    }
};
//點p在有向直線L的左側
bool OnLeft(Line L, Point p){
    return Cross(L.v, p - L.p) > 0;
}
//兩直線交點。假定交點唯一存在
Point GetIntersection(Line a, Line b){
    Vector u = a.p - b.p;
    double t = Cross(b.v, u)/Cross(a.v, b.v);
    return a.p + a.v*t;
}
//半平面交的主過程
int HalfplaneIntersection(Line* L, int n, Point* poly){
    sort(L, L + n);//按照極角排序
    int fst = 0, lst = 0;//雙端佇列的第一個元素和最後一個元素
    Point *P = new Point[n];//p[i] 為 q[i]與q[i + 1]的交點
    Line *q = new Line[n];//雙端佇列
    q[fst = lst = 0] = L[0];//初始化為只有一個半平面L[0]
    for(int i = 1; i < n; ++i){
        while(fst < lst && !OnLeft(L[i], P[lst - 1])) --lst;
        while(fst < lst && !OnLeft(L[i], P[fst])) ++fst;
        q[++lst] = L[i];
        if(sgn(Cross(q[lst].v, q[lst - 1].v)) == 0){
            //兩向量平行且同向，取內側一個
            --lst;
            if(OnLeft(q[lst], L[i].p)) q[lst] = L[i];
        }
        if(fst < lst)
            P[lst - 1] = GetIntersection(q[lst - 1], q[lst]);
    }
    while(fst < lst && !OnLeft(q[fst], P[lst - 1])) --lst;
    //刪除無用平面
    if(lst - fst <= 1) return 0;//空集
    P[lst] = GetIntersection(q[lst], q[fst]);//計算首尾兩個半平面的交點
    //從deque複製到輸出中
    int m = 0;
    for(int i = fst; i <= lst; ++i) poly[m++] = P[i];
    return m;
}

8 平面最近點對 O(nlogn)

const double eps = 1e-6;
struct Point{
    double x, y;
    Point(double x = 0, double y = 0):x(x),y(y){}
};
typedef Point Vector;
Vector operator + (Vector A, Vector B){
    return Vector(A.x+B.x, A.y+B.y);
}
Vector operator - (Point A, Point B){
    return Vector(A.x-B.x, A.y-B.y);
}
Vector operator * (Vector A, double p){
    return Vector(A.x*p, A.y*p);
}
int sgn(double x){
    if(fabs(x) < eps)
        return 0;
    if(x < 0)
        return -1;
    return 1;
}
double Dot(Vector A, Vector B){
    return A.x*B.x + A.y*B.y;
}
double Cross(Vector A, Vector B){
    return A.x*B.y-A.y*B.x;
}
double Length(Vector A){
    return sqrt(Dot(A, A));
}
struct Line{
    Point p;//直線上任意一點
    Vector v;//方向向量，它的左邊就是對應的半平面
    double ang;//極角，即從x軸正半軸旋轉到向量v所需要的角（弧度）
    Line(){}
    Line(Point p, Vector v) : p(p), v(v){
        ang = atan2(v.y, v.x);
    }
    bool operator < (const Line& L) const {//排序用的比較運算子
        return ang < L.ang;
    }
};
//點p在有向直線L的左側
bool OnLeft(Line L, Point p){
    return Cross(L.v, p - L.p) > 0;
}
//兩直線交點。假定交點唯一存在
Point GetIntersection(Line a, Line b){
    Vector u = a.p - b.p;
    double t = Cross(b.v, u)/Cross(a.v, b.v);
    return a.p + a.v*t;
}
//半平面交的主過程
int HalfplaneIntersection(Line* L, int n, Point* poly){
    sort(L, L + n);//按照極角排序
    int fst = 0, lst = 0;//雙端佇列的第一個元素和最後一個元素
    Point *P = new Point[n];//p[i] 為 q[i]與q[i + 1]的交點
    Line *q = new Line[n];//雙端佇列
    q[fst = lst = 0] = L[0];//初始化為只有一個半平面L[0]
    for(int i = 1; i < n; ++i){
        while(fst < lst && !OnLeft(L[i], P[lst - 1])) --lst;
        while(fst < lst && !OnLeft(L[i], P[fst])) ++fst;
        q[++lst] = L[i];
        if(sgn(Cross(q[lst].v, q[lst - 1].v)) == 0){
            //兩向量平行且同向，取內側一個
            --lst;
            if(OnLeft(q[lst], L[i].p)) q[lst] = L[i];
        }
        if(fst < lst)
            P[lst - 1] = GetIntersection(q[lst - 1], q[lst]);
    }
    while(fst < lst && !OnLeft(q[fst], P[lst - 1])) --lst;
    //刪除無用平面
    if(lst - fst <= 1) return 0;//空集
    P[lst] = GetIntersection(q[lst], q[fst]);//計算首尾兩個半平面的交點
    //從deque複製到輸出中
    int m = 0;
    for(int i = fst; i <= lst; ++i) poly[m++] = P[i];
    return m;
}

9 旋轉卡殼 O(nlogn)

double Dist2(Point p1, Point p2) { //計算距離的平方
    double ret = Dot(p1 - p2, p1 - p2);
    return ret;
}
double RotatingCalipers(Point* ch, int m) {//返回平面最大距離的平方
    if(m == 1) return 0.0;
    if(m == 2) return Dist2(ch[0], ch[1]);
    double ret = 0.0;
    ch[m] = ch[0];
    int j = 2;
    for(int i = 0; i < m; ++i) {
        while(Cross(ch[i + 1] - ch[i], ch[j] - ch[i]) < Cross(ch[i + 1] - ch[i], ch[j + 1] - ch[i]))
            j = (j + 1)%m;
        ret = max(ret, max(Dist2(ch[j], ch[i]), Dist2(ch[j], ch[i + 1])));
    }
    return ret;
}

10 三點確定外接圓圓心座標

struct Point {
    double x,y;
    Point(double x = 0, double y = 0):x(x),y(y){}
};
Point Excenter(Point a, Point b, Point c){
    double a1 = b.x - a.x;
    double b1 = b.y - a.y;
    double c1 = (a1*a1 + b1*b1)/2;
    double a2 = c.x - a.x;
    double b2 = c.y - a.y;
    double c2 = (a2*a2 + b2*b2)/2;
    double d = a1*b2 - a2*b1;
    return Point(a.x + (c1*b2 - c2*b1)/d, a.y + (a1*c2 - a2*c1)/d);
}