凸包算法（高级算法设计与分析实验1）

凸包问题

求解凸包问题：输入是平面上n个点的集合Q，凸包问题是要输出一个Q的凸包。其中，Q的凸包是一个凸多边形P，Q中的点或者在P上或者在P中。

实现基于枚举方法的凸包求解算法

提示：考虑Q中的任意四个点A、B、C、D，如果A处于BCD构成的三角形内部，那么A一定不属于凸包P的顶点集合。

这一方法属于暴力解法，任意枚举点集中的四个点，如果有一个点在其他三个点构成的三角形内部，则将这个点从点集中剔除。实验主要是通过python来实现的，先定义一个平面点类。代码自动生成二维平面[0,0]到[100,100]的数据集，四重循环逐一进行判断。列举关键代码如下：

class Point:
    def __init__(self, x, y):
        self.x = x
        self.y = y  
        self.angle=0

from Point import Point
import random
import matplotlib.pyplot as plt

#判断p4是否在p1,p2,p3构成的三角形中
def istriangle(a,b,c,p):

    signOfTrig = (b.x - a.x)*(c.y - a.y) - (b.y - a.y)*(c.x - a.x);
    signOfAB = (b.x - a.x)*(p.y - a.y) - (b.y - a.y)*(p.x - a.x);
    signOfCA = (a.x - c.x)*(p.y - c.y) - (a.y - c.y)*(p.x - c.x);
    signOfBC = (c.x - b.x)*(p.y - b.y) - (c.y - b.y)*(p.x - b.x);

    d1 = (signOfAB * signOfTrig > 0)
    d2 = (signOfCA * signOfTrig > 0)
    d3 = (signOfBC * signOfTrig > 0)

    return d1 and d2 and d3

def generate_data(size,data,data_o):
    for i in range(size):
        #x=random.randint(0,100)
        #y=random.randint(0,100)
        x=random.uniform(0,100)
        y=random.uniform(0,100)
        p=Point(x,y)
        data_o.append(p)
    return data,data_o

def d(A,B,P):
     ans=(B.x - A.x)*(P.y - A.y) - (B.y - A.y)*(P.x - A.x)
     return ans

def BruteForceCH(size,data,data_o):
    
    L=[]
    R=[]
    flag = [0 for x in range(size)]

    for i in range(size-3):
        if flag[i]==1:
            continue
        for j in range(i+1,size-2):
            if flag[j]==1:
                continue
            for k in range(j+1,size-1):
                if flag[k]==1:
                    continue
                for l in range(k+1,size):
                    if flag[l]==1:
                        continue
                    elif istriangle(data_o[i],data_o[j],data_o[k],data_o[l]):
                        flag[l]=1
                        continue
                    elif istriangle(data_o[l],data_o[j],data_o[k],data_o[i]):
                        flag[i]=1
                        continue
                    elif istriangle(data_o[i],data_o[l],data_o[k],data_o[j]):
                        flag[j]=1
                        continue
                    elif istriangle(data_o[i],data_o[j],data_o[l],data_o[k]):
                        flag[k]=1
    for i in range(size):
        if flag[i]==0:
            data.append(data_o[i])
    data=sorted(data,key = lambda point: (point.x,point.y))

    A=data[0]
    B=data[-1]
    del data[0]
    del data[-1]
    for P in data:
        if(d(A,B,P)>0):
            L.append(P)
        elif(d(A,B,P)<0):
            R.append(P)
    Lr=L.reverse()
    
   # '''
    for p in data_o:
        plt.scatter(p.x, p.y, c='g', marker='.')

    plt.scatter(A.x, A.y, c='r', marker='.')
    plt.plot([A.x,R[0].x],[A.y,R[0].y], color='r')
    for i in range(len(R)-1):
        plt.scatter(R[i].x, R[i].y, c='r', marker='.')
        plt.plot([R[i].x,R[i+1].x],[R[i].y,R[i+1].y], color='r')
    plt.scatter(R[-1].x, R[-1].y, c='r', marker='.')
    plt.plot([R[-1].x,B.x],[R[-1].y,B.y], color='r')
    plt.scatter(B.x, B.y, c='r', marker='.')
    plt.plot([B.x,L[0].x],[B.y,L[0].y], color='r')
    for i in range(len(L)-1):
        plt.scatter(L[i].x, L[i].y, c='r', marker='.')
        plt.plot([L[i].x,L[i+1].x],[L[i].y,L[i+1].y], color='r')
    plt.scatter(L[-1].x, L[-1].y, c='r', marker='.')
    plt.plot([L[-1].x,A.x],[L[-1].y,A.y], color='r')
    plt.show()
   # '''

实现基于Graham-Scan的凸包求解算法

