依然放上本文的github地址...

作業QwQ

先來說一下上次留下的例題.
poj這道題並沒有實數比較模式..
所以被精度勢力幹翻.
交上去WA掉竟然是因為-0.00和0.00不相等?
根據對拍結果別的地方應該沒什麽問題了OvO
下面給出並不能AC的"正確"代碼:

#include <cstdio>
#include <cmath>
const double eps=1e-8;
inline int dcmp(const double &a){
    if(fabs(a)<eps) return 0;
    return a<0?-1:1
 
;
}
int errr=0;
inline int gn(int a=0,char c=0,int f=1){
    for(;(c<48||c>57)&&c!='-';c=getchar());
    if(c=='-')f=-1,c=getchar();
    for(;c>47&&c<58;c=getchar())a=a*10+c-48;return a*f;
}
struct vec{
    double x,y;
    vec(double p=0,double q=0):x(p),y(q){}
};
double
 
 operator *(const vec &a,const vec &b){
    return a.x*b.y-a.y*b.x;
} //向量運算只需要叉乘
inline vec Cut(const vec &p0,const vec &p1,const vec &q0,const vec &q1){
    errr=0;
    double a1,b1,c1,a2,b2,c2,d;
    a1=p1.y-p0.y; b1=p0.x-p1.x; c1=p0*p1;
    a2=q1.y-q0.y; b2=q0.x-q1.x; c2=q0*q1;
    d=a1*b2-a2*b1;
    if
 
(!dcmp(d)){ //斜率相等
        if(!dcmp(c1/a1-c2/a2)) errr=1; //平行
        else errr=-1;   //重合
        return vec(0,0);
    }
    return vec((b2*c1-b1*c2)/d,(a1*c2-a2*c1)/d);
}
int main(){
    int T=gn(); vec a[4];
    puts("INTERSECTING LINES OUTPUT");
    while(T--){
        for(int i=0;i<4;++i) a[i].x=gn(),a[i].y=gn();
        vec s=Cut(a[0],a[1],a[2],a[3]);
        if(errr<0) puts("NONE"); else if(errr>0) puts("LINE");
        else printf("POINT %.2lf %.2lf\n",s.x,s.y);
    }
    puts("END OF OUTPUT");
}

新課OvO

這次學學多邊形...
多邊形就是多個點首尾順次相接圍成的圖形.

多邊形分為凹多邊形和凸多邊形.
大家見到的比較美觀的多邊形都是凸多邊形, 比如正方形一類的.
它們的特點是比較可愛 任取一條邊, 整個圖形都在這條邊所在直線的一側.

而凹多邊形則是不可愛滿足這一點的多邊形.

多邊形的形式非常的多, 我們通常用一個儲存了所有頂點坐標的數組來表示多邊形.
關於多邊形有如下常見的問題:

判斷點是否在多邊形內部

這是一個非常經典的問題, 導致出現了很多種解法.

射線法(交叉法)

以該點為起點, 做一條射線(一般方向取水平向右), 計算與多邊形交點的奇偶性.

可以看到 多邊形內部的點為起點做的射線與多邊形的交點個數是奇數, 而以外部的點為起點做的射線與多邊形的交點個數是偶數.
想法雖然很美好, 但是很顯然, 這樣做如果交點恰好是某個頂點就會出現問題.
這樣我們就需要新的定義.
一種標準的做法是, 在左側或者下邊的點屬於內部，在右側或者頂邊的點在外部.(就是計左不計右或計下不計上辣~)
這樣有公共的邊界線段的兩個多邊形就不能擁有公共點(起碼我們不視為這樣)
在這裏我們的交叉法有一個專用的邊緣交叉規則:

1. 一條向上的邊包括它的起始點, 不包括它的終止點.(記下不記上)
2. 一條向下的邊不包括它的起始點, 包括它的終止點.(記上不記下)
3. 不計所有水平的邊. (既不上也不下的特殊情況不予考慮)
4. 交點必須嚴格出現在P的右側. (右邊界屬於外部, 左邊界屬於內部)

舉個栗子, 交點是頂點的情況大約有下面幾種:

根據邊緣交叉規則, 綠色的點才是滿足條件的交點.
我們可以很容易的寫出算法流程:

