第七个专题了,初期计算几何:
(1)、基本公式
拿了白书上面的三个例题做做。。。
1、uva 11178
题意:作三角形ABC每个内角的三等分线,相交成三角形DEF,则DEF是等边三角形。已知A,B,C三点坐标,问D,E,F三点坐标。
分析:简单的求直线交点、内角等分线可以通过直线旋转角度求出。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
using namespace std;
typedef long long LL;
const double EPS = 1e-10;
const int INF = 0x3f3f3f3f;
const int N = 310;
int dcmp(double x) {return abs(x) < EPS ? 0 : x < 0 ? 1 : -1;}
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline void in() { scanf("%lf %lf", &x, &y); }
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) {return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断一个点是否在一条线段上(不包含端点).
bool Onsegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q点到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
po solve(po a, po b, po c){
po v1 = rot(c - b, angle(a - b, c - b) / 3);
po v2 = rot(b - c, -angle(a - c, b - c) / 3);
return getlineintersection(b, v1, c, v2);
}
int main() {
int t;
scanf("%d", &t);
while (t--){
po a, b, c;
a.in(); b.in(); c.in();
po d = solve(a, b, c), e = solve(b, c, a), f = solve(c, a, b);
printf("%.6f %.6f %.6f %.6f %.6f %.6f\n", d.x, d.y, e.x, e.y, f.x, f.y);
}
return 0;
}
2、uvalive 3263
题意:给定一个一笔画构成的平面图的点、求出这一笔画把整个平面分成了几个部分。
分析:根据欧拉定理可以求出面数,假设一个平面图的点数为V,边数为E,则面数F = E + 2 – V。点数可以把所有直线的交点求出后去重得出、至于边数,可以枚举每一个点,如果一个点在一条边上,那么增加一条边。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
using namespace std;
const double EPS = 1e-10;
const int INF = 0x3f3f3f3f;
const int N = 310;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline void in() { scanf("%lf %lf", &x, &y); }
inline void out() { printf("%f %f", x, y); }
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) {return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
//bool on_seg(po p1, po p2, po q) { return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0; }
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
po p[N], v[N*N];
int main() {
int n, k = 0;
while (scanf("%d", &n), n) {
for (int i = 0; i < n; i++) p[i].in();
memcpy(v, p, sizeof(po) * n);
n--;
int e = n, c = n;
for (int i = 0; i < n; i++)
for (int j = i + 2; j < n; j++)
if (SegmentProperIntersection(p[i], p[i+1], p[j], p[j+1]))
v[c++] = intersection(p[i], p[i+1], p[j], p[j+1]);
sort(v, v + c);
c = unique(v, v + c) - v;
for (int i = 0; i < c; i++)
for (int j = 0; j < n; j++)
if (OnSegment(v[i], p[j], p[j+1])) e++;
printf("Case %d: There are %d pieces.\n", ++k, e + 2 - c);
}
return 0;
}
3、uva 11796
题意:两条狗匀速分别沿着折线跑,已知同时出发,同时到达,求奔跑过程中两条狗最大、最小距离的差值。
分析:直接求解比较麻烦、因为速度未知且时间未知而且奔跑的路线是一条折线。但是我们可以首先解决奔跑路线都是一条线段时的简单子问题,以第一条狗为参照物、这样就是求点到线段的距离了。而且可以知道两条狗的速度之比为各自的距离之比,时间对于答案无影响,所以可以假设速度分别为各自的总距离。这样从最起点开始,首先判断谁先到达下一个点,然后求出这个过程当中的最大最小距离,最大距离一定在第二条狗的路线的端点处求得,最小距离即为求点到线段的距离。每次一定有一条狗到达下一个点,那么整个问题就可以分解为A + B个简单的子问题了。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
using namespace std;
const double EPS = 1e-10;
const double INF = 0x3f3f3f3f;
const double PI = acos(-1.0);
const int N = 60;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) { return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {//线段相交判定
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po p, po a, po b) {//求出点p到线段ab之间的距离
if (a == b) return (p - a).length();
po v1 = b - a, v2 = p - a, v3 = p - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
double area(int n, po* p) {//求出任意多边形的有向面积。
double ans = 0;
for (int i = 0; i < n - 1; i++) ans += det(p[i] - p[0], p[i+1] - p[0]);
return ans / 2;
}
bool isconvexhull(po* p, int n) {//判断该多边形是否为凸包
int tag = dcmp(det(p[1] - p[0], p[2] - p[1]));
for (int i = 1; i < n; i++)
if (tag * dcmp(det(p[i+1] - p[i], p[i+2] - p[i+1])) < 0) return 0;
return 1;
}
bool ispointinpolygon(po* p, int n, po c) {//判断c点是否在多边形内
bool flag = 0;
for (int i = 0, j = n - 1; i < n; j = i++) {
if (((p[i].y > c.y) != (p[j].y > c.y)) && (c.x < (p[j].x - p[i].x) * (c.y - p[i].y) / (p[j].y - p[i].y) + p[i].x))
flag = !flag;
}
return flag;
}
vector<po> convexhull(po* p, int n) {//求出顶点数组构成的凸包
sort(p, p + n);
int k = 0;
vector<po> q(n << 1);
for (int i = 0; i < n; i++) {
while (k > 1 && det(q[k-1] - q[k-2], p[i] - q[k-1]) <= 0) k--;
q[k++] = p[i];
}
for (int i = n - 2, t = k; i >= 0; i--) {
while (k > t && det(q[k-1] - q[k-2], p[i] - q[k-1]) <= 0) k--;
q[k++] = p[i];
}
q.resize(k - 1);
return q;
}
double ma, mi;
void update(po p, po a, po b) {
mi = min(mi, DistanceToSegment(p, a, b));
ma = max(ma, (p - a).length());
ma = max(ma, (p - b).length());
}
po pa[N], pb[N];
int main() {
int t, ca = 0;
scanf("%d", &t);
while (t--) {
int a, b;
scanf("%d %d", &a, &b);
for (int i = 0; i < a; i++) scanf("%lf %lf", &pa[i].x, &pa[i].y);
for (int i = 0; i < b; i++) scanf("%lf %lf", &pb[i].x, &pb[i].y);
ma = -INF, mi = INF;
int f = 0, s = 0;
double lena = 0, lenb = 0;
for (int i = 1; i < a; i++) lena += (pa[i] - pa[i-1]).length();
for (int i = 1; i < b; i++) lenb += (pb[i] - pb[i-1]).length();
while (f < a - 1 && s < b - 1) {
double la = (pa[f+1] - pa[f]).length(), lb = (pb[s+1] - pb[s]).length();
double ti = min(la / lena, lb / lenb);
po va = (pa[f+1] - pa[f]) / la * lena * ti, vb = (pb[s+1] - pb[s]) / lb * lenb * ti;
update(pa[f], pb[s], pb[s] + vb - va);
pa[f] = pa[f] + va; pb[s] = pb[s] + vb;
if (pa[f] == pa[f+1]) f++;
if (pb[s] == pb[s+1]) s++;
}
printf("Case %d: %.0f\n", ++ca, ma - mi);
}
return 0;
}
(2)、点积和叉积的运用
1、poj2031
题意:给出三维坐标系上的一些球的球心坐标和其半径,搭建通路,使得他们能够相互连通。如果两个球有重叠的部分则算为已连通,无需再搭桥。求搭建通路的最小费用(费用就是边权,就是两个球面之间的距离)。
分析:就是使得所有点联通的最小代价、明显的最小生成树,如果两个球有公共部分,则他们之间的边权就为0.直接prim求解即可。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
using namespace std;
const int N = 110;
const double EPS = 1e-6;
const double INF = 0x3f3f3f3f;
double dis[N][N], d[N];
bool vis[N];
struct po { double x, y, z, r; } p[N];
void init(int n) {
fill(d + 1, d + n, INF);
memset(vis, 0, 4 * n);
}
double prim(int n) {
init(n+2);
double res = 0;
while (1) {
int v = -1;
for (int u = 0; u < n; u++)
if (!vis[u] && (v == -1 || d[u] < d[v])) v = u;
if (v == -1) break;
vis[v] = 1; res += d[v];
for (int u = 0; u < n; u++)
d[u] = min(d[u], dis[v][u]);
}
return res;
}
int main() {
int n;
while (scanf("%d", &n), n) {
for (int i = 0; i < n; i++) scanf("%lf %lf %lf %lf", &p[i].x, &p[i].y, &p[i].z, &p[i].r);
for (int i = 0; i < n; i++) {
for (int j = i + 1; j < n; j++) {
double d = sqrt((p[i].x - p[j].x) * (p[i].x - p[j].x) + (p[i].y - p[j].y) * (p[i].y - p[j].y) + (p[i].z - p[j].z) * (p[i].z - p[j].z));
d -= p[i].r + p[j].r;
dis[i][j] = dis[j][i] = d > EPS ? d : 0;
}
}
printf("%.3f\n", prim(n));
}
return 0;
}
2、poj1039
题意:有一宽度为1的折线管道,上面顶点为(xi,yi),所对应的下面顶点为(xi,yi-1),假设管道都是不透明的,不反射的,光线从左边入口处的(x1,y1),(x1,y1-1)之间射入,向四面八方传播,求解光线最远能传播到哪里(取x坐标)或者是否能穿透整个管道.
分析:光线最优时一定会经过一个上顶点和一个下顶点,否则总可以把光线进行旋转然后使得光线传播到更远的地方,这样,我们只需要枚举上顶点和下顶点,然后求出一束光线最远可以与哪根管道相交,只需要判断管道弯曲处的垂直x轴的线段,如果光线与该线段相交,那么则证明光线一定与弯曲处的管道的上壁或者下壁相交。每次取最大值即可。
#include <cstdio>#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
using namespace std;
const double EPS = 1e-10;
const double INF = 0x3f3f3f3f;
const int N = 30;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) {return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
//bool on_seg(po p1, po p2, po q) { return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0; }
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
double area(po* p, int n) {
double ans = 0;
for (int i = 0; i < n - 1; i++) ans += det(p[i] - p[0], p[i+1] - p[0]);
return ans / 2;
}
po p[N], p1[N];
int main() {
int n;
while (scanf("%d", &n), n) {
for (int i = 1; i <= n; i++) {
double x, y;
scanf("%lf %lf", &x, &y);
p[i] = po(x, y); p1[i] = po(x, y - 1);
};
double ans = -INF;
bool flag = 1;
for (int i = 1; i <= n && flag; i++) {
for (int j = 1; j <= n && flag; j++) {
if (i != j) {
po b; int k;
for (k = 1; k <= n; k++) {
b = intersection(p[i], p1[j], p1[k], p[k]);
if (dcmp(b.y - p1[k].y) < 0 || dcmp(b.y - p[k].y) > 0) break;
}
if (k == n + 1) flag = 0;
else if (k > max(i, j)) {
b = intersection(p[i], p1[j], p[k], p[k-1]);
if (dcmp(b.x - p[k-1].x) >= 0 && dcmp(b.x - p[k].x) <= 0) ans = max(ans, b.x);
b = intersection(p[i], p1[j], p1[k], p1[k-1]);
if (dcmp(b.x - p1[k-1].x) >= 0 && dcmp(b.x - p1[k].x) <= 0) ans = max(ans, b.x);
}
}
}
}
flag ? printf("%.2f\n", ans) : puts("Through all the pipe.");
}
return 0;
}
(3)、多边形的简单算法
1、 poj1408
题意:有一个1*1的正方形,分别给出下,上,左,右边每个边上的n个点,对边对应点连线,问这些线段相交的最大的四边形面积是多少。
分析:求多边形面积可以用简单的叉积实现,这样,算出所有的多边形面积取最大值就行了。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
using namespace std;
const double EPS = 1e-10;
const int INF = 0x3f3f3f3f;
const int N = 40;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline void in() { scanf("%lf %lf", &x, &y); }
inline void out() { printf("%f %f", x, y); }
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) {return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
//bool on_seg(po p1, po p2, po q) { return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0; }
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
double area(po a, po b, po c, po d) {
return (det(b - a, c - a) + det(c - a, d - a)) / 2;
}
po a[N], b[N], c[N], d[N], t[N][N];
int main() {
int n;
while (scanf("%d", &n), n) {
double in;
for (int i = 1; i <= n; i++) scanf("%lf", &in), a[i] = po(in, 0);
for (int i = 1; i <= n; i++) scanf("%lf", &in), b[i] = po(in, 1);
for (int i = 1; i <= n; i++) scanf("%lf", &in), c[i] = po(0, in);
for (int i = 1; i <= n; i++) scanf("%lf", &in), d[i] = po(1, in);
t[0][0] = po(0, 0); t[0][n+1] = po(1, 0);
for (int i = 1; i <= n; i++) t[0][i] = a[i];
t[n+1][0] = po(0, 1); t[n+1][n+1] = po(1, 1);
for (int i = 1; i <= n; i++) t[n+1][i] = b[i];
for (int i = 1; i <= n; i++) {
t[i][0] = c[i];
for (int j = 1; j <= n; j++) {
t[i][j] = intersection(c[i], d[i], a[j], b[j]);
}
t[i][n+1] = d[i];
}
double ans = 0;
for (int i = 0; i <= n; i++)
for (int j = 0; j <= n; j++) ans = max(ans, abs(area(t[i][j], t[i+1][j], t[i+1][j+1], t[i][j+1])));
printf("%.6f\n", ans);
}
return 0;
}
2、 poj1584
题意:按照顺时针或逆时针方向输入一个n边形的顶点坐标集,先判断这个n边形是否为凸包。再给定一个圆形(圆心坐标和半径),判断这个圆是否完全在n边形内部。
分析:要判断一个多边形是否为凸包、我们只需要从第一个点开始,依次判断第一条边向量与后面一条边向量的叉积,这个叉积可以决定方向是顺时针还是逆时针、这样,如果计算时出现顺时针与逆时针则证明不是凸包,这个自己手动模拟下即可。至于判断圆是否在多边形内,我们只需要类似判断三角形内圆的方式判断、也就是圆心到每条边的距离不小于半径、不过首先圆心必须在多边形内。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
using namespace std;
const double EPS = 1e-10;
const double INF = 0x3f3f3f3f;
const double PI = acos(-1.0);
const int N = 110;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) {return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
//bool on_seg(po p1, po p2, po q) { return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0; }
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
double area(po* p, int n) {
double ans = 0;
for (int i = 0; i < n - 1; i++) ans += det(p[i] - p[0], p[i+1] - p[0]);
return ans / 2;
}
po p[N], c;
int tag;
bool isconvexhull(int n) {
tag = dcmp(det(p[1] - p[0], p[2] - p[1]));
for (int i = 1; i < n; i++)
if (tag * dcmp(det(p[i+1] - p[i], p[i+2] - p[i+1])) < 0) return 0;
return 1;
}
bool ispointinpolygon(int n) {
bool flag = 0;
for (int i = 0, j = n - 1; i < n; j = i++) {
if (((p[i].y > c.y) != (p[j].y > c.y)) && (c.x < (p[j].x - p[i].x) * (c.y - p[i].y) / (p[j].y - p[i].y) + p[i].x))
flag = !flag;
}
return flag;
}
int main() {
int n;
double r;
while (scanf("%d", &n), n >= 3) {
scanf("%lf %lf %lf", &r, &c.x, &c.y);
for (int i = 0; i < n; i++) scanf("%lf %lf", &p[i].x, &p[i].y);
p[n] = p[0], p[n+1] = p[1];
bool flag = isconvexhull(n);
if (flag) {
flag &= ispointinpolygon(n);
for (int i = 0; i < n && flag; i++) {
double s = abs(det(p[i] - c, p[i+1] - c));
if (dcmp(s / (p[i+1] - p[i]).length() - r) == -1) flag = 0;
}
puts(flag ? "PEG WILL FIT" : "PEG WILL NOT FIT");
}
else puts("HOLE IS ILL-FORMED");
}
return 0;
}
(4)、凸包
1、 poj2187
题意:给定n个点,求最远点对。
分析:枚举效率太低,但是可以知道最远点对一定出现在这n个点的凸包中,对于这n个点的凸包,两两枚举即可。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
using namespace std;
const double EPS = 1e-10;
const int INF = 0x3f3f3f3f;
const int N = 50010;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline void in() { scanf("%lf %lf", &x, &y); }
inline void out() { printf("%f %f", x, y); }
inline double length() { return sqrt(x * x + y * y); }
};
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) {return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
//bool on_seg(po p1, po p2, po q) { return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0; }
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
double area(po a, po b, po c, po d) {
return (det(b - a, c - a) + det(c - a, d - a)) / 2;
}
po p[N];
vector<po> convex_hull(po* p, int n) {
sort(p, p + n);
int k = 0;
vector<po> q(n << 1);
for (int i = 0; i < n; i++) {
while (k > 1 && det(q[k-1] - q[k-2], p[i] - q[k-1]) <= 0) k--;
q[k++] = p[i];
}
for (int i = n - 2, t = k; i >= 0; i--) {
while (k > t && det(q[k-1] - q[k-2], p[i] - q[k-1]) <= 0) k--;
q[k++] = p[i];
}
q.resize(k - 1);
return q;
}
int main() {
int n;
scanf("%d", &n);
for (int i = 0; i < n; i++) p[i].in();
vector<po> q = convex_hull(p, n);
double res = 0;
for (int i = 0; i < q.size(); i++)
for (int j = 0; j < i; j++) res = max(res, (q[i] - q[j]).length());
printf("%.0f\n", res * res);
return 0;
}
2、 poj1113
题意:给一个多边形、求一个圈包围整个多边形且这个圈离多边形的距离不能小于给定的L,求最小的圈的长度。
分析:可以知道如果圈离最外围的点的距离大于L,离里面的点的距离也会大于L,那么只需要求出包含这个多边形的凸包的圈的最小长度,手动画图加上猜想与推理可以知道圈有一些半径为L的圆的圆弧与一些线段组成且长度就为一个完整的圆的周长加上凸包的周长。
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
using namespace std;
const double EPS = 1e-10;
const double INF = 0x3f3f3f3f;
const double PI = acos(-1.0);
const int N = 110;
int dcmp(double x) { return abs(x) < EPS ? 0 : x < 0 ? -1 : 1; }
struct po {
double x, y;
po(double x = 0, double y = 0) : x(x), y(y) {}
inline double length() { return sqrt(x * x + y * y); }
};
po p[N];
inline po operator + (po a, po b) { return po(b.x + a.x, b.y + a.y); }
inline po operator - (po a, po b) { return po(a.x - b.x, a.y - b.y); }
inline po operator * (po a, double d) { return po(a.x * d, a.y * d); }
inline po operator / (po a, double d) { return po(a.x / d, a.y / d); }
inline bool operator < (const po &a, const po &b) { return a.x == b.x ? a.y < b.y : a.x < b.x; }
inline bool operator == (const po &a, const po &b) { return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0; }
inline double dot(po a, po b) { return a.x * b.x + a.y * b.y; }//点积
inline double det(po a, po b) { return a.x * b.y - a.y * b.x; }//叉积
//求出以向量a绕起点逆时针旋转a角度之后该向量的坐标
inline po rot(po a, double rad) { return po(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad)); }
double angle(po a, po b) { return acos(dot(a, b) / a.length() / b.length()); }//求向量a和b之间的角度
double area(po a, po b, po c) { return det(b - a, c - a); }//求出以ab、ac为邻边组成的三角形的面积
inline po normal(po a) {
double l = a.length();
return po(-a.y / l, a.x / l);
}
//求出直线p1+tv和直线q+tw的交点,v和w为方向向量
po getlineintersection(po p, po v, po q, po w) {
po u = p - q;
double t = det(w, u) / det(v, w);
return p + v * t;
}
//求出p1、p2点确定的直线与q1、q2点确定的直线的交点,使用确保两条直线不会平行
po intersection(po p1, po p2, po q1, po q2) { return p1 + (p2 - p1) * det(q2 - q1, q1 - p1) / det(q2 - q1, p2 - p1); }
//线段相交判定
bool SegmentProperIntersection(po a1, po a2, po b1, po b2) {
double c1 = det(a2 - a1, b1 - a1), c2 = det(a2 - a1, b2 - a1),
c3 = det(b2 - b1, a1 - b1), c4 = det(b2 - b1, a2 - b1);
return dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0;
}
//判断q点是否在p1、p2点组成的线段上
//bool on_seg(po p1, po p2, po q) { return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0; }
bool OnSegment(po p, po a1, po a2) { return dcmp(det(a1 - p, a2 - p)) == 0 && dcmp(dot(a1 - p, a2 - p)) < 0; }
double DistanceToLine(po q, po a, po b) {//求出q到直线ab的距离。
po v1 = b - a, v2 = q - a;
return abs(det(v1, v2)) / v1.length();
}
double DistanceToSegment(po q, po a, po b) {//求出点q到线段ab之间的距离
if (a == b) return (q - a).length();
po v1 = b - a, v2 = q - a, v3 = q - b;
if (dcmp(dot(v1, v2)) < 0) return v2.length();
if (dcmp(dot(v1, v3)) > 0) return v3.length();
return abs(det(v1, v2)) / v1.length();
}
po GetLinePorjection(po q, po a, po b) {//获得点q在直线ab上的投影
po v = b - a;
return a + v * dot(v, q - a) / dot(v, v);
}
//求出任意多边形的有向面积。
double area(int n) {
double ans = 0;
for (int i = 0; i < n - 1; i++) ans += det(p[i] - p[0], p[i+1] - p[0]);
return ans / 2;
}
//判断该多边形是否为凸包
bool isconvexhull(int n) {
int tag = dcmp(det(p[1] - p[0], p[2] - p[1]));
for (int i = 1; i < n; i++)
if (tag * dcmp(det(p[i+1] - p[i], p[i+2] - p[i+1])) < 0) return 0;
return 1;
}
//判断c点是否在多边形内
bool ispointinpolygon(int n, po c) {
bool flag = 0;
for (int i = 0, j = n - 1; i < n; j = i++) {
if (((p[i].y > c.y) != (p[j].y > c.y)) && (c.x < (p[j].x - p[i].x) * (c.y - p[i].y) / (p[j].y - p[i].y) + p[i].x))
flag = !flag;
}
return flag;
}
//求出顶点数组构成的凸包
vector<po> convexhull(int n) {
sort(p, p + n);
int k = 0;
vector<po> q(n << 1);
for (int i = 0; i < n; i++) {
while (k > 1 && det(q[k-1] - q[k-2], p[i] - q[k-1]) <= 0) k--;
q[k++] = p[i];
}
for (int i = n - 2, t = k; i >= 0; i--) {
while (k > t && det(q[k-1] - q[k-2], p[i] - q[k-1]) <= 0) k--;
q[k++] = p[i];
}
q.resize(k - 1);
return q;
}
int main() {
int n, l;
while (~scanf("%d %d", &n, &l)) {
for (int i = 0; i < n; i++) scanf("%lf %lf", &p[i].x, &p[i].y);
vector<po> t = convexhull(n);
int len = t.size() - 1;
double ans = 0;
for (int i = 0; i < len; i++) ans += (t[i+1] - t[i]).length();
ans += (t[len] - t[0]).length();
printf("%.0f\n", ans + 2 * PI * l);
}
return 0;
}