最開頭的話
與點,線,多邊形,圓形等各種幾何圖形相關的算法稱為計算幾何算法。
計算幾何已經稱為3D圖形,CAD,機器人等多中領域的計基礎,在程式設計競賽中也常常出現有關的問題。
計算幾何包含很多的内容,涉及範圍很廣。但是,在競賽中出現的題目大多集中于基礎性的問題上。而解決這些問題的關鍵是,把大學線性代數和高中幾何學的知識轉換成代碼。是以,這篇部落格的主要内容是如何将基礎數學理論轉化為簡潔而無異常的代碼形式。
計算幾何的工具
向量的實作
向量最直覺的表現形式為帶有方向的箭頭。向量有長度和方向兩個因素構成,是以箭頭的起點并不重要。是以,可以将向量的起點定義為坐标空間的原點,那麼就可以将向量表示成箭頭的終位置(x,y)。經過這樣的處理後,就能把向量變成帶有兩個成員變量的類。
以下的代碼實作了這種方法的vector2類函數。并且實作了多種運算符的重載。
const double PI = * acos();
struct vector2
{
double x, y;
//explicit 避免加入實數
explicit vector2(double x_, double y_):x(x_),y(y_){}
// 比較兩個向量
bool operator == (const vector2& rhs) const
{
return (x == rhs.x && y == rhs.y);
}
bool operator < (const vector2& rhs) const
{
return x != rhs.x ? x < rhs.x : y < rhs.y;
}
//向量之間的加減
vector2 operator + (const vector2& rhs) const
{
return vector2(x+rhs.x, y+rhs.y);
}
vector2 operator - (const vector2& rhs) const
{
return vector2(x-rhs.x, y-rhs.y);
}
//向量乘實數
vector2 operator * (double rhs) const
{
return vector2(x*rhs, y*rhs);
}
//傳回向量的長度
double norm(){return hypot(x, y);}
//傳回機關向量
vector2 normalize()
{
return vector2(x/norm(), y/norm());
}
//傳回從x軸正方向逆時針到達目前向量的角度
double polar() const
{
return fmod(atan2(y,x) + *PI, *PI);
}
//點乘
double dot(vector2& rhs) const
{
return x*rhs.x + y*rhs.y;
}
//叉乘
double cross(const vector2& rhs) const
{
return x * rhs.y - y * rhs.x;
}
};
其中可能存在大家不太了解的函數,麻煩自行百度。
點與直線、線段的表示方法
我們将把線段表示成以這兩個端點為終點的兩個向量。
将直線表示成包含于該直線的任意一條線段。
這表示的最直接的好處,是可以和已有的向量建立起聯系,可以使用現成的函數。
向量的內積和叉積
這兩者在計算幾何中很常用,都有着各自的幾何含義。
内積
定義
各種符号的,不好打,自行百度。
用處
- 向量的夾角
- 判斷兩個向量是否 相等
- 向量的投影
叉積
定義
略。。。
用處
-
計算面積
向量a和b的叉積絕對值等于将向量a和b用作兩邊的平行四邊形的面積。
-
判斷兩個向量的方向
兩個向量的叉積,如果為正數,就能判斷a在b的逆時針方向的180度以内;
如果為負數,則在順時針方向。
叉積的實作代碼:
//原點向量在向量a的逆時針上,傳回正數;順時針 傳回負數;平行 傳回 0;
double ccw(vector2 a, vector2 b)
{
return a.cross(b);
}
//把點p視為基準點時,傳回值情況同上;
double ccw(vector2 p, vector2 a, vector2 b)
{
return ccw(a-p, b-p);
}
相交、距離、面積
直線與直線相交
确認直線之間是否相交在幾何題中很常見,但是,編寫這樣的代碼卻不太容易,以為要考慮各種特殊的情況。
表示相交直線的最好方式為,将直線表示成一個點和一個方向向量。
<因為不知道向量符号咋打的, 是以 下文中兩個相同字元的視為向量符号加上單個字元>
例如:求解直線l1: aa + p * bb ; 直線 l2: cc + q * dd;
可得: aax + p * bbx = ccx + q * ddx; aay + p * bby = ccy + q * ddy;
是以:
p=(ccx−aax)ddy−(ccy−aay)ddybbx∗ddy−bby∗ddx
化簡得: p=(cc−aa)×ddbb×dd
最後将p的值 帶入直線l1 中即可。實作代碼:
//計算兩個直線(a,b)和直線(c,d)的交點
double lineIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& x)
{
double det = (b-a).cross(d-c);
if(fabs(det) < e-) return false;
x = a + (b-a)*((c-a).cross(d-c)/det);
return true;
}
線段和線段相交,不需要求交點
自己畫一下圖,很好了解的。
//線段(a,b)和線段(c,d) 是否有交點。
bool segmentIntersection(vector2 a, vector2 b, vector2 c, vector2 d)
{
double ab = ccw(a,b,c)*ccw(a,b,d);
double cd = ccw(c,d,a)*ccw(c,b,d);
if(ab == && cd == )
{
if(b < a) swap(a,b);
if(d < c) swap(d,c);
return !(b < c || d < a);
}
return ab <= && cd <= ;
}
凸包
先上張圖檔
這裡寫圖檔描述
先來了解一下,什麼是凸包。
點集Q的凸包(convex hull)是指一個最小凸多邊形,滿足Q中的點或者在多邊形邊上或者在其内。
有人把這些點看成在木闆上的釘子,而凸包就是一個橡皮筋套住所有釘子是的輪廓。
求解凸包的簡單方法是卷包裹法。
卷包裹算法從一個必然在凸包上的點開始向着一個方向依次選擇最外側的點
很容易了解對不對。上代碼:
//凸包(卷包裹法)
typedef vector<vector2> polygon;
polygon giftWrap(const vector<vector2>& points)
{
int n = points.size();
polygon hull;
//找出最左下角的點
vector2 pivot = *min_element(points.begin(), points.end());
hull.push_back(pivot);
while(true)
{
//找出最左邊的點next
vector2 ph = hull.back(),next = points[];
for(int i = ; i < n; i++)
{
double cross = ccw(ph, next, points[]);
double dist = (next - ph).norm() - (points[] - ph).norm();
if(cross > || (cross == && dist < ))
next = points[i];
}
if(next == pivot) break;
hull.push_back(next);
}
}