计算几何习题
Q1(Problem source : poj 2624):
Description
Input
Output
#include<cstdio> #include<cstring> struct point { double x , y; }p[4]; struct vec { double x , y; }v[3]; int main() { while(scanf("%lf%lf%lf%lf%lf%lf%lf%lf",&p[0].x,&p[0].y,&p[1].x,&p[1].y,&p[2].x,&p[2].y,&p[3].x,&p[3].y) != EOF) { double x , y; x = p[0].x + p[1].x + p[2].x + p[3].x; y = p[0].y + p[1].y + p[2].y + p[3].y; point temp; int flag; if((p[0].x == p[2].x && p[0].y == p[2].y) || (p[0].x == p[3].x && p[0].y == p[3].y)) flag = 0; else flag = 1; x -= 3*p[flag].x; y -= 3*p[flag].y; printf("%.3lf %.3lf\n",x , y); } }
从这道题目不难看出,在进行模拟运算之前最好人为得将问题进行一定程度的简化,尽量缩减编程模拟的过程。
关于半平面求交:
Q1(Problem source : poj 1474):
题目大意:给出几何图形的n个顶点,让你判断是否存在某个区域,该区域的点可以看到几何图形边上的所有点。
分析:很典型的利用半平面求交判断几何图形是否存在内核区域,半平面求交的过程中的面显然是n个顶点的几何图形,而直线则是遍历该几何图形的n条边,注意遍历过程中的方向一定要保持一致。
参考代码如下:
#include <iostream> #include <cstdio> #include <algorithm> #include <cmath> using namespace std; const double eps = 1e-8; const int maxn = 105; int dq[maxn], top, bot, pn, order[maxn], ln; struct Point { double x, y; } p[maxn]; struct Line { Point a, b; double angle; } l[maxn]; int dblcmp(double k) { if (fabs(k) < eps) return 0; return k > 0 ? 1 : -1; } double multi(Point p0, Point p1, Point p2) { return (p1.x-p0.x)*(p2.y-p0.y)-(p1.y-p0.y)*(p2.x-p0.x); } bool cmp(int u, int v) { int d = dblcmp(l[u].angle-l[v].angle); if (!d) return dblcmp(multi(l[u].a, l[v].a, l[v].b)) < 0; return d < 0; } void getIntersect(Line l1, Line l2, Point& p) { double dot1,dot2; dot1 = multi(l2.a, l1.b, l1.a); dot2 = multi(l1.b, l2.b, l1.a); p.x = (l2.a.x * dot2 + l2.b.x * dot1) / (dot2 + dot1); p.y = (l2.a.y * dot2 + l2.b.y * dot1) / (dot2 + dot1); } bool judge(Line l0, Line l1, Line l2) { Point p; getIntersect(l1, l2, p); return dblcmp(multi(p, l0.a, l0.b)) > 0; } void addLine(double x1, double y1, double x2, double y2) { l[ln].a.x = x1; l[ln].a.y = y1; l[ln].b.x = x2; l[ln].b.y = y2; l[ln].angle = atan2(y2-y1, x2-x1); order[ln] = ln; ln++; } void halfPlaneIntersection() { int i, j; sort(order, order+ln, cmp); for (i = 1, j = 0; i < ln; i++) if (dblcmp(l[order[i]].angle-l[order[j]].angle) > 0) order[++j] = order[i]; ln = j + 1; dq[0] = order[0]; dq[1] = order[1]; bot = 0; top = 1; for (i = 2; i < ln; i++) { while (bot < top && judge(l[order[i]], l[dq[top-1]], l[dq[top]])) top--; while (bot < top && judge(l[order[i]], l[dq[bot+1]], l[dq[bot]])) bot++; dq[++top] = order[i]; } while (bot < top && judge(l[dq[bot]], l[dq[top-1]], l[dq[top]])) top--; while (bot < top && judge(l[dq[top]], l[dq[bot+1]], l[dq[bot]])) bot++; } bool isThereACore() { if (top-bot > 1) return true; return false; } int main() { int t, i; int tt = 1; while (scanf ("%d", &pn) && pn) { for (i = 0; i < pn; i++) scanf ("%lf%lf", &p[i].x, &p[i].y); for (ln = i = 0; i < pn-1; i++) addLine(p[i].x, p[i].y, p[i+1].x, p[i+1].y); addLine(p[i].x, p[i].y, p[0].x, p[0].y); halfPlaneIntersection(); if (isThereACore()) printf ("Floor #%d\nSurveillance is possible.\n\n",tt++); else printf ("Floor #%d\nSurveillance is impossible.\n\n",tt++); } return 0; }
Q2(Problem source : poj 1279):
题目大意:给出平面几何图形的n个顶点,数据保证该图形有内核,求解内核的面积。
#include <iostream> #include <cstdio> #include <cmath> #include <algorithm> using namespace std; const int maxn = 1505; const double eps = 1e-8; int n, pn, dq[maxn], top, bot; struct Point { double x, y; } p[maxn]; struct Line { Point a, b; double angle; Line& operator= (Line l) { a.x = l.a.x; a.y = l.a.y; b.x = l.b.x; b.y = l.b.y; angle = l.angle; return *this; } } l[maxn]; int dblcmp(double k) { if (fabs(k) < eps) return 0; return k > 0 ? 1 : -1; } double multi(Point p0, Point p1, Point p2) { return (p1.x-p0.x)*(p2.y-p0.y)-(p1.y-p0.y)*(p2.x-p0.x); } bool cmp(const Line& l1, const Line& l2) { int d = dblcmp(l1.angle-l2.angle); if (!d) return dblcmp(multi(l1.a, l2.a, l2.b)) < 0; //大于0取半平面的左半,小于0取右半 return d < 0; } void addLine(Line& l, double x1, double y1, double x2, double y2) { l.a.x = x1; l.a.y = y1; l.b.x = x2; l.b.y = y2; l.angle = atan2(y2-y1, x2-x1); } void getIntersect(Line l1, Line l2, Point& p) { double A1 = l1.b.y - l1.a.y; double B1 = l1.a.x - l1.b.x; double C1 = (l1.b.x - l1.a.x) * l1.a.y - (l1.b.y - l1.a.y) * l1.a.x; double A2 = l2.b.y - l2.a.y; double B2 = l2.a.x - l2.b.x; double C2 = (l2.b.x - l2.a.x) * l2.a.y - (l2.b.y - l2.a.y) * l2.a.x; p.x = (C2 * B1 - C1 * B2) / (A1 * B2 - A2 * B1); p.y = (C1 * A2 - C2 * A1) / (A1 * B2 - A2 * B1); } bool judge(Line l0, Line l1, Line l2) { Point p; getIntersect(l1, l2, p); return dblcmp(multi(p, l0.a, l0.b)) > 0; //与上面的注释处的大于小于符号相反,大于0,是p在向量l0.a->l0.b的左边,小于0是在右边,当p不在半平面l0内时,返回true } void HalfPlaneIntersect( ){ int i, j; sort(l, l+n, cmp); for (i = 0, j = 0; i < n; i++) if (dblcmp(l[i].angle-l[j].angle) > 0) l[++j] = l[i]; n = j + 1; dq[0] = 0; dq[1] = 1; top = 1; bot = 0; for (i = 2; i < n; i++) { while (top > bot && judge(l[i], l[dq[top]], l[dq[top-1]])) top--; while (top > bot && judge(l[i], l[dq[bot]], l[dq[bot+1]])) bot++; dq[++top] = i; } while (top > bot && judge(l[dq[bot]], l[dq[top]], l[dq[top-1]])) top--; while (top > bot && judge(l[dq[top]], l[dq[bot]], l[dq[bot+1]])) bot++; dq[++top] = dq[bot]; for (pn = 0, i = bot; i < top; i++, pn++) getIntersect(l[dq[i+1]], l[dq[i]], p[pn]); } double getArea() { if (pn < 3) return 0; double area = 0; for (int i = 1; i < pn-1; i++) area += multi(p[0], p[i], p[i+1]); if (area < 0) area = -area; return area/2; } int main() { int t, i; scanf ("%d", &t); while (t--) { scanf ("%d", &n); for (i = 0; i < n; i++) scanf ("%lf%lf", &p[i].x, &p[i].y); for (i = 0; i < n-1; i++) addLine(l[i], p[i].x, p[i].y, p[i+1].x, p[i+1].y); addLine(l[i], p[i].x, p[i].y, p[0].x, p[0].y); HalfPlaneIntersect(); printf ("%.2lf\n", getArea()); } return 0; }
其实几何类的问题也可以分成两大类,不仅仅所有的几何问题都是叉积和设计繁琐的穷举,也会出现一些很偏几何的题目,这种时候,进行充分的数理分析是非常重要的。
Q3:
修水管问题:
话说最近柴小俊迷上了一款战略经营游戏《部落冲突》,他在修建筑的时候面临这样一个问题,他的军营被夹在两条河流之间,现在它想在两条河流的边上修筑供水点,并且需要用水管将两个供水点分别与军营连接起来(两个供水点之间也要相连),那么身为Acmer的你,能否当一次参谋,告诉柴小俊他最少需要准备多长的水管呢?
输入:我们以一条河流为x轴建系数,两条河流在(0,0)交汇。
点A的坐标(另一条河流上的一点,假定其在第一象限)。
点B的坐标(数据保证它落在两条直线的内部)。
输出:
输出柴小俊修建水管的最小长度。
分析:我们首先要做的应该是将这个实际问题抽象成数学(几何)问题,很容易看到,这道问题描述的是如下一个图形。
给出了直线OA,和夹在x、OA之间的B,那么现在需要你求解三角形BCD的最小周长,其中C、D必须分别落在OA和x轴上。
首先我们先进行最优解的分析,也就是我们得找到这个最优三角形。我们任意画一个三角形BCD,很容易看到,这里我们分别做B关于OA、x轴的对称点,记作E、F,那么根据中垂线的性质,我们发现三角形的周长变成了DC+CD+DF,那么为了使得和最小,应该令E、C、D、F共线,也就是说,连接E、F,与OA、x的交点就是最优三角形的顶点C、D。
那么现在问题的关键变成了怎么求一个点关于一条直线的对称点(点A关于直线l的对称点B)。
容易看到,我们可以把这个过程分为两个过程:
(1) 作点关于A在直线l上的投影,记作C。
(2) 基于投影,我们进行向量运算,点C + 向量CA。(想一想,为什么)
那么我们的问题再次简化:如何求A在l上的投影?
这里就需要进行一些简单的运算。
有了这些铺垫,我们就不难求得两个对称点,它们之间的距离也就是这道问题的最终解。
Q4:
等周问题.
话说柴小俊成功修完水管之后,实力大增(节省了材料),成功攻占了地方阵营,经过谈判,在敌方阵营中(一个三角形区域),柴小俊能够用长度为L的绳子圈出一块区域作为自己的领地,那么机智的你,能否告诉柴小俊他最多能够圈多大的面积呢?
输入:
给出三角形三个顶点的坐标、 绳长L。
输出:
柴小俊能够圈出的最大面积。
分析:首先我们脱离这个具体问题本身,来想清楚另外一个问题:周长L的所有平面几何图形中,面积最大的是谁?
其实这就是著名的等周问题,答案是圆。其最严谨的证明在数学分析的课本中往往会介绍变分法来证明,但是这里我们利用一个初等的技巧来简单的证明一下。
证明:
首先我们能够看到,最大面积的图形必然是“凸型曲线”,因为如果是凹陷的,我们能够将凹陷处向外翻折,周长没变但是显然面积增大了。
其次,我们应该看到,最大面积图形中,任取一点A,再做另外一点B,使得曲线AB的长度是L的一半,直线AB左右两侧的面积应该是相同的。因为如果不同,我们可以做那个面积较大的部分关于AB的对称,依然是周长没变,面积变大。
最后,基于上边的要点,我们发现,只要找到了不闭合曲线L/2和直线AB围成图形的最大面积,那么我们通过对称就可以得到最终的最大面积图形。这里在长度为L/2的不闭合曲线上取一个点C,假设这里弦AC及其对应的弧围成的面积是最大的、BC及其对应的弧围成的面积是最大的,那么现在最优解的问题又转化了:一直三角形的两边AC、BC,△ABC的最大面积是多少?很显然,由三角形的面积公式S=1/2absinC可知,C是直角的时候面积最大。那么在这段弧上处处取这样的点,都要满足这个结论,最终我们会得到一个半圆。
最后我们就得到了一个圆。
证毕。
那么我们现在再次回到这个问题,会出现如下三种情况:
(1)L小于三角形内切圆的周长,那么这表明L可以在三角形内部圈出一个完整的圆。
(2)L大于三角形的周长,那么圈出的面积就是这个三角形。
(3)介于二者之间。
能够看到,对于(1)(2)两种情况非常好处理,难处理的是(3)。
看这样一个图。
首先我们要搞明白这种情况下的曲线的性质。
这里首先确定这样一件事情:图中直线部分的长度和曲线部分长度分别都是定值,然后对于曲线部分,利用等周定理,我们应该拼出一个圆才能够使得面积最大。
也就是我们得到这样要给结论,这种情况下最大面积的图形,三段曲线是一个整圆的三部分。
然后依据性质做具体而有效的计算。
S(白) = S(大三角)-S(小三角)+S(小内切圆)
结合初等数学相似比的方法,便可求解。
简单的参考代码如下:
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;
const double eps=1e-9;
const double pi=acos(-1.0);
double a, b, c, l;
bool cmp(double a , double b)
{
return a < b ? 1 : 0;
}
int main()
{
double len, hl, cosc, area, res, r, k;
int tt=1;
double edge[3];
while(scanf("%lf%lf%lf%lf", &edge[0], &edge[1], &edge[2], &l) != EOF)
{
if(fabs(edge[0])<eps && fabs(edge[1])<eps && fabs(edge[2])<eps && fabs(l)<eps) break;
sort(edge,edge + 3,cmp);
a = edge[0],b = edge[1] , c = edge[2];
len=a+b+c;
cosc=(a*a + b*b - c*c)/(2*a*b);
if(fabs(cosc)<eps) area=a*b/2 , r=(a+b-c)*0.5;
else area=a*b*sqrt(1-cosc*cosc)*0.5, r=area*2.0/len;
if(2*pi*r>=l)
{
r = l/(2*pi);
res = pi*r*r;
}
else if(l>=len) res=area;
else
{
k=(len-l)/(len-2*pi*r);
r*=k;
res=area-area*k*k+pi*r*r;
}
printf("Case %d: %.2lf\n", tt++, res);
}
return 0;
}
Q5(hdu5476):
给出一个等腰三角形ABC,AB=AC,M是BC中点。求解三角形ABC内部点P的轨迹长度。其中P满足使得max(∠APC + ∠BPM , ∠APB + ∠MPC)最小。
分析:
容易将约束条件等价转化成∠APB + ∠MPC =∠APC + ∠BPM=180°.然后满足要求的P的轨迹就是下图AM+弧BC。
Q6(hdu 6206):
题目顺序给出四个点,判断第四个点是否在前三个点构成的圆内。
分析:直观上这道题目很好理解,但是这种题目往往卡精度,导致我们转化推导过程,用原始的三个点的横纵坐标来表示结果(传统的外心公式带来的误差就很大)。因此对于外心坐标的推到,我们就要用更加代数化的方法,而要避开几何求法。
下面给出不共线三个点外心坐标推导过程:
那么现在对于判断第四个点是否在该圆范围内,我们仅仅需要对比半径平方R^2和距离平方D^2的分母即可(消除除法对精度的影响)。由于这道题坐标的范围没有给出,数值上限会非常大,用java高精度解决。
总结一下,在处理一些简单的几何问题中,卡的并不是解题思路而是精度,这要求我们把几何化的、基于三角函数的推导转化成基于基本四则运算(还要避开除法)的代数化推导(往往求解线性方程组会变得很常见)。