1. 設點的坐標為\((x_0,y_0)\),
2. 遍歷多邊形的每一條邊e
3. 看線段e是否跨立直線\(y=y_0\).
4. 根據向上還是向下分別判斷是否符合#1 #2 #3
5. 求出e與直線\(y=y_0\)交點的橫坐標\(x\), 跟\(x_0\)比較, 如果?\(x>x_0\)(滿足#4), 則?\(ans++\).
6. 判斷\(ans\)的奇偶性, 奇數則在內部, 偶數則在外部.

一道本子……啊呸 板子題.. 稀有的hdu中文題目...
雖然題目描述並不是太嚴謹, 但是數據水啊, 但是題意很精煉.
剛才分析那麽多了直接貼代碼吧.(這個題1A了就很爽..

#include <cmath> 
#include <cstdio>
#include <cstring>
const double eps=1e-9;
int dcmp(const double &a){
    if(fabs(a)<eps) return 0;
    return a<0?-1:1;
}
struct point{
    double x,y;
}poly[202],pt;
double operator *(const point &a,const point &b){return a.x*b.y-a.y*b.x;}
//求p0p1兩點確定的直線和直線y=y0的交點 
double Cutx(const double &y0,const point &p0,const point &p1){
    return ((y0-p0.y)*p1.x+(p1.y-y0)*p0.x)/(p1.y-p0.y);
}
//交叉法判斷點是否在多邊形內
//參數p表示點 pts表示多邊形的頂點集 pcnt表示頂點數 
bool cn_ptInPoly(const point &p,point *pts,int pcnt){
    pts[pcnt]=pts[0]; int cnt=0;
    for(int i=1;i<=pcnt;++i)
        if(((dcmp(pts[i].y-p.y)>0&&dcmp(pts[i-1].y-p.y)<1)      //#1
            ||(dcmp(pts[i].y-p.y)<1&&dcmp(pts[i-1].y-p.y>0)))   //#2(#3)
            &&(dcmp(p.x-Cutx(p.y,pts[i],pts[i-1]))<0))          //#4
            ++cnt; //交點計數一次 
    return cnt&1; //交點奇數時在內, 偶數時在外.
}
int main(){
    int n,m;
    while(~scanf("%d",&n)){
        memset(poly,0,sizeof(poly));
        for(int i=0;i<n;++i) scanf("%lf%lf",&poly[i].x,&poly[i].y);
        scanf("%d",&m);
        for(int i=0;i<m;++i) {
            scanf("%lf%lf",&pt.x,&pt.y);
            puts(cn_ptInPoly(pt,poly,n)?"Yes":"No");
        } 
    }
}

這裏的話唯一沒有說的就是求交點了.
首先經歷了一波特判我們知道一定是有且只有一個交點的..
而且縱坐標知道了只需要求橫坐標, 就是解析幾何化式子的過程了...
我們這裏讓點\(p\)的坐標為\((x_p,y_p)\), 線段兩端點分別是\((x_0,y_0),(x_1,y_1)\).
因為是求交點所以日常聯立.. 根據必修二所學知識, 我們直線方程可以設兩點式.
\[ \left\{\begin{matrix} y=y_p\\ \frac{x_p-x_0}{x_1-x_0}=\frac{y-y_0}{y_1-y_0} \end{matrix}\right. \]
然後將①代入②, 隨便(???)化一波式子得到
\[ x_p=x=\frac{(y_p-y_0)x_1+(y_1-y_p)x_0}{y_1-y_0} \]
或者什麽其他式子也可以..只要能算就行..然後就可以了..

這個算法比較簡單好用, 而且速度也不錯. 時間復雜度是\(O(n)\)的.. 也比較常用.

繞數法

這個射線法已經是比較好用多了, 但是有一個缺陷..
多邊形還有一種分類是簡單和非簡單...
非簡單多邊形長得就比凹多邊形更為難受...
它們竟然會自我重疊...
就像這樣...(圖中的數字表示頂點的順序)

這個多邊形就很"不簡單".. 因為黃色的區域自我重疊了..
我知道這樣非常的不好看, 但是就是有這樣的圖形我也很絕望啊←_←
然後再用射線法的時候就會出現一些問題..

你看這樣一統計.. 交點個數是偶數 所以在外面?! 這樣就出問題了...
那我們該怎麽辦?? 辦法總比問題多, 我們可以再想有沒有別的方法嘛...
有一種神奇的操作, 我們可以計算多邊形繞某個點繞過的角度...

雖然畫的有點亂,但是還是可以跟著紫色箭頭看一下的..
先從這個點向每個頂點做射線, 然後從1繞一圈繞回1. 看轉過的角度是不是0.
這個點在外部, 當且僅當繞過的角度是0.(其實就是沒繞嘛, 這麽就好理解了)
很明顯圖中的角度是\(2\pi\)而不是0, 所以在內部.
在外部的點就是0了.像這樣:

可以看到, 黃色區域的點也可以被判斷到內部(不過繞了兩圈)

那我們如何統計呢?
因為都是向量的夾角, 我們又想到了點積. 我們最後要求的就是
\[ \sum_{i=0}^{n-1}arccos(\frac{\vec{PV_i}\vec{PV_{i+1}}}{|\vec{PV_1}||\vec{PV_{i+1}}|}) \]
但是這裏我們用到了一個昂貴的函數arccos, 這樣運行會變慢, 而且精度會受到影響.
或許最初這也就是繞數法不如射線法常用的原因?(據說可能會慢20倍?)
所以我們要考慮優化.
我也不知道誰想出了一種天才的優化方式..
這種方式借鑒了射線法, 依然做一條射線, 經過上行邊時, 計數器+1; 經過下行邊時, 計數器-1, 最後看計數器是否為0, 為0則點在多邊形外, 否則點在多邊形內.

舉個栗子:

圖上數得很清楚了 自己看看就能懂..
但是還有一些細節問題..
比如邊界問題.. 見上面的邊緣交叉規則即可..
在這裏我們沒有必要求出交點的坐標, 我們只需要判斷交點會不會在定點的右邊,
所以我們匯總一下一共有這幾種情況:

對於上行邊來說, 邊的起始點指向定點的向量要在邊的逆時針方向; 對於下行邊來說, 則要在順時針方向, 這樣我們一點積就出來了..
還是給出hdu1756的代碼(又是1A哈哈哈哈

#include <cmath> 
#include <cstdio>
#include <cstring>
const double eps=1e-9;
int dcmp(const double &a){
    if(fabs(a)<eps) return 0;
    return a<0?-1:1;
}
struct point{
    double x,y;
    point(double X=0,double Y=0):x(X),y(Y){}
}poly[202],pt;
point operator -(const point &a,const point &b){return point(a.x-b.x,a.y-b.y);}
double operator *(const point &a,const point &b){return a.x*b.y-a.y*b.x;}
//繞數法確定點是否在多邊形內 
bool wn_ptInPoly(const point &p,point *pts,int pcnt){
    pts[pcnt]=pts[0]; int cnt=0;
    for(int i=1;i<=pcnt;++i)
        if((dcmp(pts[i].y-p.y)>0&&dcmp(pts[i-1].y-p.y)<1)       //#1(#3) 
            &&dcmp((p-pts[i-1])*(pts[i]-pts[i-1]))<0)           //#4逆時針 
            ++cnt; 
        else if((dcmp(pts[i].y-p.y)<1&&dcmp(pts[i-1].y-p.y)>0)  //#2(#3) 
            &&dcmp((p-pts[i-1])*(pts[i]-pts[i-1]))>0)           //#4順時針
            --cnt;
    return cnt; //只要不為0 
}
int main(){
    int n,m;
    while(~scanf("%d",&n)){
        memset(poly,0,sizeof(poly));
        for(int i=0;i<n;++i) scanf("%lf%lf",&poly[i].x,&poly[i].y);
        scanf("%d",&m);
        for(int i=0;i<m;++i) {
            scanf("%lf%lf",&pt.x,&pt.y);
            puts(wn_ptInPoly(pt,poly,n)?"Yes":"No");
        } 
    }
}

不過竟然跑了15ms 不知道為什麽跟射線法的0ms還有一定差距..
所以簡單的多邊形還是跑射線法似乎靠譜一點...
不過出於正確性的考慮, 如果不保證多邊形簡單的話, 還是要用繞數的..

二分法

我們之前說過, 凸多邊形非常的美觀.
美觀的圖形就一定要有美觀的性質.
於是我們有一種特殊的判斷點在凸多邊形內的方法.

我們來看這張圖, 第一步, 我們可以做一個"快速排斥實驗".
(P.S. 以下均假設點按照逆時針方向給出, 順時針請自求多福做個對稱..)
判斷一下\(\vec {PV_0}\)與\(\vec{PV_1}\)和\(\vec{PV_8}\)的順逆時針關系.
如果在\(\angle V_1V_0V_8\)之外的話就不需要考慮了(因為是凸多邊形).
然後我們在\(1\sim8\)之間二分, 每次判斷順逆時針關系.
就可以在\(O(logn)\)時間內確定是在\(\angle V_iV_0V_{i+1}\)內了.
然後我們再\(O(1)\)判斷\(P\)點在\(\vec {V_iV_{i+1}}\)左側還是右側(用叉積)就好了~
(沒找到裸題TAT 那就只能幹寫了 正確性就不敢說了喲~)

//二分確定點是否在凸多邊形內
//P.S. 假設點以逆時針給出 
bool ptInConvexPoly(const point &p,point *pts,int pcnt){
    point po=p-pts[0]; //比較常用,所以存一下 
    if(dcmp(po*(pts[1]-pts[0]))>0
        ||dcmp(po*(pts[pcnt-1]-pts[0]))<0)
        return false;  //"快速排斥"一波
    int l=1,r=pcnt-1;
    while(r-l>1){
        int mid=(l+r)>>1;
        if(!dcmp(po*(pts[mid]-pts[0]))<0) //逆時針方向
            l=mid;
        else r=mid;
    }
    return dcmp((p-pts[l])*(pts[l+1]-pts[l]))<0; //在邊的左側 
} 

然後完整的代碼依然去github找好了...

最後發現自己其實一天只學了判斷點在多邊形內..
心情復雜...

【筆記篇】最良心的計算幾何學習筆記（二）