计算几何
基础知识
1.acos是cos的反函数,acos(-1) = pi
2.正弦定理和余弦定理
3.比较两个数的大小时,用eps规避误差
int cmp(db x, db y)
{
if(fabs(x - y) < eps) return 0;
else if(x < y) return -1;
return 1;
}
4.向量 坐标的相加减
5.内积(点积):一个向量在另外一个向量上投影的长度乘上这个向量的长度。
double dot(node a, node b)
{
return a.x * b.x + a.y * b.y;
}
6.外积(叉积):由向量A和向量B构成的三角形的面积的两倍。
double cross(node a, node b)
{
return a.x * b.y - b.x * a.y;
}
逆时针的时候正数,顺时针的时候是负数。(我猜测是因为三角函数的图像中sin的角度是按照逆时针方向来计算的)
用右手螺旋定则,右手指向起点向量的方向,朝中点防线旋转,朝向外面的就是正,朝向里面的就是负。
当然,一般使用的会直接取绝对值
7.计算向量的夹角
double get_angle(node a, node b)
{
return acos(dot(a, b) / get_length(a) / get_length(b));
}
8.计算由三个已知点构成的唯一一个平行四边形的面积。
double area(node a, node b, node c)
{
return cross(b - a, c - a);
}
9.求一个向量顺时针旋转某个角度所得到的向量。
node rotate(node a, double angle)//顺时针旋转angle的角度
{
return (node){a.x * cos(angle) + a.y * sin(angle), a.y * cos(angle) - a.x * sin(angle); }
}
10.直线的表示
ax + by + c = 0
y = kx + b
a + tv
11.如何判断某点是不是在某条直线上面:叉积结果是0的话就在上面,否则就不在,可以画一个图来证明。
node get_line_intersection(node p, vector v, node q, vector w)
{
if(cross(v, w) == 0) return (node){inf, inf};
vector u = p - q;
double t = cross(u, w) / cross(v, w);
return p + t * v;
}
12.点到直线的距离
该点所对应的向量与直线所对应的向量的叉积除以该直线对应向量的模长。
13.点到线段的距离需要分别判断一下线段两端端点 和该点 和该线段 形成的叫是否都是锐角。如果是钝角直接返回该点到那个端点的长度就可以了,否则直接求点到直线的距离。
14.点在某条直线上面的投影。
15.点是否在线段上。
向量(a, p) 和向量(b, p)的叉积等于0 并且 点积小于等于0
16.判断两条线的线段是否相交。直接看两条线的两个点是否互相跨立在对方的两侧。
直接用叉积判断 就是一条直线的一个端点到另外一条直线两个端点的向量的和该直线向量的叉积之积 <= 0代表端点相交也算 < 0不算端点。
「POJ2318」TOYS
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 5005;
int n, m, ans[N], a[N], b[N];
struct node
{
int x, y;
};
int x_1, y_1, x_2, y_2;
ll cross(node a, node b)
{
return (ll)a.x * b.y - (ll)a.y * b.x;
}
ll area(node a, node b, node c)
{
return cross((node){b.x - a.x, b.y - a.y}, (node){c.x - a.x, c.y - a.y});
}
int find_pos(node p)
{
int l = 0, r = n - 1, mid;
while(l <= r)
{
mid = l + r >> 1;
if(area(p, (node){b[mid], y_2}, (node){a[mid], y_1}) > 0) r = mid - 1;
else l = mid + 1;
}
return r + 1;
}
int main()
{
bool pd = true;
while(scanf("%d", &n) != EOF && n)
{
scanf("%d %d %d %d %d", &m, &x_1, &y_1, &x_2, &y_2);
for (int i = 0; i < n; ++ i)
{
ans[i] = 0;
scanf("%d %d", &a[i], &b[i]);
}
ans[n] = 0;
for(int i = 0; i < m; ++ i)
{
node p;
scanf("%d %d", &p.x, &p.y);
ans[find_pos(p)] ++;
}
if(pd == true) pd = false;
else printf("\n");
for (int i = 0; i <= n; ++ i) printf("%d: %d\n", i, ans[i]);
}
return 0;
}
[POJ3304]线段
#include <bits/stdc++.h>
using namespace std;
const int N = 405;
typedef double db;
const db eps = 1e-8;
int n, cnt = 0;
struct node
{
db x, y;
}a[N], b[N], p[N];
int cmp(db u, db v)
{
if(fabs(u - v) < eps) return 0;
else if(u < v) return -1;
return 1;
}
db cross(db x_1, db y_1, db x_2, db y_2)
{
return x_1 * y_2 - x_2 * y_1;
}
db area(node a, node b, node c)
{
return cross(b.x - a.x, b.y - a.y, c.x - a.x, c.y - a.y);
}
bool check()
{
for (int i = 1; i <= cnt; ++ i)
{
for (int j = 1; j <= cnt; ++ j)
{
if(!cmp(p[i].x, p[j].x) && !cmp(p[i].y, p[j].y)) continue;
bool pd = true;
for (int k = 1; k <= n; ++ k)
{
if(cmp(area(p[i], p[j], a[k]), 0) * cmp(area(p[i], p[j], b[k]), 0) > 0) {pd = false; break; }
}
if(pd) return true;
}
}
return false;
}
int main()
{
int T;
scanf("%d", &T);
while(T --)
{
scanf("%d", &n);
cnt = 0;
for (int i = 1; i <= n; ++ i)
{
scanf("%lf %lf %lf %lf", &a[i].x, &a[i].y, &b[i].x, &b[i].y);
p[++ cnt] = a[i], p[++ cnt] = b[i];
}
if(check()) printf("Yes!\n");
else printf("No!\n");
}
return 0;
}
凸包
平面上给我们很多点,用一个橡皮筋把所有点围住,会形成一个凸多边形,这个凸多边形就成为这个图形的凸包。
Andrew算法:
将第一关键字从小到大排序,第二关键字随便排序。
先从左到右做上半部分,再从右到左做下半部分。
然后再用单调栈维护,如果该点是在之前向量的顺时针方向(叉乘为正数),那么栈的当前点需要删除,否则的话直接加进去。至于在同一条直线上面的点到底删不删除,是要根据题目进行区分的。
#include <bits/stdc++.h>
using namespace std;
const int N = 10005;
typedef double db;
const db eps = 1e-8;
bool used[N];
struct node
{
db x, y;
friend node operator - (node a, node b) {return (node){a.x - b.x, a.y - b.y}; }
}p[N];
int stk[N], top = 0, n;
bool sort_cmp(node a, node b)
{
if(a.x == b.x) return a.y > b.y;
return a.x < b.x;
}
db get_dis(node a, node b)
{
db dx = a.x - b.x;
db dy = a.y - b.y;
return sqrt(dx * dx + dy * dy);
}
int cmp(db a, db b)
{
if(fabs(a - b) < eps) return 0;
else if(a < b) return -1;
return 1;
}
db cross(node a, node b)
{
return a.x * b.y - a.y * b.x;
}
db area(node a, node b, node c)
{
return cross(b - a, c - a);
}
db Andrew()
{
sort(p + 1, p + n + 1, sort_cmp);
stk[++ top] = 1;
for (int i = 1; i <= n; ++ i)
{
while(top >= 2 && area(p[stk[top - 1]], p[stk[top]], p[i]) >= 0)
{
used[stk[top]] = 0;
top --;
}
stk[++ top] = i, used[i] = 1;
}
used[1] = 0;
for (int i = n; i >= 1; -- i)
{
if(used[i]) continue;
while(top >= 2 && area(p[stk[top - 1]], p[stk[top]], p[i]) >= 0) top --;
stk[++ top] = i;
}
db res = 0;
for (int i = 1; i <= top - 1; ++ i) res += get_dis(p[stk[i]], p[stk[i + 1]]);
return res;
}
int main()
{
scanf("%d", &n);
for (int i = 1; i <= n; ++ i) scanf("%lf %lf", &p[i].x, &p[i].y);
printf("%.2lf", Andrew());
}
「BZOJ2829」信用卡凸包
调了很久,因为周长公式写成面积公式了。
#include <bits/stdc++.h>
using namespace std;
const int N = 40005;
typedef double db;
const db eps = 1e-8;
int n, cnt = 0, top = 0, s[N];
db a, b, r, res = 0;
bool used[N];
struct node
{
db x, y;
friend node operator - (node a, node b) {return (node){a.x - b.x, a.y - b.y}; }
}p[N];
node turn(db x, db y, db nx, db ny, db theta)
{
return (node){(nx - x) * cos(theta) - (ny - y) * sin(theta) + x, (ny - y) * cos(theta) + (nx - x) * sin(theta) + y};
}
int cmp(db a, db b)
{
if(fabs(a - b) < eps) return 0;
else if(a < b) return -1;
return 1;
}
bool sort_cmp(node a, node b)
{
if(cmp(a.x, b.x) == 0) return a.y < b.y;
return a.x < b.x;
}
db cross(node a, node b)
{
return a.x * b.y - a.y * b.x;
}
db area(node a, node b, node c)
{
return cross(b - a, c - a);
}
db get_dis(node a, node b)
{
return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
db andrew()
{
sort(p + 1, p + cnt + 1, sort_cmp);
for (int i = 1; i <= cnt; ++ i)
{
while(top >= 2 && area(p[s[top - 1]], p[s[top]], p[i]) >= 0)
{
used[s[top]] = 0, top --;
}
s[++ top] = i, used[i] = 1;
}
used[1] = 0;
for (int i = cnt; i >= 1; -- i)
{
if(used[i]) continue;
while(top >= 2 && area(p[s[top - 1]], p[s[top]], p[i]) >= 0) top --;
s[++ top] = i;
}
for (int i = 1; i <= top - 1; ++ i) res += get_dis(p[s[i]], p[s[i + 1]]);
}
int main()
{
scanf("%d", &n);
scanf("%lf %lf %lf", &a, &b, &r);
for (int i = 1; i <= n; ++ i)
{
db x, y, theta;
scanf("%lf %lf %lf", &x, &y, &theta);
db nx, ny;
nx = x + b / 2 - r, ny = y + a / 2 - r;
p[++ cnt] = turn(x, y, nx, ny, theta);
nx = x + b / 2 - r, ny = y - a / 2 + r;
p[++ cnt] = turn(x, y, nx, ny, theta);
nx = x - b / 2 + r, ny = y + a / 2 - r;
p[++ cnt] = turn(x, y, nx, ny, theta);
nx = x - b / 2 + r, ny = y - a / 2 + r;
p[++ cnt] = turn(x, y, nx, ny, theta);
}
res = acos(-1) * r * 2;
andrew();
printf("%.2lf", res);
return 0;
}
今日小trick。
所有函数记得写return,不然没有返回值的话容易RE。
所有==可以把常数写在前面。
半平面交
「模板」凸多边形/半平面交模板题
#include <bits/stdc++.h>
using namespace std;
const int N = 505;
typedef double db;
const db eps = 1e-8;
int n, m, cnt = 0;
struct node
{
db x, y;
friend node operator - (node a, node b) {return (node){a.x - b.x, a.y - b.y}; }
}pg[N];
struct line
{
node st, ed;
}p[N];
db get_angle(line u)
{
return atan2(u.ed.y - u.st.y, u.ed.x - u.st.x);
}
int cmp(db a, db b)
{
if(fabs(a - b) < eps) return 0;
else if(a < b) return -1;
return 1;
}
db cross(node a, node b)
{
return a.x * b.y - a.y * b.x;
}
db area(node a, node b, node c)
{
return cross(b - a, c - a);
}
bool sort_cmp(line u, line v)
{
db gu = get_angle(u), gv = get_angle(v);
if(0 == cmp(gu, gv)) return area(u.st, u.ed, v.ed) < 0;
else return gu < gv;
}
int tt = 1, hh = 0, q[N];
node get_intersection_(node p, node v, node q, node w)
{
node u = p - q;
db t = cross(w, u) / cross(v, w);
return (node){p.x + t * v.x, p.y + t * v.y};
}
node get_intersection(line a, line b)
{
return get_intersection_(a.st, a.ed - a.st, b.st, b.ed - b.st);
}
bool on_right(line a, line b, line c)
{
node o = get_intersection(b, c);
return area(a.st, a.ed, o) <= 0;
}
db get_half_intersection()
{
sort(p + 1, p + cnt + 1, sort_cmp);
tt = 1, hh = 0;
for (int i = 1; i <= cnt; ++ i)
{
if(1 != i && 0 == cmp(get_angle(p[i]), get_angle(p[i - 1]))) continue;
while(tt + 1 <= hh && on_right(p[i], p[q[hh - 1]], p[q[hh]])) hh --;
while(tt + 1 <= hh && on_right(p[i], p[q[tt]], p[q[tt + 1]])) tt ++;
q[++ hh] = i;
}
while(tt <= hh + 1 && on_right(p[q[tt]], p[q[hh - 1]], p[q[hh]])) hh --;
// while(tt <= hh + 1 && on_right(p[q[hh]], p[q[tt]], p[q[tt + 1]])) tt ++;
q[++ hh] = q[tt];
int k = 0;
for (int i = tt; i <= hh - 1; ++ i) pg[k ++] = get_intersection(p[q[i]], p[q[i + 1]]);
db res = 0.0;
for (int i = 1; i < k - 1; ++ i)
{
res = res + area(pg[0], pg[i], pg[i + 1]);
}
res = res / 2.0;
return res;
}
int main()
{
scanf("%d", &n);
for (int i = 1; i <= n; ++ i)
{
scanf("%d", &m);
for (int j = 1; j <= m; ++ j) scanf("%lf %lf", &pg[j].x, &pg[j].y);
for (int j = 1; j <= m - 1; ++ j) p[++ cnt].st = pg[j], p[cnt].ed = pg[j + 1];
p[++ cnt].st = pg[m], p[cnt].ed = pg[1];
}
printf("%.3lf", get_half_intersection());
return 0;
}
JLOI2013赛车
#include <bits/stdc++.h>
using namespace std;
const int N = 10005;
typedef long double db;
typedef long long ll;
const db eps = 1e-18;
int n, cnt = 0, ans[N];
struct node
{
db x, y;
friend node operator - (node a, node b) {return (node){a.x - b.x, a.y - b.y}; }
};
struct Line
{
node st, ed;
vector<int> id;
}line[N];
struct car
{
int k, b;
bool operator < (const car &a) const
{
if(k == a.k) return b < a.b;
return k < a.k;
}
}cr[N];
map<car, int> mp;
db get_angle(Line a)
{
return atan2(a.ed.y - a.st.y, a.ed.x - a.st.x);
}
int cmp(db a, db b)
{
if(fabs(a - b) < eps) return 0;
else if(a < b) return -1;
return 1;
}
db cross(node a, node b)
{
return a.x * b.y - a.y * b.x;
}
db area(node a, node b, node c)
{
return cross(b - a, c - a);
}
bool sort_cmp(Line a, Line b)
{
db ga = get_angle(a), gb = get_angle(b);
if(cmp(ga, gb) == 0) return area(a.st, a.ed, b.st) < 0;
return ga < gb;
}
int tt = 1, hh = 0, q[N];
node get_intersection_(node p, node v, node q, node w)
{
node u = p - q;
db t = cross(w, u) / cross(v, w);
return (node){p.x + t * v.x, p.y + t * v.y};
}
node get_intersection(Line a, Line b)
{
return get_intersection_(a.st, a.ed - a.st, b.st, b.ed - b.st);
}
bool on_right(Line a, Line b, Line c)
{
node o = get_intersection(b, c);
return cmp(area(a.st, a.ed, o), 0) == -1;
}
void get_half_intersection()
{
sort(line + 1, line + cnt + 1, sort_cmp);
for (int i = 1; i <= cnt; ++ i)
{
if(i != 1 && cmp(get_angle(line[i]), get_angle(line[i - 1])) == 0) continue;
while(tt + 1 <= hh && on_right(line[i], line[q[hh - 1]], line[q[hh]])) hh --;
while(tt + 1 <= hh && on_right(line[i], line[q[tt]], line[q[tt + 1]])) tt ++;
q[++ hh] = i;
}
while(tt + 1 <= hh && on_right(line[q[tt]], line[q[hh - 1]], line[q[hh]])) hh --;
while(tt + 1 <= hh && on_right(line[q[hh]], line[q[tt]], line[q[tt + 1]])) tt ++;
int k = 0;
for (int i = tt; i <= hh; ++ i)
{
int pos = q[i];
for (int j = 0; j < line[pos].id.size(); ++ j) ans[++ k] = line[pos].id[j];
}
sort(ans + 1, ans + k + 1);
printf("%d\n", k);
for (int i = 1; i <= k; ++ i) printf("%d ", ans[i]);
return ;
}
int main()
{
scanf("%d", &n);
ll s = 0;
for (int i = 1; i <= n; ++ i) scanf("%d", &cr[i].b), s += cr[i].b;
for (int i = 1; i <= n; ++ i) scanf("%d", &cr[i].k);
if(s == 0)
{
printf("%d\n", n);
for (int i = 1; i <= n; ++ i) printf("%d ", i);
return 0;
}
line[++ cnt] = (Line){(node){0, 1}, (node){0, 0}};
line[++ cnt] = (Line){(node){0, 0}, (node){1, 0}};
for (int i = 1; i <= n; ++ i)
{
if(!mp.count(cr[i])) mp[cr[i]] = ++ cnt, line[cnt] = (Line){(node){0.0, cr[i].b * 1.0}, (node){1.0, 1.0 * cr[i].b + 1.0 * cr[i].k}};
int pos = mp[cr[i]];
line[pos].id.push_back(i);
}
get_half_intersection();
return 0;
}

浙公网安备 33010602011771号