scipy.spatial.Delaunay 関数を使用して、Python の Matlab delaunayn 関数によって実行される N 次元 Delaunay 三角形分割を複製しようとしています。ただし、Matlab 関数は私が望んで期待する結果をもたらしますが、scipy は別の結果をもたらします。どちらも QHull ライブラリのラッパーであることを考えると、これは奇妙だと思います。Matlab は、その呼び出しで暗黙的にさまざまなパラメーターを設定していると思います。私がそれらの 2 つの間で再現しようとしている状況は、Matlab のドキュメントにあります。
セットアップは、以下のように中心にポイントを持つ立方体を持つことです。形状を視覚化するために提供した青い線ですが、この問題には何の目的も意味もありません。
この結果から予想される三角形分割は、12 個のシンプリス (Matlab の例にリストされています) になり、次のようになります。
ただし、この python に相当するものは「余分な」単純化を生成します。
x = np.array([[-1,-1,-1],[-1,-1,1],[-1,1,-1],[1,-1,-1],[1,1,1],[1,1,-1],[1,-1,1],[-1,1,1],[0,0,0]])
simp = scipy.spatial.Delaunay(x).simplices
返される変数simpは、M x N 配列である必要があります。ここで、M は見つかったシンプレックスの数 (私の場合は 12 である必要があります) であり、N はシンプレックス内のポイントの数です。この場合、各シンプレックスは N が 4 であることを意味する四面体である必要があります。
私が見つけたのは、M が実際には 18 であり、余分な 6 つのシンプリックスは四面体ではなく、立方体の 6 つの面であるということです。
何が起きてる？返されるシンプリックスを四面体のみに制限するにはどうすればよいですか? この単純なケースを使用して問題を示したので、この問題に合わせていないソリューションが必要です。
編集
Amro の回答のおかげで、これを理解することができ、Matlab と Scipy の間で簡単に一致させることができました。これには 2 つの要因がありました。まず、指摘したように、Matlab と Scipy は異なる QHull オプションを使用します。次に、QHull はボリュームがゼロのシンプリスを返します。Matlab はこれらを削除しますが、Scipy は削除しません。上記の例では、6 つの余分なシンプリスがすべて立方体のゼロ ボリュームの同一平面上にあるため、これは明らかです。これらは、次のコードを使用して、N 次元で削除できます。
N = 3 # The dimensions of our points
options = 'Qt Qbb Qc' if N <= 3 else 'Qt Qbb Qc Qx' # Set the QHull options
tri = scipy.spatial.Delaunay(points, qhull_options = options).simplices
keep = np.ones(len(tri), dtype = bool)
for i, t in enumerate(tri):
if abs(np.linalg.det(np.hstack((points[t], np.ones([1,N+1]).T)))) < 1E-15:
keep[i] = False # Point is coplanar, we don't want to keep it
tri = tri[keep]
他の条件に対処する必要があると思いますが、ポイントには既に重複が含まれていないことが保証されており、方向条件は識別できる出力に影響を与えていないようです.