计算几何全家桶

最新推荐文章于 2023-02-09 13:41:59 发布
最新推荐文章于 2023-02-09 13:41:59 发布
阅读量209 收藏 2
点赞数 2
分类专栏： 算法 文章标签： 算法
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/messywind/article/details/119373554
版权

二维计算几何

#include <bits/stdc++.h>
#define x first
#define y second
using namespace std;
const int N = 1e5 + 5;
const double pi = acos(-1), eps = 1e-8;
typedef pair<int, int> PII;
typedef pair<double, double> PDD;
int sign(double x) {
    if (fabs(x) < eps) return 0;
    if (x < 0) return -1;
    return 1;
}
int cmp(double x, double y) {
    if (fabs(x - y) < eps) return 0;
    if (x < y) return -1;
    return 1;
}
struct Point {
    double x, y;
    Point(double x = 0, double y = 0) : x(x), y(y) {}
    bool operator < (const Point &p) {
        if (x == p.x) return y < p.y;
        return x < p.x;
    }
    Point operator + (Point a) {
        return Point(x + a.x, y + a.y);
    }
    Point operator - (Point a) {
        return Point(x - a.x, y - a.y);
    }
    Point operator * (double t) {
        return Point(x * t, y * t);
    }
    Point operator / (double t) {
        return Point(x / t, y / t);
    }
    bool operator == (Point a) {
        return x == a.x && y == a.y;
    }
    double lenth() { //模
        return sqrt(x * x + y * y);
    }
};
struct Line {
    Point st, ed;
    Line() {}
    Line(Point a, Point b){
        this->st = a;
        this->ed = b;
    }
};
struct Circle {
    Point p;
    double r;
    Circle() {}
    Circle(Point p, double r){
        this->p = p;
        this->r = r;
    }
};
double get_dist(Point a, Point b) { //两点距离
    double dx = a.x - b.x;
    double dy = a.y - b.y;
    return sqrt(dx * dx + dy * dy);
}
double get_angle_line(const Line &a) { //线的角度
    return atan2(a.ed.y - a.st.y, a.ed.x - a.st.x);
}
double dot(Point a, Point b) { //点积
    return a.x * b.x + a.y * b.y;
}
double cross(Point a, Point b) { //叉积
    return a.x * b.y - b.x * a.y;
}
double get_angle(Point a, Point b) { //两向量夹角
    return acos(dot(a, b) / a.lenth() / b.lenth());
}
double area(Point a, Point b, Point c) {
    return cross(b - a, c - a);
}
Point norm(Point a) { //单位向量
    return a / a.lenth();
}
Point rotate(Point a, double angle) { //点的旋转
    return Point(a.x * cos(angle) + a.y * sin(angle), -a.x * sin(angle) + a.y * cos(angle));
}
Point get_line_intersection(Point p, Point v, Point q, Point w) { //两直线交点
    Point u = p - q, res;
    double t = cross(w, u) / cross(v, w);
    return {p.x + v.x * t, p.y + v.y * t};
}
Point get_line_intersection(Line a, Line b) { //两直线交点
    return get_line_intersection(a.st, a.ed - a.st, b.st, b.ed - b.st);
}
bool on_right(Line &a, Line &b, Line &c) { //bc的交点是否在a的右侧
    auto o = get_line_intersection(b, c);
    return sign(area(a.st, a.ed, o)) <= 0;
}
double distance_to_line(Point p, Point a, Point b) { //点到线的距离
    Point v1 = b - a, v2 = p - a;
    return fabs(cross(v1, v2) / v1.lenth());
}
double distance_to_segment(Point p, Point a, Point b) { //点到线段的距离
    if (a == b) return (p - a).lenth();
    Point v1 = b - a, v2 = p - a, v3 = p - b;
    if (sign(dot(v1, v2)) < 0) return v2.lenth();
    if (sign(dot(v1, v3)) > 0) return v3.lenth();
    return distance_to_line(p, a, b);
}
Point get_line_projection(Point p, Point a, Point b) { //点在线上的投影
    Point v = b - a;
    return a + v * (dot(v, p - a) / dot(v, v));
}
double project(Point a, Point b, Point c) { //ac在ab上的投影
    return dot(b - a, c - a) / (b - a).lenth();
}
bool on_segment(Point p, Point a, Point b) { //点是否在线段上
    return sign(cross(p - a, p - b)) == 0 && sign(dot(p - a, p - b)) <= 0;
}
bool segment_intersection(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, a2 - b1), c4 = cross(b2 - b1, a1 - b1);
    return sign(c1) * sign(c2) <= 0 && sign(c3) * sign(c4) <= 0;
}
bool cmp_line(const Line &a, const Line &b) { //线段按角度排序
    double A = get_angle_line(a), B = get_angle_line(b);
    if (!cmp(A, B)) return area(a.st, a.ed, b.ed) < 0;
    return A < B;
}
Point a[N], ans[N], r;
Line line[N];
Circle cir[N];
PII A[N];
PDD D[N];
double res, R;
int stk[N], top, n, q[N], cnt, m;
bool used[N];
double polygon_area() { //多边形面积
    double res = 0;
    for (int i = 0; i < n; i ++) res += cross(a[i], a[(i + 1) % n]);
    return fabs(0.5 * res);
}
void andrew() { //求凸包
    sort(a, a + n);
    for (int i = 0; i < n; i ++) {
        while (top >= 2 && area(a[stk[top - 1]], a[stk[top]], a[i]) <= 0) {
            if (area(a[stk[top - 1]], a[stk[top]], a[i]) < 0) {
                used[stk[top --]] = 0;
            } else {
                top --;
            }
        }
        stk[++ top] = i;
        used[i] = 1;
    }
    used[0] = 0;
    for (int i = n - 1; ~i; i --){
        if (used[i]) continue;
        while (top >= 2 && area(a[stk[top - 1]], a[stk[top]], a[i]) <= 0) top -- ;
        stk[++ top] = i;
    }
}
double rotating_calipers() { //旋转卡壳
    if (top <= 3) return get_dist(a[0], a[n - 1]);
    double res = 0;
    for (int i = 1, j = 3; i <= top; i ++) {
        auto d = a[stk[i]], e = a[stk[i + 1]];
        while (area(d, e, a[stk[j]]) < area(d, e, a[stk[j + 1]])) {
            j = (j + 1) % top;
            if (!j) j = 1;
        }
        res = max(res, max(get_dist(d, a[stk[j]]), get_dist(e, a[stk[j]])));
    }
    return res;
}
double half_plane_intersection() { //半平面交
    sort(line, line + cnt, cmp_line);
    int hh = 0, tt = -1;
    for (int i = 0; i < cnt; i ++) {
        if (i && !cmp(get_angle_line(line[i]), get_angle_line(line[i - 1]))) continue;
        while (hh + 1 <= tt && on_right(line[i], line[q[tt - 1]], line[q[tt]])) tt -- ;
        while (hh + 1 <= tt && on_right(line[i], line[q[hh]], line[q[hh + 1]])) hh ++ ;
        q[++ tt] = i;
    }
    while (hh + 1 <= tt && on_right(line[q[hh]], line[q[tt - 1]], line[q[tt]])) tt -- ;
    while (hh + 1 <= tt && on_right(line[q[tt]], line[q[hh]], line[q[hh + 1]])) hh ++ ;
    q[++ tt] = q[hh];
    int k = 0;
    for (int i = hh; i < tt; i ++) ans[k ++] = get_line_intersection(line[q[i]], line[q[i + 1]]);
    double res = 0;
    for (int i = 1; i + 1 < k; i ++) res += area(ans[0], ans[i], ans[i + 1]);
    return res / 2;
}
Line get_line(Point a, Point b) {
    return Line((a + b) / 2, rotate(b - a, pi / 2));
}
double get_circle_line_intersection(Point a, Point b, Point &pa, Point &pb) { //圆与直线交点
    auto e = get_line_intersection(a, b - a, r, rotate(b - a, pi / 2));
    auto mind = get_dist(r, e);
    if (!on_segment(e, a, b)) mind = min(get_dist(r, a), get_dist(r, b));
    if (cmp(R, mind) <= 0) return mind;
    auto len = sqrt(R * R - get_dist(r, e) * get_dist(r, e));
    pa = e + norm(a - b) * len;
    pb = e + norm(b - a) * len;
    return mind;
}
double get_sector(Point a, Point b) { //扇形面积
    auto angle = acos(dot(a, b) / a.lenth() / b.lenth());
    if (sign(cross(a, b)) < 0) angle = -angle;
    return R * R * angle / 2;
}
double get_circle_triangle_area(Point a, Point b) { //三角剖分 圆与三角面积交
    auto da = get_dist(r, a), db = get_dist(r, b);
    if (cmp(R, da) >= 0 && cmp(R, db) >= 0) return cross(a, b) / 2;
    if (!sign(cross(a, b))) return 0;
    Point pa, pb;
    auto mind = get_circle_line_intersection(a, b, pa, pb);
    if (cmp(R, mind) <= 0) return get_sector(a, b);
    if (cmp(R, da) >= 0) return cross(a, pb) / 2 + get_sector(pb, b);
    if (cmp(R, db) >= 0) return get_sector(a, pa) + cross(pa, b) / 2;
    return get_sector(a, pa) + cross(pa, pb) / 2 + get_sector(pb, b);
}
Circle get_circle(Point a, Point b, Point c) { //三点求圆
    auto u = get_line(a, b), v = get_line(a, c);
    auto p = get_line_intersection(u.st, u.ed, v.st, v.ed);
    return Circle(p, get_dist(p, a));
}
Circle min_circle_cover() { //最小圆覆盖
    random_shuffle(a, a + n);
    Circle c = Circle(a[0], 0);
    for (int i = 1; i < n; i ++) {
        if (cmp(c.r, get_dist(c.p, a[i])) < 0) {
            c = {a[i], 0};
            for (int j = 0; j < i; j ++) {
                if (cmp(c.r, get_dist(c.p, a[j])) < 0) {
                    c = {(a[i] + a[j]) / 2, get_dist(a[i], a[j]) / 2};
                    for (int k = 0; k < j; k ++) {
                        if (cmp(c.r, get_dist(c.p, a[k])) < 0) {
                            c = get_circle(a[i], a[j], a[k]);
                        }
                    }
                }
            }
        }
    }
    return c;
}
double f(double x) { //辛普森积分函数(圆的面积并)
    int cnt = 0;
    for (int i = 0; i < n; i ++) {
        auto X = fabs(x - cir[i].p.x), R = cir[i].r;
        if (cmp(X, R) < 0) {
            auto Y = sqrt(R * R - X * X);
            a[cnt ++] = {cir[i].p.y - Y, cir[i].p.y + Y};
        }
    }
    if (!cnt) return 0;
    sort(a, a + cnt);
    double res = 0, st = a[0].x, ed = a[0].y;
    for (int i = 1; i < cnt; i ++) {
        if (a[i].x <= ed) {
            ed = max(ed, a[i].y);
        } else {
            res += ed - st;
            st = a[i].x, ed = a[i].y;
        }
    }
    return res + ed - st;
}
double simpson(double l, double r) { //辛普森积分
    auto mid = (l + r) / 2;
    return (r - l) * (f(l) + 4 * f(mid) + f(r)) / 6;
}
double asr(double l, double r, double s) { //自适应
    auto mid = (l + r) / 2;
    auto left = simpson(l, mid), right = simpson(mid, r);
    if (fabs(s - left - right) < eps) return left + right;
    return asr(l, mid, left) + asr(mid, r, right);
}
signed main() {
    
}