具体实现方法就是，先找到点集中纵坐标最小的点，再此基础上再找横坐标最小的点，把这个点作为基点，这里需要遍历一遍点集，需要$O(n)$的时间复杂度。之后再遍历一遍点集，计算每个点与基点的逆时针角度，按照逆时针角度对点集进行排序，这里需要$O(nlogn)$的时间开销。之后依次按照逆时针方向计算是否为凸包。

from Point import Point
import random
import matplotlib.pyplot as plt
import math

def generate_data(size,data):
    for i in range(size):
        #x=random.randint(0,100)
        #y=random.randint(0,100)
        x=random.uniform(0,100)
        y=random.uniform(0,100)
        p=Point(x,y)
        data.append(p)
    return data
def judge(a,b,c):
    ans=(c.x - a.x)*(b.y - a.y) - (c.y - a.y)*(b.x - a.x)
    if ans>=0:
        return True
    else:
        return False
def GrahamScanCH(size,data):
    p0=Point(101,101)
    for p in data:
        if p.y<p0.y:
            p0=p
        elif p.y==p0.y and p.x<p0.x:
            p0=p
    data.remove(p0)
    for p in data:
        if p.x==p0.x:
            p.angle=math.pi/2
        elif p.x>p0.x:
            p.angle=math.atan((p.y-p0.y)/(p.x-p0.x))
        else:
            p.angle=math.atan((p0.x-p.x)/(p.y-p0.y))+math.pi/2
    data=sorted(data,key = lambda point: (point.angle))
    stack=[]
    stack.append(p0)
    stack.append(data[0])
    stack.append(data[1])
    for i in range(2,len(data)):
        while judge(stack[-2],stack[-1],data[i]):
            stack.pop(-1)
        stack.append(data[i])
    #'''
    plt.scatter(p0.x,p0.y,c='g',marker='.')
    for p in data:
        plt.scatter(p.x, p.y, c='g', marker='.')
    for i in range(len(stack)-1):
        plt.scatter(stack[i].x, stack[i].y, c='r', marker='.')
        plt.plot([stack[i].x,stack[i+1].x],[stack[i].y,stack[i+1].y], color='r')
    plt.scatter(stack[-1].x, stack[-1].y, c='r', marker='.')
    plt.plot([stack[-1].x,p0.x],[stack[-1].y,p0.y], color='r')
    plt.show()
    #'''

实现基于分治思想的凸包求解算法

分治算法主要是解决如何划分和合并的问题，划分采用方法是选择点集中横坐标的中位数，小于中位数的点集和大于中位数的点集分别进行子问题的计算。划分至点集小于等于三个点。合并的过程是基于Graham-Scan算法。

from Point import Point
import random
import matplotlib.pyplot as plt
import math

def generate_data(size,data):
    for i in range(size):
        #x=random.randint(0,10)
        #y=random.randint(0,10)
        x=random.uniform(0,100)
        y=random.uniform(0,100)
        p=Point(x,y)
        data.append(p)
    return data

def judge(a,b,c):
    ans=(c.x - a.x)*(b.y - a.y) - (c.y - a.y)*(b.x - a.x)
    if ans>=0:
        return True
    else:
        return False

def divide(data):
    maxn=-1
    minn=101
    data1=[]
    data2=[]
    for p in data:
        if p.x<minn:
            minn=p.x
        if p.x>maxn:
            maxn=p.x
    minddle=(maxn+minn)/2
    for p in data:
        if p.x<=minddle:
            data1.append(p)
        else:
            data2.append(p)
    return data1,data2


