我的計算幾何之路
計算幾何之路
計算幾何 Part.1---點，線，面，形基本關係，點積叉積的理解
計算幾何 Part.2---凸包問題
計算幾何 Part.3---面積公式矩形切割
計算幾何 Part.4---半平面交
計算幾何 Part.5---計算幾何背景，實際上解題的關鍵是其他問題（資料結構、組合數學，或者是列舉思想）若干道經典的離散化＋掃描線的題目
計算幾何 Part.6---解析幾何
計算幾何 Part.7---旋轉卡殼
計算幾何 Part.8---極角排序
END：推薦BLOG
計算幾何題的特點與做題要領：
1.大部分不會很難，少部分題目思路很巧妙
2.做計算幾何題目，模板很重要，模板必須高度可靠。
3.要注意程式碼的組織，因為計算幾何的題目很容易上兩百行程式碼，裡面大部分是模板。如果程式碼一片混亂，那麼會嚴重影響做題正確率。
4.注意精度控制。
5.能用整數的地方儘量用整數，要想到擴大資料的方法（擴大一倍，或擴大sqrt2）。因為整數不用考慮浮點誤差，而且運算比浮點快。
計算幾何 Part.1---點，線，面，形基本關係，點積叉積的理解
1.
題意：給了m個點，落在n+1個區域中，問各個區域有多少個點。
就是利用叉積去判斷點線上段的哪一側，可以二分去做，比較快。
程式碼：
#include <iostream>
#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <queue>
#include <map>
#include <vector>
#include <set>
#include <string>
#include <math.h>
using namespace std;
struct Point
{
int x,y;
Point(){}
Point(int _x,int _y)
{
x = _x;y = _y;
}
Point operator -(const Point &b)const
{
return Point(x - b.x,y - b.y);
}
int operator *(const Point &b)const
{
return x*b.x + y*b.y;
}
int operator ^(const Point &b)const
{
return x*b.y - y*b.x;
}
};
struct Line
{
Point s,e;
Line(){}
Line(Point _s,Point _e)
{
s = _s;e = _e;
}
};
int xmult(Point p0,Point p1,Point p2) //計算p0p1 X p0p2
{
return (p1-p0)^(p2-p0);
}
const int MAXN = 5050;
Line line[MAXN];
int ans[MAXN];
int main()
{
//freopen("in.txt","r",stdin);
//freopen("out.txt","w",stdout);
int n,m,x1,y1,x2,y2;
bool first = true;
while(scanf("%d",&n) == 1 && n)
{
if(first)first = false;
else printf("\n");
scanf("%d%d%d%d%d",&m,&x1,&y1,&x2,&y2);
int Ui,Li;
for(int i = 0;i < n;i++)
{
scanf("%d%d",&Ui,&Li);
line[i] = Line(Point(Ui,y1),Point(Li,y2));
}
line[n] = Line(Point(x2,y1),Point(x2,y2));
int x,y;
Point p;
memset(ans,0,sizeof(ans));
while( m-- )
{
scanf("%d%d",&x,&y);
p = Point(x,y);
int l = 0,r = n;
int tmp;
while( l <= r)
{
int mid = (l + r)/2;
if(xmult(p,line[mid].s,line[mid].e) < 0)
{
tmp = mid;
r = mid - 1;
}
else l = mid + 1;
}
ans[tmp]++;
}
for(int i = 0; i <= n;i++)
printf("%d: %d\n",i,ans[i]);
}
return 0;
}
先判斷兩條直線是不是同線，不是的話再判斷是否平行，再不是的話就只能是相交的，求出交點。
如何判斷是否同線？由叉積的原理知道如果p1，p2，p3共線的話那麼(p2-p1)X(p3-p1)=0。因此如果p1，p2，p3共線，p1，p2，p4共線，那麼兩條直線共線。direction（）求叉積，叉積為0說明共線。
如何判斷是否平行？由向量可以判斷出兩直線是否平行。如果兩直線平行，那麼向量p1p2、p3p4也是平等的。即((p1.x-p2.x)*(p3.y-p4.y)-(p1.y-p2.y)*(p3.x-p4.x))==0說明向量平等。
如何求出交點？這裡也用到叉積的原理。假設交點為p0(x0,y0)。則有：
(p1-p0)X(p2-p0)=0
(p3-p0)X(p2-p0)=0
展開後即是
(y1-y2)x0+(x2-x1)y0+x1y2-x2y1=0
(y3-y4)x0+(x4-x3)y0+x3y4-x4y3=0
將x0,y0作為變數求解二元一次方程組。
假設有二元一次方程組
a1x+b1y+c1=0;
a2x+b2y+c2=0
那麼
x=(c1*b2-c2*b1)/(a2*b1-a1*b2);
y=(a2*c1-a1*c2)/(a1*b2-a2*b1);
因為此處兩直線不會平行，所以分母不會為0。
#include<iostream>
#include<stdio.h>
using namespace std;
const double epx=1e-10;
struct Point
{
double x;
double y;
};
//求解二元一次方程
Point solve(double a1,double b1,doublec1,double a2,double b2,double c2)
{
Point p;
p.x=(c1*b2-c2*b1)/(a2*b1-a1*b2);
p.y=(a2*c1-a1*c2)/(a1*b2-a2*b1);
return p;
}
//p1p3,p1p2的叉積
double direction(Point p1,Point p2,Pointp3)
{
return (p3.x-p1.x)*(p2.y-p1.y)-(p2.x-p1.x)*(p3.y-p1.y);
}
//判斷兩直線的關係
/*
在這裡有三種關係：1共線 2平行 3相交
1 共線可通過叉積來判斷
2 平行通過向量來判斷
3 通了上面2種情況的其他情況
求交點可通過叉積及解二元一次方程來求解
*/
int N;
Point p1,p2,p3,p4;
Point p0;//交點
double a1,b1,c1,a2,b2,c2;
int main()
{
scanf("%d",&N);
int i;
printf("INTERSECTING LINES OUTPUT\n");
for(i=0;i<N;++i)
{
scanf("%lf%lf%lf%lf%lf%lf%lf%lf", &p1.x, &p1.y,&p2.x, &p2.y, &p3.x, &p3.y, &p4.x, &p4.y);
if(direction(p3,p4,p1)==0 && direction(p3,p4,p2)==0)//共線
printf("LINE\n");
else
{
if( ((p1.x-p2.x)*(p3.y-p4.y)-(p1.y-p2.y)*(p3.x-p4.x))==0 )//平行
printf("NONE\n");
else
{
a1=p1.y-p2.y;b1=p2.x-p1.x;c1=p1.x*p2.y-p2.x*p1.y;
a2=p3.y-p4.y;b2=p4.x-p3.x;c2=p3.x*p4.y-p4.x*p3.y;
p0=solve(a1,b1,c1,a2,b2,c2);
printf("POINT %.2f %.2f\n",p0.x,p0.y);
}
}
}
printf("END OF OUTPUT\n");
return 0;
}
題目大意：給出n條線段兩個端點的座標，問所有線段投影到一條直線上，如果這些所有投影至少相交於一點就輸出Yes！，否則輸出No！。
解題思路：如果有存在這樣的直線，過投影相交區域作直線的垂線，該垂線必定與每條線段相交，問題轉化為問是否存在一條線和所有線段相交
直線肯定經過兩個端點。
列舉端點，判斷直線和線段是否相交。
細節要注意，判斷重合點。
還有就是加入只有一條線段的話，剛好直線是過同一條直線的。
所以保險的做法是列舉所有的兩個端點，包括同一條直線的。
#include<cstdio>
#include<string>
#include<cstring>
#include<iostream>
#include<cmath>
#include<algorithm>
#include<vector>
using namespace std;
#define all(x) (x).begin(), (x).end()
#define for0(a, n) for (int (a) = 0; (a)< (n); (a)++)
#define for1(a, n) for (int (a) = 1; (a)<= (n); (a)++)
#define mes(a,x,s) memset(a,x,(s)*sizeof a[0])
#define mem(a,x) memset(a,x,sizeof a)
#define ysk(x) (1<<(x))
typedef long long ll;
typedef pair<int, int> pii;
const int INF =0x3f3f3f3f;
const int maxn= 100 ;
int n;
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 -(Vector A,Vector 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); }
Vector operator -(Vector A) {return Vector(-A.x,-A.y);}
const double eps=1e-8;
int dcmp(double x)
{
if(fabs(x)<eps) return 0;
else return x<0?-1:1;
}
bool operator ==(Point A,Point B) {returndcmp(A.x-B.x)==0&&dcmp(A.y-B.y)==0;}
double Cross(Vector A,Vector B)//叉乘
{
return A.x*B.y-A.y*B.x;
}
struct Line
{
Point p,p2;
Vector v;
Line(){}
Line(Point p,Vector v):p(p),v(v){}//點線式
void twoPointIntial(Point p,Point p2)//兩點式
{
this->p=p;
this->p2=p2;
v= p2-p;
}
}seg[maxn+5];
typedef Line Seg;
bool intersect(Line a,Seg b)//判斷直線和線段是否相交
{
Point A=b.p;
Point B=b.p2;
Vector v1=a.p-A;
Vector v2=a.p2-A;
Vector v3=a.p-B;
Vector v4=a.p2-B;
return dcmp(Cross(v1,v2))*dcmp(Cross(v3,v4))<=0; //如果線段的端點在直線上,肯定相交
}
bool solve(Point A,Point B)
{
Line L;L.twoPointIntial(A,B);
if( A==B ) return false;
for(int i=0;i<n;i++)
{
if(!intersect(L,seg[i]) ) return false;
}
return true;
}
int main()
{
std::ios::sync_with_stdio(false);
int T;cin>>T;
while(T--)
{
cin>>n;
for0(i,n)
{
cin>>seg[i].p.x>>seg[i].p.y>>seg[i].p2.x>>seg[i].p2.y;
}
if(n==1) {puts("Yes!");continue;}
bool find=false;
for(int i=0;i<n&&!find;i++)
{
for(int j=i+1;j<n;j++)
{
find|= solve(seg[i].p ,seg[j].p);
find|= solve(seg[i].p ,seg[j].p2);
find|= solve(seg[i].p2 ,seg[j].p);
find|= solve(seg[i].p2 ,seg[j].p2);
if(find) break;
}
}
puts(find?"Yes!":"No!");
}
return 0;
}
題目大意：有n根木條，一根一根的往一個座標系上丟(給出木條兩點的座標)，問最後不被覆蓋的木條有哪些，即丟的木條如果和前面丟的木條交叉的話，就會覆蓋前面那根木條解題思路：用到線段交叉的判斷，我用佇列模擬是否被覆蓋，沒被覆蓋的入隊，並且將當前火柴和前面沒被覆蓋的木條比較如果交叉則前面那根出隊，這樣的暴力做法在時間上不知道怎麼優化了，不過比較了下，我的記憶體還是很低的
#include <iostream>
#include <queue>
#include <stdio.h>
using namespace std;
#define min(a,b) (a>b?b:a)
#define max(a,b) (a<b?b:a)
typedef struct
{
intcnt;
doublex1,y1,x2,y2;
}Line;
int n;
double direction(double x,double y,doublex1,double y1,double x2,double y2)
{ //叉積
doublea1=x1-x;
doubleb1=y1-y;
doublea2=x2-x;
doubleb2=y2-y;
returna1*b2-a2*b1;
}
int on_segment(double x1,double y1,doublex2,double y2,double x,double y)
{ //判斷共線
if((min(x1,x2)<=x&& x<=max(x1,x2)) && (min(y1,y2)<=y &&y<=max(y1,y2)))
return1;
return0;
}
int cross(Line v,Line t)
{ //是否交叉，根據演算法導論寫的
doubled1,d2,d3,d4;
d1=direction(t.x1,t.y1,t.x2,t.y2,v.x1,v.y1); //算叉積
d2=direction(t.x1,t.y1,t.x2,t.y2,v.x2,v.y2);
d3=direction(v.x1,v.y1,v.x2,v.y2,t.x1,t.y1);
d4=direction(v.x1,v.y1,v.x2,v.y2,t.x2,t.y2);
if(d1*d2<0&& d3*d4<0) return 1; //直接和0比較的話時間是625MS
if(!d1&& on_segment(t.x1,t.y1,t.x2,t.y2,v.x1,v.y1)) return 1;
if(!d2&& on_segment(t.x1,t.y1,t.x2,t.y2,v.x2,v.y2)) return 1;
if(!d3&& on_segment(v.x1,v.y1,v.x2,v.y2,t.x1,t.y1)) return 1;
if(!d4&& on_segment(v.x1,v.y1,v.x2,v.y2,t.x2,t.y2)) return 1;
return0;
}
void input()
{
queue<Line>q;
Linev,t;
inti;
scanf("%lf%lf%lf%lf",&v.x1,&v.y1,&v.x2,&v.y2);
v.cnt=1;
q.push(v);
for(i=2;i<=n;i++)
{
scanf("%lf%lf%lf%lf",&t.x1,&t.y1,&t.x2,&t.y2);
t.cnt=i;
q.push(t); //用當前的作為隊尾
while(!q.empty())
{
v=q.front();q.pop(); //一個一個出隊，直到隊尾
if(t.cnt==v.cnt)
{
q.push(t);
break;
}
if(!cross(v,t))q.push(v); //如果不交叉繼續入隊
}
}
v=q.front();q.pop();
printf("Topsticks: %d",v.cnt);
while(!q.empty())
{
v=q.front();
q.pop();
printf(",%d",v.cnt);
}
printf(".\n");
}
int main()
{
while(scanf("%d",&n)&& n)
{
input();
}
return0;
}
2.
問題描述
為了尋找失去的愛人Cupid，Psyche需要完成Venus的最後一項任務：前往冥界，收集一盒冥界皇后Prosperina的美貌。
冥界有nn個神殿，可以表示為平面上的nn個整點。Psyche想要找到這nn座神殿中，最近的兩座神殿之間的距離。傳說那就是通往主殿的密碼。
但是冥界神祕莫測，在不同的時刻，這nn座神殿中的某一座會消失。
Psyche想要知道，對於nn座神殿中的任意一座消失的情況，最近的兩座神殿之間的距離。你只需要輸出它們的和。
為避免精度誤差，定義兩點(x_1, y_1), (x_2, y_2)(x1,y1),(x2,y2)間的距離為d = (x_1 - x_2) ^ 2 + (y_1 - y_2) ^ 2d=(x1−x2)2+(y1−y2)2。
輸入描述
第一行，一個整數TT(1 \le T \le 5)(1≤T≤5)，代表資料組數。
對於每組資料，第一行，一個整數nn(3 \le n \le 10 ^ 5)(3≤n≤105)，代表神殿個數。
下面nn行，每行兩個整數x, yx,y(-10 ^ 5 \le x,y \le 10 ^ 5)(−105≤x,y≤105)，代表神殿的位置在(x, y)(x,y)。
注意可能存在兩座神殿坐落在同一位置。
輸出描述
輸出TT行，對於每組資料，輸出nn座神殿中的任意一座消失的情況，最近兩座神殿之間的距離的和。
輸入樣例
1
3
0 0
1 1
2 2
輸出樣例
12
Hint
神殿(0,0)(0,0)消失時，d = (1-2) ^ 2 + (1 - 2) ^ 2 = 2d=(1−2)2+(1−2)2=2；
神殿(1,1)(1,1)消失時，d = (0-2) ^ 2 + (0 - 2) ^ 2 = 8d=(0−2)2+(0−2)2=8；
神殿(2,2)(2,2)消失時，d = (0-1) ^ 2 + (0-1) ^ 2 = 2d=(0−1)2+(0−1)2=2；
故答案為2 + 8 + 2 = 122+8+2=12。
解題思路：
感覺就是很裸的平面最近點對的模板題= =題解上的說的也很明確了，只有三種情況
一種情況是刪的點不屬於最近點對裡面的點
剩下兩種情況各刪的最近點對裡面的其中一點
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
#define INF 1e18
#define Min(a, b) ((a) <= (b) ? (a) :(b))
typedef long long LL;
typedef struct node{
LL x, y;
LL index;
}Coord;
const int maxn = 1e5 + 5;
LL p1, p2, ans;
Coord p[maxn], ppp[maxn], tmp[maxn];
inline LL f(LL x){
return (x >= 0 ? x : -x);
}
inline bool cmp(Coord a, Coord b){
return a.x < b.x;
}
inline bool cmp2(Coord a, Coord b){
return a.y < b.y;
}
inline LL Dist(Coord a, Coord b){
return (a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y);
}
LL Cal_Closest(LL left, LL right, intstep){
LL d = INF;
if(left == right){
return d;
}
if(left + 1 == right){
LL tmp2 = Dist(p[left], p[right]);
if(step == 0){
if(step == 0 && tmp2< ans){
p1 = left;
p2 = right;
ans = tmp2;
}
}
return tmp2;
}
LL mid = (left + right) >> 1;
LL d1 = Cal_Closest(left, mid, step);
LL d2 = Cal_Closest(mid + 1, right, step);
d= Min(d1, d2);
LL k = 0;
for(LL i = left; i <= right; ++i){
if(f(p[mid].x - p[i].x) <= d){
tmp[k++] = p[i];
}
}
sort(tmp, tmp + k, cmp2);
for(LL i = 0; i < k; ++i){
for(LL j = i + 1; j < k && tmp[j].y - tmp[i].y < d; ++j){
if(d > Dist(tmp[i], tmp[j])){
d = Dist(tmp[i], tmp[j]);
if(step == 0 && d <ans){
p1 = tmp[i].index;
p2 = tmp[j].index;
ans = d;
}
}
}
}
return d;
}
int main()
{
LL t, n;
scanf("%lld", &t);
while(t--){
scanf("%lld", &n);
memset(p, 0, sizeof(p));
memset(ppp, 0, sizeof(ppp));
for(LL i = 0; i < n; ++i){
scanf("%lld%lld", &p[i].x, &p[i].y);
}
sort(p, p + n, cmp);
for(LL i = 0; i < n; ++i) p[i].index = i;
for(LL i = 0; i < n; ++i) ppp[i] = p[i];
p1 = 0; p2 = 0;
ans = INF;
ans = Cal_Closest(0, n - 1, 0);
ans = ans * (n - 2);
LL t = 0;
for(LL i = 0; i < n; ++i){
if(ppp[i].index == p1) continue;
p[t++] = ppp[i];
}
sort(p, p + t, cmp);
ans += Cal_Closest(0, t - 1, 1);
t = 0;
for(LL i = 0; i < n; ++i){
if(ppp[i].index == p2) continue;
p[t++] = ppp[i];
}
sort(p, p + t, cmp);
ans += Cal_Closest(0, t - 1, 1);
printf("%lld\n", ans);
}
return 0;
}
計算幾何 Part.2---凸包問題
關於凸包的嚴格定義，這裡不打算寫出來，大家可以自行Google或者百度，因為嚴格的數學定義反而不太好理解，用最通俗的話來解釋凸包：給定二維平面上的點集，凸包就是將最外層的點連線起來構成的凸多邊型，它能包含點集中所有的點。
1.
acm.hdu.edu.cn/showproblem.php?pid=3847
/*
求出多邊形最窄的地段長度
列舉邊，求出所有點中到邊的距離最大的值
這些值中最小的就是答案
*/
#include <cstdio>
#include <cstring>
#include <cmath>
#include <cstdlib>
const int MAXN = 109;
const double eps = 1e-4;
struct point{
double x,y;
}p[MAXN],h[MAXN];
inline double distance(const point &p1,const point &p2){
return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
}// 求兩點之間的距離
inline double multiply(const point &sp,const point &ep,const point &op){
return ((sp.x-op.x)*(ep.y-op.y)-(ep.x-op.x)*(sp.y-op.y));
}//判斷sp,ep，op是否滿足左轉
int cmp(const void *a,const void *b){//按極角排序
point *p1 = (point *)a;
point *p2 = (point *)b;
double t = multiply(*p2,*p1,p[0]);
if(t>eps) return 1;
else if(fabs(t)<=eps)
{
if(distance(p[0],*p1)>distance(p[0],*p2)) return 1;
else return -1;
}
else return -1;
}
void anglesort(point p[],int n){//找到最左下方的點
int i,k=0;
point temp;
for(i=1;i<n;i++)
if(p[i].x<p[k].x ||( fabs(p[i].x-p[k].x)<eps && (p[i].y<p[k].y)))
k=i;
temp=p[0],p[0]=p[k],p[k]=temp;
qsort(p+1,n-1,sizeof(point),cmp);
}
void Graham_scan(point p[],point ch[],int n,int &len){//建立凸包
int i,top=2;
anglesort(p,n);
if(n<3){
for(i=0,len=n;i<n;i++) ch[i]=p[i];
return;
}
ch[0]=p[0],ch[1]=p[1],ch[2]=p[2];
p[n]=p[0];
for(i=3;i<n;i++){
while(multiply(p[i],ch[top],ch[top-1])>=eps) top--;
ch[++top]=p[i];
}
len=top+1;
}
double judge(point _x,point _y,int len)
{
double res=0;
for(int i=0;i<len;i++)
{
double tmp=fabs(multiply(h[i],_x,_y))/distance(_x,_y);//面積法求出距離
if(tmp>res)
res=tmp;
}
return res;
}
double solve(int len)
{
h[len]=h[0];
double ans=1<<30;
for(int i=0;i<len;i++)
{
double tmp=judge(h[i],h[i+1],len);//列舉邊
if(tmp<ans)
ans=tmp;
}
return ans;
}
int main(){
int n,len;
int d=0;
while(scanf("%d",&n),n)
{
for(int i=0;i<n;i++) scanf("%lf %lf",&p[i].x,&p[i].y);
Graham_scan(p,h,n,len);
printf("Case %d: %.2lf\n",++d,solve(len)+0.005);
}
return 0;
}
2.
大致題意：
給定平面上的一些散點集，求最遠兩點距離的平方值。
解題思路：
別想著暴力列舉任意亮點距離找最大，行不通，想想三點共線吧！
平面上的散點集的最遠的兩點距離必然在這個散點集的凸包的某兩個頂點上出現。
那麼先求凸包，再列舉頂點距離就OK了。
別看是3000ms就想用簡單的捲包裹，這題資料規模極大，捲包裹鐵超（我一開始就是這麼做的。。。） 萬般無奈不得不用GrahamScanAlgorithm。。。。O(nlogn)用來做這題還是相當可觀的。
GrahamScan理解是不困難的，同學們百度搜搜就知道基本思想了：資料結構：棧。而且每個點最多入棧出棧一次。
問題是用棧構造凸包之前要對散點集進行排序，水平序我不說了，我沒看懂，聽說很簡單。我用的是極座標排序。
利用叉積的旋向 配合 比較排序 （如插入排序，冒泡，快排等）可以對極座標排序，推薦用qsort，不要用冒泡之類的，太慢了，用Graham都是想快而已。
qsort（）難點在於 比較規則，（我的程式為cmp函式），必須把“qsort對cmp返回值的處理、cmp返回值、叉積旋向返回值”三者有機結合，否則一定排序出錯，詳見我的程式，有詳細解釋。
Cmp比較規則為“按極角大小逆時針排序，極角相等則近極座標的點優先”。在網上看到有些童鞋在“極角相同”的情況下，利用第二關鍵字“距離”計算排序時，是通過計算兩點的橫座標到 原點橫座標之差 作為判斷“距離”依據。乍一看好像沒錯，也能AC，其實那是POJ資料庫不強大，試試多點與原點的連線垂直於x軸的情況吧！
最後注意的：在使用qsort之前，必須先找到 散點集中 最左下角的點，把它放在散點集陣列的最後一個元素位，作為 極座標原點（快排的參考點），否則排序也會出錯。
其實只要qsort的cmp函式能處理好了，基本這題就過了，難度不大。
#include<iostream>
#include<algorithm>
using namespace std;
const int inf=10001;
typedef class location
{
public:
int x,y;
}node;
/*點A到點B的距離的平方*/
int distsquare(node A,node B)
{
return (B.x-A.x)*(B.x-A.x)+(B.y-A.y)*(B.y-A.y);
}
/*計算叉積*/
int det(int x1,int y1,int x2,int y2)
{
return x1*y2-x2*y1;
}
int cross(node A,node B,node C,node D)
{
return det(B.x-A.x , B.y-A.y , D.x-C.x , D.y-C.y);
}
/*qsort比較規則*/
node* p; //極座標原點
int cmp(const void* a,const void* b)
{
node* s=(node*)a;
node* t=(node*)b;
int temp=cross(*p,*s,*p,*t); //叉乘ps X pt
if(temp>0) //說明pt向量的極角小於 ps向量的極角
return -1; //從逆時針排序角度，t點位置在s點前面，即t<s ，根據qsort規則返回-1
else if(temp==0) //說明pt向量的極角 等於 ps向量的極角
return distsquare(*p,*t)-distsquare(*p,*s); //距離原點近的點優先排序，用減法能夠實現3出口：- 0 +
//注意，網上有些程式在這裡不是比較 距離之差，而是比較橫座標絕對值 之差
//這是欠缺考慮 多點與原點連線垂直於x軸，不完善，之所以能AC是因為POJ的資料庫不足夠大而已
else
return 1; //pt向量的極角 大於 ps向量的極角
}
int main()
{
int n,i,j;
while(cin>>n)
{
node* farm=new node[n+1];
int* conbag=new int[n+2]; //conbag[]順序記錄輸入的點中構造成凸包的頂點集的各點在farm[]中的下標
/*Input & search the first vertex*/
int min_x=inf;
int fi=0;
for(i=1;i<=n;i++)
{
cin>>farm[i].x>>farm[i].y;
if(min_x > farm[i].x)
{
min_x = farm[i].x; //先以橫座標為第一關鍵字搜尋最左下角的點
fi=i;
}
else if(min_x == farm[i].x)
{
if(farm[fi].y >farm[i].y) //若第一關鍵字相同，則以縱座標作為第二關鍵搜尋
fi=i;
}
}
/*Quicksort Point Sets*/
farm[0]=farm[n]; //這三步非常關鍵，是使用qsort前的準備工作
farm[n]=farm[fi]; //目的是把前面找到的最左下角的點作為 極座標原點
farm[fi]=farm[0]; //即把第fi個點移到farm[]最後的位置,qsort則會把它作為排序的參考點
p=&farm[n]; //極座標原點傳參
qsort(farm+1,n,sizeof(node),cmp); //farm[]散點集排序
/*Graham Scan Algorithm*/
int pc=1; //conbag[]指標
conbag[1]=n; //(約定)極座標原點為凸包第1個頂點
conbag[++pc]=1; //(在前面基礎上,)有序點集farm[]的第一個點 (必)為凸包第2個頂點
conbag[0]=2; //凸包頂點計數器
for(i=2;i<=n;)
if(cross(farm[ conbag[pc-1] ],farm[ conbag[pc] ],farm[ conbag[pc]],farm[i]) >= 0)
{ //檢查向量(前一點pc-1,當前點pc) 與 向量(當前點pc,待入棧點i) 的旋轉關係
conbag[++pc]=i++; //入棧
conbag[0]++;
}
else
{
pc--; //當前點pc出棧
conbag[0]--;
}
/*Search The Max distant*/
int maxdist=0;
for(i=1;i<=conbag[0]-2;i++) //散點集的兩點最大距離必定出現在該散點集的凸包的某兩個頂點之間
for(j=i+1;j<=conbag[0]-1;j++)
{
int temp=distsquare(farm[conbag[i] ],farm[ conbag[j] ]);
if(maxdist < temp)
maxdist = temp;
}
/*Output*/
cout<<maxdist<<endl;
delete farm;
delete conbag;
}
return 0;
}
計算幾何 Part.3---面積公式矩形切割
題意：以原點為圓心，給出圓上的一點，要求兩位兩點，是的這三個點的距離和最大，很容易想到這是一個等邊三角形，而事實上，經過對題目給出樣例的測試也證明了這確實是一個等邊三角形
思路：幾何水題
我們可以得到方程組
x^2+y^2 = r^2
(a-x)^2+(b-y)^2=3r^2
解方程組得到的兩點即為三角形的另外兩點
#include <stdio.h>
#include <math.h>
int main()
{
int t;
double x,y,x2,y2,r;
double ax,ay,bx,by,k,m,l,A,B,C;
scanf("%d",&t);
while(t--)
{
scanf("%lf%lf",&x,&y);
r = x*x+y*y;
A = r;
B = y*r;
C = r*r/4-r*x*x;
ay = (-B-sqrt(B*B-4*A*C))/(2*A);
by = (-B+sqrt(B*B-4*A*C))/(2*A);
if(fabs(x-0)<1e-7)//防止除數出現0的情況
{
ax=-sqrt(r-ay*ay);
bx=sqrt(r-by*by);
}
else
{
ax=(-r/2-ay*y)/x;//由於ay必定小於by，所以ax也必定小於bx，所以無需進行大小判定
bx=(-r/2-by*y)/x;
}
printf("%.3lf %.3lf %.3lf %.3lf\n",ax,ay,bx,by);
}
return 0;
}
計算幾何 Part.4---半平面交
1.
題目：鐵人三項，每個人在某一項中有確定的速度，裁判可以決定某一項比賽的路程為多少，問對於某個人，是否存在一種安排能使他拿到第一，而且不能是並列。
我們假設三項的路程分別人X,Y,Z。
比較其中的兩個人。A的時間為X / U1+Y / V1+Z /W1 B的時間為X / U2 +Y / V2 +Z / W2
如果A想要獲勝媽，X / U1+Y / V1+Z / W1- X / U2 +Y / V2 +Z / W2 < 0
由於我寫的是順時針的，所以把不等式變個號。這裡一定要注意，小於0大於0和順時針逆時針的對應關係
這個不等式還是有3個未知數。顯然用三維就太麻煩了。由於我們最終不需要求出X,Y,Z具體為多少，而且Z>0，所以把不等式兩邊同時除以Z，則把X/Z看成一個未知量，Y/Z看成另外一個。
問題轉化成一系列的不等式是否為解，而且注意條件X>0Y>0
那麼可以設立一個初始範圍，(0,0)(0,inf)(inf,inf)(inf,0)
然後通過兩個人的引數，求出AX+BY+C>0，通過半平面交解決，最終判斷面積是否為0
有幾個地方需要注意：
這題要求的精度很高，大家一味的說需要1e-16，其實不然，1e-8也過了。主要是中間的處理細節，對於1/U1-1/U2，普通的處理是需要兩次除法，精度嚴重受損，可以改成(U2-U1)/(U1*U2)。
另外就是特判，題目要求是不能並列，所以最終結果是大於0才行。而且如果遇到A==0&&B==0&&C<=0說明不等式無解，直接返回
#include <stdio.h>
#include <string.h>
#include <iostream>
#include <algorithm>
#include <vector>
#include <queue>
#include <set>
#include <map>
#include <string>
#include <math.h>
#include <stdlib.h>
#include <time.h>
using namespace std;
const double eps = 1e-18;
int sgn(double x)
{
if(fabs(x) < eps)return 0;
if(x < 0)return -1;
elsereturn 1;
}
struct Point
{
double x,y;
Point(){}
Point(double _x,double _y)
{
x = _x; y = _y;
}
Point operator -(const Point &b)const
{
return Point(x - b.x, y - b.y);
}
double operator ^(const Point &b)const
{
return x*b.y - y*b.x;
}
double operator *(const Point &b)const
{
return x*b.x + y*b.y;
}
};
//計算多邊形面積
double CalcArea(Point p[],intn)
{
double res = 0;
for(int i = 0;i < n;i++)
res += (p[i]^p[(i+1)%n]);
return fabs(res/2);
}
//通過兩點，確定直線方程
void Get_equation(Pointp1,Point p2,double &a,double &b,double &c)
{
a = p2.y - p1.y;
b = p1.x - p2.x;
c = p2.x*p1.y - p1.x*p2.y;
}
//求交點
Point Intersection(Pointp1,Point p2,double a,double b,double c)
{
double u = fabs(a*p1.x + b*p1.y + c);
double v = fabs(a*p2.x + b*p2.y + c);
Point t;
t.x = (p1.x*v + p2.x*u)/(u+v);
t.y = (p1.y*v + p2.y*u)/(u+v);
return t;
}
Point tp[110];
void Cut(double a,doubleb,double c,Point p[],int &cnt)
{
int tmp = 0;
for(int i = 1;i <= cnt;i++)
{
//當前點在左側，逆時針的點
if(a*p[i].x + b*p[i].y + c <eps)tp[++tmp] = p[i];
else
{
if(a*p[i-1].x + b*p[i-1].y + c <-eps)
tp[++tmp] =Intersection(p[i-1],p[i],a,b,c);
if(a*p[i+1].x + b*p[i+1].y + c <-eps)
tp[++tmp] =Intersection(p[i],p[i+1],a,b,c);
}
}
for(int i = 1;i <= tmp;i++)
p[i] = tp[i];
p[0] = p[tmp];
p[tmp+1] = p[1];
cnt = tmp;
}
double V[110],U[110],W[110];
int n;
const double INF = 100000000000.0;
Point p[110];
bool solve(int id)
{
p[1] = Point(0,0);
p[2] = Point(INF,0);
p[3] = Point(INF,INF);
p[4] = Point(0,INF);
p[0] = p[4];
p[5] = p[1];
int cnt = 4;
for(int i = 0;i < n;i++)
if(i != id)
{
double a = (V[i] -V[id])/(V[i]*V[id]);
double b = (U[i] -U[id])/(U[i]*U[id]);
double c = (W[i] -W[id])/(W[i]*W[id]);
if(sgn(a) == 0 && sgn(b) ==0)
{
if(sgn(c) >= 0)return false;
else continue;
}
Cut(a,b,c,p,cnt);
}
if(sgn(CalcArea(p,cnt)) == 0)return false;
else return true;
}
int main()
{
//freopen("in.txt","r",stdin);
//freopen("out.txt","w",stdout);
while(scanf("%d",&n) == 1)
{
for(int i = 0;i < n;i++)
scanf("%lf%lf%lf",&V[i],&U[i],&W[i]);
for(int i = 0;i < n;i++)
{
if(solve(i))printf("Yes\n");
else printf("No\n");
}
}
return 0;
}
2.
對於多邊形核不懂的，自行查資料了。
第一道半平面交，只會寫N^2。
這裡預設是順時針的順序，否則就要調整一下。
將每條邊化作一個不等式，ax+by+c>0，所以要固定順序，方便求解。
半平面交其實就是對一系列的不等式組進行求解可行解。
如果某點在直線右側，說明那個點在區域內，否則出現在左邊，就可能會有交點，將交點求出加入。
#include<iostream>
#include<fstream>
#include<iomanip>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cstdlib>
#include<cmath>
#include<set>
#include<map>
#include<queue>
#include<stack>
#include<string>
#include<vector>
#include<sstream>
#include<cassert>
#define LL long long
#define eps 1e-7
#define inf 1<<30
using namespace std;
struct Point{
double x,y;
}p[1505],tp[1505],q[1505];
//叉積
double xmul(Point p0,Point p1,Point p2){
return(p1.x-p0.x)*(p2.y-p0.y)-(p1.y-p0.y)*(p2.x-p0.x);
}
//通過兩點，確定直線方程
double Get_equation(Point p1,Point p2,double &a,double&b,double &c){
a=p2.y-p1.y;
b=p1.x-p2.x;
c=p2.x*p1.y-p1.x*p2.y;
}
//求交點
Point Intersection(Point p1,Point p2,double a,double b,double c){
doubleu=fabs(a*p1.x+b*p1.y+c);
doublev=fabs(a*p2.x+b*p2.y+c);
Point t;
t.x=(p1.x*v+p2.x*u)/(u+v);t.y=(p1.y*v+p2.y*u)/(u+v);
return t;
}
//求面積，正為順時針，和叉積寫法有關
double Get_area(Point p[],int n){
double area=0;
for(int i=2;i<n;i++)
area+=xmul(p[1],p[i],p[i+1]);
return -area/2.0;
}
//改變順序
double Change_dir(Point p[],int n){
for(int i=1;i<=n/2;i++)
swap(p[i],p[n+1-i]);
}
//加入一條邊，切割
void Cut(double a,double b,double c,Point p[],int &cnt){
int tmp=0;
for(inti=1;i<=cnt;i++){
//當前點就在右側
if(a*p[i].x+b*p[i].y+c>-eps) tp[++tmp]=p[i];
else{
//前一個點在右側，產生交點
if(a*p[i-1].x+b*p[i-1].y+c>eps)
tp[++tmp]=Intersection(p[i-1],p[i],a,b,c);
//同理
if(a*p[i+1].x+b*p[i+1].y+c>eps)
tp[++tmp]=Intersection(p[i],p[i+1],a,b,c);
}
}
for(int i=1;i<=tmp;i++)
p[i]=tp[i];
p[0]=p[tmp];p[tmp+1]=p[1];
cnt=tmp;
}
void slove(Point p[],int n){
//預設順時針，通過面積判斷一下
if(Get_area(p,n)<eps)Change_dir(p,n);
p[0]=p[n];p[n+1]=p[1];
//原來的點要備份一遍，查了好久
for(int i=0;i<=n+1;i++)q[i]=p[i];
int cnt=n;
for(int i=1;i<=n;i++){
double a,b,c;
Get_equation(q[i],q[i+1],a,b,c);
Cut(a,b,c,p,cnt);
}
printf("%.2f\n",fabs(Get_area(p,cnt)));
}
int main(){
int t,n;
scanf("%d",&t);
while(t--){
scanf("%d",&n);
for(inti=1;i<=n;i++)
scanf("%lf%lf",&p[i].x,&p[i].y);
slove(p,n);
}
return 0;
}
3.
題意：在（0，10000）*（0，10000）的座標系上，給定n個半平面，求出它們圍成的圖形的面積
每個半平面由兩點（x1，y1）（x2，y2）確定的直線確定，規定半平面為直線的左邊，即存在一點（x，y）
使得(x – x1) *(y – y2) – (x – x2) * (y – y1) = (x1 – x) * (y2 – y) – (x2 – x) * (y1 – y)>=0成立
//模板題直接copy網上模板
#include <iostream>
#include <cstdio>
#include <cmath>
#include <algorithm>
using namespace std;
const double eps = 1e-8;
int n, pn, dq[20005], top, bot;
struct Point {
double x, y;
} p[20005];
struct Line {
Point a, b;
double angle;
Line& operator= (Linel) {
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[20005];
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) returndblcmp(multi(l1.a, l2.a, l2.b)) > 0; //極角相同時，將更靠半平面裡面的放在前面
return d < 0;
}
void addLine(Line& l, double x1, double y1, double x2, doubley2) {
l.a.x = x1; l.a.y = y1;
l.b.x = x2; l.b.y = y2;
/*atan2(y,x)返回值的取值範圍為-PI到PI，實際上就是根據向量（x,y）確定極角，
若向量在1，2象限，則值大於0；3，4象限，則值小於0
所以atan2（1，1）與atan（-1，-1）值不同
*/
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;
}
void HalfPlaneIntersect( ){
int i, j;
/*排序是在滿足所有半平面A*x+B*y+C>0或（<,<=,>=）,
也就是所有半平面的符號均相同的情況下對極角進行排序，
此題中題目對半平面的定義就等價於這個條件
*/
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的條件後，
它又會影響到底部的兩條直線，當底部的兩條直線的交點，在當前的半平面外部時，底部彈棧
*/
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 i;
double x1, y1, x2, y2;
while (scanf("%d", &n) != EOF) {
//輸入半平面，由一條線段確定，半平面為線上段的哪邊題目已經給出
for (i = 0; i < n;i++) {
scanf("%lf%lf%lf%lf", &x1, &y1, &x2, &y2);
addLine(l[i], x1,y1, x2, y2);
}
//這四個半平面的新增要滿足題目給出的半平面的定義
addLine(l[n++], 0, 0,10000, 0);
addLine(l[n++], 10000,0, 10000, 10000);
addLine(l[n++], 10000,10000, 0, 10000);
addLine(l[n++], 0,10000, 0, 0);
HalfPlaneIntersect();
printf("%.1lf\n", getArea());
}
return 0;
}
計算幾何 Part.5---計算幾何背景，實際上解題的關鍵是其他問題（資料結構、組合數學，或者是列舉思想）若干道經典的離散化＋掃描線的題目
1.
離散化: 將所有的x軸座標存在一個數組裡..排序.當進入一條線段時..通過二分的方式確定其左右點對應的離散值...
掃描線..可以看成一根平行於x軸的直線..至y=0開始往上掃..直到掃出最後一條平行於x軸的邊..但是真正在做的時候..不需要完全模擬這個過程..掃描線的做法是從最下面的邊開始掃到最上面的邊.
線段樹: 本題用於動態維護掃描線在往上走時..x哪些區域是有合法面積的..
/*
POJ 1151 Atlantis
求矩形並的面積（線段樹+離散化）
*/
#include<stdio.h>
#include<iostream>
#include<algorithm>
using namespace std;
#define MAXN 201
struct Node
{
int l,r;//線段樹的左右整點
int c;//c用來記錄重疊情況
double cnt,lf,rf;//
//cnt用來計算實在的長度，rf,lf分別是對應的左右真實的浮點數端點
}segTree[MAXN*3];
struct Line
{
double x,y1,y2;
int f;
}line[MAXN];
//把一段段平行於y軸的線段表示成陣列 ，
//x是線段的x座標，y1,y2線段對應的下端點和上端點的座標
//一個矩形，左邊的那條邊f為1，右邊的為-1，
//用來記錄重疊情況，可以根據這個來計算，nod節點中的c
bool cmp(Line a,Line b)//sort排序的函式
{
return a.x < b.x;
}
double y[MAXN];//記錄y座標的陣列
void Buil