三维计算几何

#include <bits/stdc++.h>
using namespace std;
const int N = 1e4 + 5;
const double eps = 1e-12;
int n, m;
bool g[N][N];
double res;
double rand_eps() {
    return ((double)rand() / RAND_MAX - 0.5) * eps;
}
struct Point {
    double x, y, z;
    void shake() {
        x += rand_eps(), y += rand_eps(), z += rand_eps();
    }
    Point operator- (Point t) {
        return {x - t.x, y - t.y, z - t.z};
    }
    double operator& (Point t) {
        return x * t.x + y * t.y + z * t.z;
    }
    Point operator* (Point t) {
        return {y * t.z - t.y * z, z * t.x - x * t.z, x * t.y - y * t.x};
    }
    double len() {
        return sqrt(x * x + y * y + z * z);
    }
};
Point q[N];
struct Plane {
    int v[3];
    Point norm() { // 法向量
        return (q[v[1]] - q[v[0]]) * (q[v[2]] - q[v[0]]);
    }
    double area() { // 求面积
        return norm().len() / 2;
    }
    bool above(Point a) {
        return ((a - q[v[0]]) & norm()) >= 0;
    }
};
Plane plane[N], np[N];
void get_convex_3d() { //三维凸包
    plane[m ++] = {0, 1, 2};
    plane[m ++] = {2, 1, 0};
    for (int i = 3; i < n; i ++) {
        int cnt = 0;
        for (int j = 0; j < m; j ++) {
            bool t = plane[j].above(q[i]);
            if (!t) np[cnt ++] = plane[j];
            for (int k = 0; k < 3; k ++) g[plane[j].v[k]][plane[j].v[(k + 1) % 3]] = t;
        }
        for (int j = 0; j < m; j ++) {
            for (int k = 0; k < 3; k ++) {
                int a = plane[j].v[k], b = plane[j].v[(k + 1) % 3];
                if (g[a][b] && !g[b][a]) np[cnt ++] = {a, b, i};
            }
        }
        m = cnt;
        for (int j = 0; j < m; j ++) plane[j] = np[j];
    }
}
signed main() {
    cin >> n;
    for (int i = 0; i < n; i ++) {
        cin >> q[i].x >> q[i].y >> q[i].z;
        q[i].shake();
    }
    get_convex_3d();
    for (int i = 0; i < m; i ++) res += plane[i].area();
    printf("%lf", res);
}
凌乱之风
关注
  • 2
    点赞
  • 踩
  • 2
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
【计算几何全家桶】计算几何二维板子套路详解
SC_Linno的博客
11-16 758
title : 计算几何（二维） tags : ACM,数学,计算几何 date : 2021-10-23 author : Linno 向量 点乘与叉乘 叉乘一般用于求面积或者判断两点的相对位置：(x1,y1)×(x2,y2)=x1*y2-x2*y1; 点乘一般用于判断直线是否正交：(x1,y1)·(x2,y2)=x1*x2+y1*y2; inline double Len(Vt a){return sqrt(Dot(a,a));}//模长 inline Vt operator +(Vt a,Vt .
计算几何算法与应用(pdf)
04-19
这是一本关于计算几何的经典书籍。如果你想判断点是否在平面多边形内部，或者想知道离散点的三角网格划分的具体实现，从这本书中将得到答案。
二维几何基础大合集！《计算几何全家桶（一）》（基础运算、点、线、多边形、圆、网格）
繁凡さん的博客
09-27 6626
整理的算法模板合集： ACM模板 目录基本运算点与直线多边形 基本运算 #include<bits/stdc++.h> using namespace std; const int N = 5007, M = 50007, INF = 0x3f3f3f3f; const double DINF = 12345678910, eps = 1e-10; struct Point{ double x, y; Point(double x = 0, double y = 0):x(x
二维计算几何全家桶
最新发布
zcxxn的博客
02-09 557
计算几何——从入门到入土
多项式全家桶
olahiuj的博客
01-15 613
施工现场 多项式 在数学中，由若干个单项式相加组成的代数式叫做多项式（若有减法：减一个数等于加上它的相反数）。多项式中的每个单项式叫做多项式的项，这些单项式中的最高项次数，就是这个多项式的次数。其中多项式中不含字母的项叫做常数项。 为了方便接下来瞎逼逼，我们定义n次多项式A(x)=∑i=1nai⋅xiA(x)=\sum\limits_{i=1}^n a_i\cdot x^iA(x)=i=1∑n​...
三维几何基础大合集（三维点积叉积、点线面、凸包）《计算几何全家桶（三）》
繁凡さん的博客
10-03 2620
整理的算法模板合集： ACM模板 目录
计算几何模板
w750636248的专栏
07-26 685
#include #define MAXN 1000 #define offset 10000 #define eps 1e-8 #define PI acos(-1.0)//3.14159265358979323846 //判断一个数是否为0,是则返回true,否则返回false #define zero(x)(((x)>0?(x):-(x))<eps) //返回一个数的符号，正数返回1，负数返
繁凡的ACM算法全家桶（全新的模板整合计划）
热门推荐
繁凡さん的博客
01-21 2万+
整理的算法模板合集： ACM模板 点我看算法全家桶系列！！！ 实际上是一个全新的模板整合计划 目录一、多项式全家桶 算法全家桶！！！ 在搞完多项式全家桶之后就想着要不要把我学过的数论，组合数学，图论，数据结构，线性代数什么的都搞一个全家桶，实际上就是一个模板hhh，只不过是整合到一块的，打算把我之前搞的ACM模板全部再整理一下，完善一下，正好也全部复习一遍，把之前冗余的模板全部精炼一下。 大概需要10天的时间，明天先从数论全家桶开始 现在只有一个多项式全家桶撑场面hhh 一、多项式全家桶 点我看多项式
总结-计算几何
ghy1953
04-05 976
计算几何 做的题目很少, 而且模版也不熟. 所以还有这几天的时间, 把几个经典的算法弄熟弄懂, 模版能打出来就行了吧. 1. 凸包 2. 旋转卡壳 3. 半平面交 题目: 1. [codevs 1249] 多边形的面积 求多边形面积, 要理解叉积的意义. 2. [codevs 1298] 凸包周长 [codevs 3201] 奶牛代理商 XI 凸包周长 3. [codevs 1
【学习笔记】线性代数全家桶（在编程竞赛中的应用）
繁凡さん的博客
01-22 4710
两万字的算法竞赛线性代数全家桶！！！ OI / ACM均可食用，味美价廉
计算几何 - 你绝对找不到比这更好的计算几何
Antimonysbguy的博客
03-07 289
dzy哥哥回来给我们上课 杜老师好帅啊qwq 计算几何这些知识你需要自己推一遍向量的运算凸包半平面交 重新写计算几何 这些知识你需要自己推一遍 点 线段，及点在线段上的表示（两种向量法-一个叉乘一个等和线，一种距离法） 直线及其表示，用点斜法存储（学习斜率优化） struct node3{ double x,y,k; }line[NNN]; 直线引申为半平面 角度 double ang...
【夯实算法基础】 并查集
小新的博客
08-09 904
相信自己一定能学会
算法基础课贪心学习笔记
weixin_51368613的博客
04-07 818
一、区间问题 1.1 区间选点 自定义排序；resize；a[i].x > maxn 不是 a[i].y #include <iostream> #include <vector> #include <algorithm> #define x first #define y second using namespace std; typedef pair<int, int> PII; const int N = 1e5 + 7; int n, a,
最短路问题
m0_53203911的博客
05-08 320
最短路难在建图 单源最短路：一个起点 多源汇最短路：起点和终点不确定 n代表点数，m代表边数，m~n^2 叫稠密图，m~n叫稀疏图 求不超过k条边问题的话用bellman-ford算法 朴素dijkstra 算法，是稠密图，用邻接矩阵存 第一步表示初始化 第二步：以此枚举所有点x，从不在s中的所有点中找到离X的最近的点t，让t进入s，用t来更新其它点的距离 acwing板子题 #include<bits/stdc++.h> using namespace std; const int N=51
理解#define定义宏
learner_pu的博客
04-18 757
深入理解#define X（宏） 在写这种东西的时候,是否有困惑，它到底是什么东西？ #ifndef __TEST_H__ #define __TEST_H__ 首先回到#define的定义 #define A B : 将标识符A定义B的别名。 那么根据这个定义我们来看这个宏，它其实是说把定义的东西赋值成空，因为后面啥也没有，如果它的意思是定成空，就是一个宏，它存在在代码中也跟不存在一样。 下面的语句是可以编译的，因为它就是int a = 200; int __FACILITYOBJ_PM_H__
详解C语言中的#define、#undef、#indef、#ifndef、#else、#endif,#if,#elif
qq_46213710的博客
04-06 3741
广度优先搜索(BFS）---入门篇
qq_47982709的博客
04-06 694
填涂颜色 - 洛谷 这个题就是给我们一个图，图内有一个圈，我们将圈内的0 -> 2 即可； 我们可以将图内0看成两种，一种是圈内的0，一种是圈外的0； 有没有一种豁然开朗的感觉，我们考虑区分出这两种0就行了； 那显然是处理圈外的圆比较简单； 这个题唯一一点就是我们不能只搜一个点，因为可能会在边界存在一个突出的部分属于圈外的零； 又因为圈一定在图内，那我们就包围式搜索，从边界往里搜； #include<iostream> #include<cmath>...
2017CCPC 秦皇岛 L
ad4627的博客
09-12 117
water #include <bits/stdc++.h> #define pb push_back #define mp make_pair #define x first #define y second #define lson l,mid,rt<<1 #define rson mid+1,r,rt<<1|1 ...
#define几点注意事项
wangbinyantai的专栏
01-15 1395
1、 用无参宏定义一个简单的常量 #define LEN 12 这个是最常见的用法，但也会出错。 比如下面几个知识点你会吗？可以看下： （1）#define NAME "zhangyuncong" 程序中有"NAME"则，它会不会被替换呢？ （2）#define 0x abcd 可以吗？也就是说，可不可以用不是标识符的字母替
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值