GJK算法
1. 预备知识
  1. 闵可夫斯基差
  2. 支持向量与支持函数
2. GJK算法实现
  1. 参考资料

GJK算法

Gilbert–Johnson–Keerthi 算法(GJK算法)是碰撞检测Narrow-Phase常用的一种算法。可以计算二维或者三维的碰撞。以下介绍二维情况，且可以被比较简单地推广到三维。

预备知识

闵可夫斯基差

两个集合的闵可夫斯基差(Minkowski Difference) $A\ominus B = \{x-y|x\in A\space \land\space y \in B\}$

则$A\cap B \neq \emptyset \Leftrightarrow O \in A\ominus B$

可以证明，如果$A,B$都是凸多边形，那么$A\ominus B$也是凸多边形，且$A\ominus B$的顶点的来源是$A,B$的顶点。

如此，我们可以考虑如下算法：

bool bruteForce(Polygon A,Polygon B)
{
    vector<Vector3> points;
	for(auto vertexa : A.vertex){
        for(auto vertexb : B.vertex){
        	points.push_back(vertexa - vertexB);
        }
    }
    Polygon minkowskiDiff = FindConvexHull(points);
    return minkowskiDiff.hasOrigin();
}

上述算法的复杂度是$O(n^2logn)$

支持向量与支持函数

对于集合$A$定义支持函数(support function)如下：

$support(A,\vec d) = \underset{a \in A}{argmax}\space \vec a \cdot\vec d $

即集合$A$中所有点中与$\vec d$点乘最大的点，也即在$\vec d$方向上最远的点，称$\vec d$为支持向量。

可以证明，若$A$是一个凸多边形，则当$\vec d$旋转一周后，$support(A,\vec d)$正好遍历每个顶点一次。

于是不难发现有$support(A\ominus B,\vec d) = support(A,\vec d)-support(B,-\vec d)$

GJK算法实现

求出整个$A\ominus B$需要$O(n^2logn)$的时间，这是我们不能接受的。于是我们考虑维护一个单纯形(即确定一块区域的最简图形，二维为三角形，三维为四面体)，如果这个三角形包含远点，则$O\in A\ominus B$。如果它不包含，我们通过求解逐步调整这个三角形使他逼近远点，直到不能继续为止(详见参考资料/ppt)。

直接上代码：

public class GJK : MonoBehaviour
{
    private static Vector3 findFurthest(List<Vector3> vertexs,Vector3 dir){
        //求解support(A,dir)
        Vector3 res = vertexs[0];
        float maxVal = Vector3.Dot(res,dir);
        foreach(var vertex in vertexs){//遍历所有顶点
            if(Vector3.Dot(vertex,dir) > maxVal){
                maxVal = Vector3.Dot(vertex,dir);
                res = vertex;
            }
        }
        return res;
    }
    private static Vector3 support(List<Vector3> polygonA,List<Vector3> polygonB,Vector3 dir){
    //求解support(A-B,dir)
        return findFurthest(polygonA,dir) - findFurthest(polygonB,-dir);
    }
    public static bool nextSimplex(List<Vector3> simplex,ref Vector3 d){
    //求下一个单纯形
        if(simplex.Count == 2){
            return lineCase(simplex,ref d); //两个点的情况
        }
        return triangleCase(simplex,ref d);//三角形的情况
    }
    public static bool lineCase(List<Vector3> simplex,ref Vector3 d){
        Vector3 B = simplex[0],A = simplex[1];
        Vector3 AB = B - A,AO = -A;
        //与AB垂直，指向原点的向量
        Vector3 ABp = Vector3.Cross(Vector3.Cross(AB,AO),AB);
        d = ABp;
        return false;
    }
    public static bool triangleCase(List<Vector3> simplex,ref Vector3 d){
    	//按C,B,A顺序加入单纯形的三个点
        Vector3 C = simplex[0],B = simplex[1],A = simplex[2];
        Vector3 AB = B - A,AC = C - A,AO = -A;
       	//与AB垂直的向量
        Vector3 ABp = Vector3.Cross(Vector3.Cross(AC,AB),AB);
        //与AC垂直的向量
        Vector3 ACp = Vector3.Cross(Vector3.Cross(AB,AC),AC);
        //在R_AB区域
        if(Vector3.Dot(ABp,AO) > 0){
            simplex.Remove(C);
            d = ABp;
            return false;
        }
        //在R_AC区域
        if(Vector3.Dot(ACp,AO) > 0){
            simplex.Remove(B);
            d = ACp;
            return false;
        }
        //在三角形内部
        return true;
    }
    public static bool gjk(List<Vector3> polygonA,List<Vector3> polygonB){
        Vector3 d = support(polygonA,polygonB,Vector3.right);//初始支持向量是任意方向
        List<Vector3> simplex = new List<Vector3>();
        simplex.Add(d);
        d = -d;//调整支持向量指向远点
        while(true){
            Vector3 A = support(polygonA,polygonB,d);
            if(Vector3.Dot(A,d) < 0){
                return false; //不可能包含原点
            }
            simplex.Add(A);
            if(nextSimplex(simplex,ref d)){//求下一个单纯形
                return true;
            }
        }
    }
}