def conquer(data1,data2):
    if len(data1)==0:
        return data2
    if len(data2)==0:
        return data1
    #求纵坐标最小的点
    if data1[0].y<=data2[0].y:
        p0=data1[0]
        data1.remove(p0)
        flag=1
    else:
        p0=data2[0]
        data2.remove(p0)
        flag=0
    #计算每一个点相对于p0的极角
    p0.angle=0
    for p in data1:
        if p.x==p0.x:
            p.angle=math.pi/2
        elif p.x>p0.x:
            p.angle=math.atan((p.y-p0.y)/(p.x-p0.x))
        else:
            p.angle=math.atan((p0.x-p.x)/(p.y-p0.y))+math.pi/2
    for p in data2:
        if p.x==p0.x:
            p.angle=math.pi/2
        elif p.x>p0.x:
            p.angle=math.atan((p.y-p0.y)/(p.x-p0.x))
        else:
            p.angle=math.atan((p0.x-p.x)/(p.y-p0.y))+math.pi/2
    #利用data1和data2逆时针排序的特点，在线性时间内将所有点按极角大小排好序
    min1=10
    max1=-1
    min2=10
    max2=-1
    if flag==0:
        minindex=-1
        maxindex=-1
        for i in range(len(data1)):
            if(data1[i].angle>max1):
               max1=data1[i].angle
               maxindex=i
            if(data1[i].angle<min1):
               min1=data1[i].angle
               minindex=i
        if minindex<=maxindex:
            data3=data1[minindex:maxindex+1]
            data4=[]
            if minindex!=0:
                data4=data4+data1[minindex-1::-1]
            if maxindex!=len(data1)-1:
                data4=data4+data1[-1:maxindex:-1]
        else:
            data3=data1[minindex:maxindex:-1]
            data4=[]
            data4=data4+data1[minindex+1:]
            data4=data4+data1[0:maxindex+1]
        data5=[]
        i=0
        j=0
        while i<len(data3) and j<len(data4):
            if data3[i].angle<=data4[j].angle:
                data5.append(data3[i])
                i=i+1
            else:
                data5.append(data4[j])
                j=j+1
        while i<len(data3):
            data5.append(data3[i])
            i=i+1
        while j<len(data4):
            data5.append(data4[j])
            j=j+1
        data6=[]
        i=0
        j=0
        while i<len(data2) and j<len(data5):
            if data2[i].angle<=data5[j].angle:
                data6.append(data2[i])
                i=i+1
            else:
                data6.append(data5[j])
                j=j+1
        while i<len(data2):
            data6.append(data2[i])
            i=i+1
        while j<len(data5):
            data6.append(data5[j])
            j=j+1
    else:
        maxindex=-1
        minindex=-1
        for i in range(len(data2)):
            if(data2[i].angle>max2):
               max2=data2[i].angle
               maxindex=i
            if(data2[i].angle<min2):
               min2=data2[i].angle
               minindex=i
        if minindex<=maxindex:
            data3=data2[minindex:maxindex+1]
            data4=[]
            if minindex!=0:
                data4=data4+data2[minindex-1::-1]
            if maxindex!=len(data2)-1:
                data4=data4+data2[-1:maxindex:-1]
        else:
            data3=data2[minindex:maxindex:-1]
            data4=[]
            data4=data4+data2[minindex+1:]
            data4=data4+data2[0:maxindex+1]
        data5=[]

        i=0
        j=0
        while i<len(data3) and j<len(data4):
            if data3[i].angle<=data4[j].angle:
                data5.append(data3[i])
                i=i+1
            else:
                data5.append(data4[j])
                j=j+1
        while i<len(data3):
            data5.append(data3[i])
            i=i+1
        while j<len(data4):
            data5.append(data4[j])
            j=j+1
        data6=[]
        i=0
        j=0
        while i<len(data1) and j<len(data5):
            if data1[i].angle<=data5[j].angle:
                data6.append(data1[i])
                i=i+1
            else:
                data6.append(data5[j])
                j=j+1
        while i<len(data1):
            data6.append(data1[i])
            i=i+1
        while j<len(data5):
            data6.append(data5[j])
            j=j+1
    #洗牌，生成新的凸包
    if len(data6)<=1:
        return [p0]+data6
    stack=[]
    stack.append(p0)
    stack.append(data6[0])
    stack.append(data6[1])
    if len(data6)==2:
        return stack
    for i in range(2,len(data6)):
        while judge(stack[-2],stack[-1],data6[i]):
            stack.pop(-1)
        stack.append(data6[i])
    return stack

def divide_and_conquer(data):
    if len(data)<=1:
        return data
    if len(data)==2:
        return sorted(data,key = lambda point: (point.y,point.x))
    elif len(data)==3:
        p0=Point(101,101)
        for p in data:
            if p.y<p0.y:
                p0=p
            elif p.y==p0.y and p.x<p0.x:
                p0=p
        data.remove(p0)
        for p in data:
            if p.x==p0.x:
                p.angle=math.pi/2
            elif p.x>p0.x:
                p.angle=math.atan((p.y-p0.y)/(p.x-p0.x))
            else:
                p.angle=math.atan((p0.x-p.x)/(p.y-p0.y))+math.pi/2
        res=[]
        res.append(p0)
        if data[0].angle<=data[1].angle:
            res.append(data[0])
            res.append(data[1])
        else:
            res.append(data[1])
            res.append(data[0])
        return res
    data1,data2=divide(data)
    data1=divide_and_conquer(data1)
    data2=divide_and_conquer(data2)
    data=conquer(data1,data2)
    return data

def DivideCH(size,data):

    
    
    #plt.show()
    stack=divide_and_conquer(data)
   
    #'''
    for p in data:
        plt.scatter(p.x, p.y, c='g', marker='.')
    for i in range(len(stack)-1):
        plt.scatter(stack[i].x, stack[i].y, c='r', marker='.')
        plt.plot([stack[i].x,stack[i+1].x],[stack[i].y,stack[i+1].y], color='r')
    plt.scatter(stack[-1].x, stack[-1].y, c='r', marker='.')
    plt.plot([stack[-1].x,stack[0].x],[stack[-1].y,stack[0].y], color='r')
    plt.show()
    #'''

算法实现效果：