整理的演算法模板合集: ACM模板
目錄
- 1.平面掃描
- 2.凸包
- Andrew演算法
- 直線兩點式轉為一般式使用公式O(1)計算點到直線距離
- 判斷點是否在凸包內
- 3.旋轉卡殼
- 4.半平面交
- 有向直線
- O ( n l o g n ) O(nlogn) O(nlogn)的半平面交
- 二分 + 半平面交
- 5.閔可夫斯基和
- 6.平面區域
1.平面掃描
平面掃描是用來降低演算法的復雜度的方法,通過掃描線在平面上按照給定的軌跡移動的同時,不斷根據掃描線掃過的部分更新資訊,從而得到整體所要求的結果,既可以從左向右平移與y軸平行的直線,也可以固定射線的端點逆時針轉動,
平面上有N NN個兩兩沒有公共點的圓,i ii號圓的圓心在 ( x i , y i ) (x_i,y_i) (xi?,yi?),半徑為 r i r_i ri?,求所有最外層的,即不包含與其他圓內部的圓,
利用垂直于x軸的線去從左往右掃描,觀察點為圓的左右端點,因為圓是沒有公共點的,這使得包含的判斷變得簡單,如果某個圓被包含那么一定是被最近的兩個外圈的其中一個圓包含(y軸上),假設這個圓已經被包含,如果遠的圓如果想包含這個圓必須把兩個圓都包住,那么第一層包裹已經不再是外層,所以如果存在被包括在某個外層里,一定是被包裹在最近的兩個外層的其中一個,
typedef long long ll;
typedef pair<double, int>PDL;
const int N = 50007, M = 5000007, INF = 0x3f3f3f3f;
set<PDL>out;
vector<int>ans;
int n, m;
int cnt, num;
PDL p[N * 2];
double x[N], y[N], r[N];
bool inside(int i, int j)
{
double dx = x[i] - x[j];
double dy = y[i] - y[j];
return dx * dx + dy * dy <= r[j] * r[j];
}
int main()
{
scanf("%d", &n);
for(int i = 1; i <= n; ++ i){
scanf("%lf%lf%lf", &r[i], &x[i], &y[i]);
p[++num].first = x[i] - r[i];//左端點
p[num].second = i;
p[++num].first = x[i] + r[i];//右端點
p[num].second = i + n;
}
sort(p + 1, p + 1 + num);
for(int i = 1; i <= num; ++ i){
int id = p[i].second;
if(id <= n){//左端點
set<PDL> :: iterator it = out.lower_bound(make_pair(y[id], id));
if(it != out.end() && inside(id, it->second))continue;//被包含,不是答案,跳過
if(it != out.begin() && inside(id, (--it)->second))continue;
ans.push_back(id);
out.insert(make_pair(y[id], id));//外圈
}
else {//到了這個圓的右端點該洗掉了,已經沒用了,留著會影響答案
id -= n;
out.erase(make_pair(y[id], id));
}
}
printf("%d\n", ans.size());
sort(ans.begin(),ans.end());
for(int i = 0; i < ans.size(); ++i)
printf("%d ", ans[i]);
return 0;
}
2.凸包
過多邊形的任意一邊做一條直線,如果其他各個頂點都在這條直線的同側,則把這個多邊形叫做凸多邊形,
凸包求解演算法的基礎便是凸多邊形的定義與性質,
假設平面上n個點,過某些點作一個多邊形,使這個多邊形能把所有點都“包”起來,當這個多邊形是凸多邊形的時候,我們就叫它“凸包”,
假設每種顏料都擁有(R,G,B)三種屬性,表示該種顏料紅色,綠色,與藍色的化學成分所占的比重
給你若干種已有的不限量的顏料,問是否能夠勾兌出目標顏料(R0,G0,B0)
將每一種顏料映射為二維歐氏空間中的一個點,我們可以將已經給定的顏料與目的顏料在空間中標定出來
經過觀察與思考,我們可以發現,一個顏料能夠被勾兌出來當且僅當該顏料對應的點,位于以給定顏料所構成的凸包之中
Andrew演算法
判斷這樣連是否為凸多邊形的一部分:
如果p2-p3的斜率小于p1-p3,那么我們就沒必要選擇p2這個點,而是直接走p1~p3即可,也就是說我們把已經放進去的p2踢掉
//計算凸包,輸入點陣列p,個數為p輸出點陣列ch,函式回傳凸包頂點個數,
//輸入不能有重復的點,函式執行完后的輸入點的順序將被破壞(因為要排序,可以加一個陣列存原來的id)
//如果不希望在凸包邊上有輸入點,把兩個<=改成<即可
int ConvexHull(Point* p, int n, Point* ch)
{
sort(p, p + n);
int m = 0;
for(int i = 0; i < n; ++ i){//下凸包
//如果叉積<=0說明新邊斜率小說明已經不是凸包邊了,趕緊踢走
while(m > 1 && Cross(ch[m - 1] - ch[m - 2], p[i] - ch[m - 2]) <= 0)m -- ;
ch[m ++ ] = p[i];
}
int k = m;
for(int i = n - 2; i >= 0; -- i){//上凸包
while(m > k && Cross(ch[m - 1] - ch[m - 2], p[i] - ch[m - 2]) <= 0)m -- ;
ch[m ++ ] = p[i];
}
if(n > 1) m -- ;
return m;
}
直線兩點式轉為一般式使用公式O(1)計算點到直線距離
平面上有n個點,求一條直線,使得這n個點都在這條直線上或同側,且每個點到該直線的距離之和盡量小,
首先,一條直線不分割這n個點當且僅當不分割這n個點的凸包,并且為了使這n個點到該直線的距離最小,這條直線應是凸包上某一條邊所在直線,
由幾何知識可得,對于點
P
(
x
0
,
y
0
)
P(x_0,y_0)
P(x0?,y0?)和直線
A
x
+
B
y
+
C
=
0
Ax+By+C=0
Ax+By+C=0之間的距離
d i s = ∣ A x 0 + B y 0 + C ∣ A 2 + B 2 dis=\frac{|Ax_0+By_0+C|}{\sqrt{A^2+B^2}} dis=A2+B2 ?∣Ax0?+By0?+C∣?
又因為這n個點在這條直線的同側,所以對于所有的點 A x 0 + B y 0 + C Ax_0+By_0+C Ax0?+By0?+C符號相同,所以預處理出所有點x坐標和y坐標的和就可以 Θ ( 1 ) \Theta(1) Θ(1)算出每個點到該直線的距離之和,
還有一個問題,如何將直線的兩點式
l
P
Q
,
P
(
x
P
,
y
P
)
,
Q
(
x
Q
,
y
Q
)
l_{PQ},P(x_P,y_P),Q(x_Q,y_Q)
lPQ?,P(xP?,yP?),Q(xQ?,yQ?)轉化為一般式
A
x
+
B
y
+
C
=
0
Ax+By+C=0
Ax+By+C=0
簡單計算可得一個解為
A
=
y
Q
?
y
P
A=y_Q-y_P
A=yQ??yP?
B
=
x
P
?
x
Q
B=x_P-x_Q
B=xP??xQ?
C
=
?
x
P
?
A
?
y
P
?
B
C=-x_P * A - y_P * B
C=?xP??A?yP??B
注意千萬不能用除法,否則會計算結果會特別大或NaN導致WA.
//Ax + By + C = 0;
double get(double A, double B, double C) {
double k = fabs(A*sumx + B*sumy + n*C);
double v = sqrt(A*A + B*B);//v != 0;
return k/v;
}
double getDist(const point &a, const point &b) {
double A = a.y-b.y;
double B = b.x-a.x;
double C = a.x*b.y - a.y*b.x;
return get(A, B, C);
}
int main()
{
int counter = 0;
int T;
point t;
scanf("%d", &T);
while(T--) {
init();
scanf("%d", &n);
for(int i = 0; i < n; i++) {
scanf("%lf%lf", &t.x, &t.y);
sumx += t.x; sumy += t.y;
tmp.push_back(t);
}
polygon_convex tres = convex_hull(tmp);
int Size = (int)tres.P.size();
printf("Case #%d: ", ++counter);
if(Size == 2 || Size == 1) { //剛開始wa一次,看了看題目n>0.又把n=1考慮一下,
printf("0.000\n");
continue;
}
for(int i = 0; i < Size; i++) {
double temp = getDist(tres.P[i], tres.P[(i+1)%Size]);
ans = min(ans, temp);
}
printf("%.3lf\n", ans/n);
}
return 0;
}
判斷點是否在凸包內
先判定和邊界的關系
然后找到與其極角相鄰的兩點,憑此判斷
須保證A[1]=(0,0)
ll in(Node a)
{
if(a*A[1]>0||A[tot]*a>0) return 0;
ll ps=lower_bound(A+1,A+tot+1,a,cmp2)-A-1;
return (a-A[ps])*(A[ps%tot+1]-A[ps])<=0;
}
3.旋轉卡殼
利用凸包上一些奇妙的單調性,求解
- 多邊形直徑
- 多邊形寬度
- 最小面積矩形覆寫
復雜度 O ( n ) O(n) O(n)
for(int i=1,p=1;i<=n;i++)
{
//這份代碼只卡了一個點,還可以同時卡多個點
while(H(A[i],A[i+1],A[p%n+1])>=H(A[i],A[i+1],A[p])) p=p%n+1;
Ans=max(Ans,max((A[p]-A[i]).dis(),(A[p]-A[i+1]).dis()));
}
給定平面上 n 個點,求凸包直徑,
第一行一個正整數 n,接下來 n 行,每行兩個整數 x,y,表示一個點的坐標,輸出一行一個整數,表示答案的平方,
#include<bits/stdc++.h>
using namespace std;
const int N = 50007, M = 50007, INF = 0x3f3f3f3f;
const double DINF = 12345678910, eps = 1e-10, PI = acos(-1);
struct Point{
int 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 - (Point A, Point 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);}
bool operator < (const Point& a, Point& b) {return a.x < b.x || (a.x == b.x && a.y < b. y);}
int dcmp(double x){
if(fabs(x) < eps) return 0;
else return x < 0 ? -1 : 1;
}
bool operator == (const Point& a, const Point& b){return !dcmp(a.x - b.x) && !dcmp(a.y - b.y);}
double Polar_angle(Vector A){return atan2(A.y, A.x);}
inline double D_to_R(double D)//角度轉弧度
{ return PI/180*D; }
double Cross(Vector A, Vector B){return A.x * B.y - B.x * A.y;}
Vector Rotate(Vector A, double rad){
return Vector(A.x * cos(rad) - A.y * sin(rad), A.x * sin(rad) + A.y * cos(rad));
}
Point Get_line_intersection(Point P,Vector v,Point Q,Vector w)
{
Vector u = P - Q;
double t = Cross(w, u) / Cross(v, w);
return P + v * t;
}
double convex_polygon_area(Point* p, int n)
{
double area = 0;
for(int i = 1; i <= n - 2; ++ i)
area += Cross(p[i] - p[0], p[i + 1] - p[0]);
return area / 2;
}
int ConvexHull(Point* p, int n, Point* ch)
{
sort(p, p + n);
int m = 0;
for(int i = 0; i < n; ++ i){//下凸包
while(m > 1 && Cross(ch[m - 1] - ch[m - 2], p[i] - ch[m - 2]) <= 0)m -- ;
ch[m ++ ] = p[i];
}
int k = m;
for(int i = n - 2; i >= 0; -- i){//上凸包
while(m > k && Cross(ch[m - 1] - ch[m - 2], p[i] - ch[m - 2]) <= 0)m -- ;
ch[m ++ ] = p[i];
}
if(n > 1) m -- ;
return m;
}
int get_dist (const Point &x){return x.x * x.x + x.y * x.y;}
Point p[N], con[N];
int con_num;
double Rotating_calipers()
{
int op = 1, ans = 0;
for(int i = 0; i < con_num; ++ i){
while(Cross((con[i] - con[op]), (con[i + 1] - con[i])) < Cross((con[i] - con[op + 1]), (con[i + 1] - con[i])))
//(寫成<=會被兩個點的資料卡掉,所以必須寫成<)
op = (op + 1) % con_num;
ans = max(ans, max(get_dist(con[i] - con[op]), get_dist(con[i + 1] - con[op])));
}
printf("%d\n", ans);
return ans;
}
int n ;
int main()
{
scanf("%d", &n);
for(int i = 0; i < n; ++ i){
scanf("%d%d", &p[i].x, &p[i].y);
}
con_num = ConvexHull(p, n, con);
double res = Rotating_calipers();
return 0;
}
4.半平面交
定義一個半平面為一個向量的左側,半平面交即為若干個半平面的交集,
有向直線
//有向直線,它的左邊就是對應的半平面
struct Line
{
Point P;//直線上任意一點
Vector v;//方向向量,左邊就是對應的半平面
double deg;//極角
Line(){}
Line(Point P, Vector v):P(P), v(v){deg = atan2(v.y, v.x);}
bool operator < (const Line& L)const {//排序時使用的比較運算子
return deg < L.deg;
}
};
O ( n l o g n ) O(nlogn) O(nlogn)的半平面交
//半平面交一般是一個凸多邊形,但是有時候會是一個無界多邊形
//甚至會是一個直線、線段、點,但是結果一定是凸的
//點p在有向直線L的左邊(線上的不算)(叉積大于0a在b的左側,小于0在右側[sin夾角])
bool on_left(Line L, Point P){return Cross(L.v, P - L.P) > 0;}
//兩個有向直線的交點/假定交點唯一存在
Point get_intersection(Line a, Line b)
{
Vector u = a.P - b.P;
double t = Cross(b.v, u) / Cross(a.v, b.v);
return a.P + a.v * t;
}
//半平面交的主程序
//L陣列存所有有向直線,n為有向直線個數,半平面交的交點存在poly陣列中
int half_plane_intersection(Line* L, int n, Point* poly)
{
sort(L, L + n);//按照極角排序
int first, last;//雙端佇列
Point *p = new Point[n];//p[i]為q[i]和q[i + 1]的交點
Line *q = new Line[n];//手寫的Line型別的雙端佇列(陣列)
q[first = last = 0] = L[0];//雙端佇列初始化的時候只有一個半平面L[0]
for(int i = 1; i < n; ++ i){
while(first < last && !on_left(L[i], p[last - 1]))last -- ;
while(first < last && !on_left(L[i], p[first]))first ++ ;
q[++ last] = L[i];//新的點是一定要放進去的
if(fabs(Cross(q[last].v, q[last - 1].v)) < eps){
//相鄰的兩個向量平行且同向,取內側的那一個
last -- ;//如果新的向量上的一個點在老的向量的左側就取新的
if(on_left(q[last], L[i].P))q[last] = L[i];
}
if(first < last)p[last - 1] = get_intersection(q[last - 1], q[last]);
}
while(first < last && !on_left(q[first], p[last - 1]))last -- ;
//洗掉無用的平面
if(last - first <= 1)return 0;//空集
p[last] = get_intersection(q[last], q[first]);//計算首尾兩個半平面交(環狀)
//從手寫deque中復制答案到輸出陣列中
int m = 0;
for(int i = first ; i <= last; ++ i)poly[m ++ ] = p[i];
return m;
}
逆時針給出n個凸多邊形的頂點坐標,求它們交的面積,
二分 + 半平面交
在大海的中央,有一個凸n邊形的小島,求出島上離海邊最遠的一個點到海邊的距離,保留6位小數,
//有向直線,它的左邊就是對應的半平面
struct Line
{
Point P;//直線上任意一點
Vector v;//方向向量,左邊就是對應的半平面
double deg;//極角
Line(){}
Line(Point P, Vector v):P(P), v(v){deg = atan2(v.y, v.x);}
bool operator < (const Line& L)const {//排序時使用的比較運算子
return deg < L.deg;
}
};
bool on_left(Line L, Point P){return Cross(L.v, P - L.P) > 0;}
//兩個有向直線的交點/假定交點唯一存在
Point get_intersection(Line a, Line b)
{
Vector u = a.P - b.P;
double t = Cross(b.v, u) / Cross(a.v, b.v);
return a.P + a.v * t;
}
//半平面交的主程序
int half_plane_intersection(Line* L, int n, Point* poly)
{
sort(L, L + n);//按照極角排序
int first, last;//雙端佇列
Point *p = new Point[n];//p[i]為q[i]和q[i + 1]的交點
Line *q = new Line[n];//手寫的Line型別的雙端佇列(陣列)
q[first = last = 0] = L[0];//雙端佇列初始化的時候只有一個半平面L[0]
for(int i = 1; i < n; ++ i){
while(first < last && !on_left(L[i], p[last - 1]))last -- ;
while(first < last && !on_left(L[i], p[first]))first ++ ;
q[++ last] = L[i];//新的點是一定要放進去的
if(fabs(Cross(q[last].v, q[last - 1].v)) < eps){
//相鄰的兩個向量平行且同向,取內側的那一個
last -- ;//如果新的向量上的一個點在老的向量的左側就取新的
if(on_left(q[last], L[i].P))q[last] = L[i];
}
if(first < last)p[last - 1] = get_intersection(q[last - 1], q[last]);
}
while(first < last && !on_left(q[first], p[last - 1]))last -- ;
//洗掉無用的平面
if(last - first <= 1)return 0;//空集
p[last] = get_intersection(q[last], q[first]);//計算首尾兩個半平面交(環狀)
//從手寫deque中復制答案到輸出陣列中
int m = 0;
for(int i = first ; i <= last; ++ i)poly[m ++ ] = p[i];
return m;
}
Point p[N], poly[N];
Line L[N];
Vector v[M], v2[M];
int n;
int main()
{
while(scanf("%d", &n) != EOF && n)
{
int m, x, y;
for(int i = 0; i < n; ++i){
scanf("%d%d", &x, &y);
p[i] = Point(x, y);
}
for(int i = 0; i < n; ++ i){
v[i] = p[(i + 1) % n] - p[i];
v2[i] = Normal(v[i]);//單位法向量
}
double l = 0, r = 20000;
while(r - l > eps){
double mid = l + (r - l) / 2;
for(int i = 0; i < n; ++ i)
L[i] = Line(p[i] + v2[i] * mid, v[i]);
//收縮/放大整個凸多邊形
//點+單位向量*長度=平移后的點
m = half_plane_intersection(L, n, poly);
if(!m)r = mid;
else l = mid;
}
printf("%.6f\n", l);
}
return 0;
}
5.閔可夫斯基和
定義p+q=(p.x+q.x,p.y+q.y),給定兩個點集,求{pi+qj}的凸包(凸殼)的問題
Vector V1[N],V2[N];
inline int Mincowski(Point *P1,Re n,Point *P2,Re m,Vector *V){//【閔可夫斯基和】求兩個凸包{P1},{P2}的向量集合{V}={P1+P2}構成的凸包
for(Re i=1;i<=n;++i)V1[i]=P1[i<n?i+1:1]-P1[i];
for(Re i=1;i<=m;++i)V2[i]=P2[i<m?i+1:1]-P2[i];
Re t=0,i=1,j=1;V[++t]=P1[1]+P2[1];
while(i<=n&&j<=m)++t,V[t]=V[t-1]+(dcmp(Cro(V1[i],V2[j]))>0?V1[i++]:V2[j++]);
while(i<=n)++t,V[t]=V[t-1]+V1[i++];
while(j<=m)++t,V[t]=V[t-1]+V2[j++];
return t;
}
6.平面區域
當平面上有很多線段時,組成的圖形往往不止一個多邊形,而是一個平面直線圖(PSLG),它代表一個平面區域劃分,其中一個區域是一個多邊形,
如果只有點和邊的資訊,如何找出所有區域呢?為方便起見,我們把每條邊u-v拆成兩條半邊 u-v 和 v-u.并且每條半邊只與它左邊的面相鄰,接下來,我們從一條半邊出發遍歷,每次像卷包裹演算法那樣找一個逆時針轉的盡量多的邊作為下一條邊,直到回到出發的那條半邊,
程式實作上可以把邊表擴大一倍,讓編號為i的半邊的反向邊的編號為i^1.
首先看一下邊的存盤結構定義:
//邊
struct Edge
{
int from, to; //起點,終點
double ang; //極角
};
平面直線圖:
//平面直線圖,
struct PSLG
{
int n, m; //頂點數 邊數
int face_cnt; //面數
double x[maxn]; //頂點
double y[maxn];
vector<Edge> edges;
vector<int> G[maxn];
bool vis[maxn*2]; // 每條邊是否被訪問過
int left[maxn*2]; //左面的編號
int prev[maxn*2]; //相同起點的上一條邊,即順時針旋轉碰到的下一條邊的編號
vector<Polygon> faces; //面
double areas[maxn]; //每個polygon的面積
void init(int n)
{
this->n = n;
edges.clear();
faces.clear();
for(int i = 0; i < n; i++)
G[i].clear();
}
//有向線段from-to的極角
double getAngle(int from, int to)
{
return atan2(y[to]-y[from], x[to]-x[from]);
}
void AddEdge(int from, int to)
{
edges.push_back((Edge){from, to, getAngle(from, to)});
edges.push_back((Edge){to, from, getAngle(to, from)});
m = edges.size();
G[from].push_back(m-2);
G[to].push_back(m-1);
}
//找出faces并計算面積
void Build()
{
int u, i, j, d, from, e;
//給從u開始個各條邊按照極角排序
for(u = 0; u < n; u++)
{
d = G[u].size();
for(i = 0; i < d; i++)
for(j = i+1; j < d; j++) //這里假設從u出發的線段不會太多
if(edges[G[u][i]].ang > edges[G[u][j]].ang)
swap(G[u][i], G[u][j]);
for(i = 0; i < d; i++)
prev[G[u][(i+1)%d]] = G[u][i];
}
face_cnt = 0;
memset(vis, false, sizeof(vis));
for(u = 0; u < n; u++)
for(i = 0; i < G[u].size(); i++)
{
e = G[u][i];
if(!vis[e]) //逆時針找圈
{
face_cnt++; //找到一個新的面
Polygon poly;
for(;;)
{
vis[e] = true;
left[e] = face_cnt;
from = edges[e].from;
poly.push_back(Point(x[from], y[from])); //把新的點加入面中
e = prev[e ^ 1]; //e^1為反向邊,然后prev就是需要走的下一條邊
if(e == G[u][i]) //回到了起始邊
break;
}
faces.push_back(poly);
}
}
//對于連通圖,惟一一個面積小于零的面是無限面
//對于內部區域來說,無限面多邊形的各個頂點是順時針的
for(i = 0; i < faces.size(); i++) //計算各個面的面積
areas[i] = PolygonArea(faces[i]);
}
};
轉載請註明出處,本文鏈接:https://www.uj5u.com/qita/155703.html
標籤:其他
