Delaunay

scipy.spatial.Delaunay

class scipy.spatial.Delaunay(points, furthest_site=False, incremental=False, qhull_options=None)

%matplotlib inline


#! coding:utf-8
import numpy as np
from scipy.spatial import Delaunay
import matplotlib.pyplot as plt

ポイントの作成

# メッシュ状にポイントを設置
nx, ny = (5, 5)
x = np.linspace(0, 1, nx)
y = np.linspace(0, 1, ny)
xv, yv = np.meshgrid(x, y)

# ポイントをndarrayに変換(N x 2)
points = []
for i in range(ny):
    for j in range(nx):
        points.append([xv[i, j], yv[i, j]])
points = np.asarray(points)

print points

[[ 0.    0.  ]
 [ 0.25  0.  ]
 [ 0.5   0.  ]
 [ 0.75  0.  ]
 [ 1.    0.  ]
 [ 0.    0.25]
 [ 0.25  0.25]
 [ 0.5   0.25]
 [ 0.75  0.25]
 [ 1.    0.25]
 [ 0.    0.5 ]
 [ 0.25  0.5 ]
 [ 0.5   0.5 ]
 [ 0.75  0.5 ]
 [ 1.    0.5 ]
 [ 0.    0.75]
 [ 0.25  0.75]
 [ 0.5   0.75]
 [ 0.75  0.75]
 [ 1.    0.75]
 [ 0.    1.  ]
 [ 0.25  1.  ]
 [ 0.5   1.  ]
 [ 0.75  1.  ]
 [ 1.    1.  ]]

ドロネー関数を実行し格子分割

tri = Delaunay(points)

print(tri)

&lt;scipy.spatial.qhull.Delaunay object at 0x000000000A84B588&gt;

Attributes

print tri.transform

[[[ 4.    0.  ]
  [ 0.    4.  ]
  [ 0.    0.  ]]

 [[-4.   -0.  ]
  [-0.   -4.  ]
  [ 0.25  0.25]]

....

 [[-0.    4.  ]
  [-4.    0.  ]
  [ 0.75  0.75]]

 [[ 0.   -4.  ]
  [ 4.   -0.  ]
  [ 0.5   1.  ]]]

# 各節点が含まれる要素の番号
print tri.vertex_to_simplex

[ 0  0  5  2  2  0  1  4  3  2 13 10 11 18 21  8  9 24 17 16  8  8 29 16 16]

# 境界要素の節点番号
print tri.convex_hull

[[ 5  0]
 [ 1  0]
 [ 3  4]
 [ 9  4]
 [ 1  2]
 [ 3  2]
 [21 20]
 [15 20]
 [ 5 10]
 [15 10]
 [19 24]
 [23 24]
 [ 9 14]
 [19 14]
 [21 22]
 [23 22]]

# 分からない
print tri.vertex_neighbor_vertices

(array([  0,   2,   7,  10,  15,  17,  22,  26,  34,  38,  43,  46,  54,
        58,  66,  69,  74,  78,  86,  90,  95,  97, 102, 105, 110, 112]), 
 array([ 1,  5,  5,  0,  6,  7,  2,  7,  1,  3,  9,  4,  8,  7,  2,  9,  3,
        1,  0,  6, 11, 10,  5,  1,  7, 11,  1,  6,  2,  3,  8, 11, 12, 13,
        3,  9,  7, 13,  3,  4,  8, 13, 14,  5, 11, 15,  7,  6, 12,  5, 10,
       15, 16, 17, 11,  7, 13, 17,  7, 12,  8,  9, 14, 19, 18, 17, 13,  9,
       19, 21, 20, 16, 11, 10, 21, 15, 11, 17, 11, 16, 12, 13, 18, 21, 22,
       23, 19, 23, 13, 17, 23, 24, 18, 13, 14, 15, 21, 15, 20, 16, 17, 22,
       21, 17, 23, 19, 24, 18, 17, 22, 23, 19]))

# 節点の座標を表示
print tri.points

[[ 0.    0.  ]
 [ 0.25  0.  ]
 [ 0.5   0.  ]
 [ 0.75  0.  ]
 [ 1.    0.  ]
 [ 0.    0.25]
 [ 0.25  0.25]
 [ 0.5   0.25]
 [ 0.75  0.25]
 [ 1.    0.25]
 [ 0.    0.5 ]
 [ 0.25  0.5 ]
 [ 0.5   0.5 ]
 [ 0.75  0.5 ]
 [ 1.    0.5 ]
 [ 0.    0.75]
 [ 0.25  0.75]
 [ 0.5   0.75]
 [ 0.75  0.75]
 [ 1.    0.75]
 [ 0.    1.  ]
 [ 0.25  1.  ]
 [ 0.5   1.  ]
 [ 0.75  1.  ]
 [ 1.    1.  ]]

# 要素の節点番号を表示
print tri.simplices

[[ 1  5  0]
 [ 5  1  6]
 [ 9  3  4]
 [ 3  9  8]
 [ 1  7  6]
 [ 7  1  2]
 [ 7  3  8]
 [ 3  7  2]
 [15 21 20]
 [21 15 16]
 [ 7 11  6]
 [11  7 12]
 [11  5  6]
 [ 5 11 10]
 [15 11 16]
 [11 15 10]
 [23 19 24]
 [19 23 18]
 [ 7 13 12]
 [13  7  8]
 [ 9 13  8]
 [13  9 14]
 [13 19 18]
 [19 13 14]
 [11 17 16]
 [17 11 12]
 [17 13 18]
 [13 17 12]
 [17 21 16]
 [21 17 22]
 [23 17 18]
 [17 23 22]]

# 隣接要素番号を表示(-1は境界)
print tri.neighbors

[[-1 -1  1]
 [ 4 12  0]
 [-1 -1  3]
 [20  6  2]
 [10  1  5]
 [-1  7  4]
 [ 3 19  7]
 [ 5 -1  6]
 [-1 -1  9]
 [14 28  8]
 [12  4 11]
 [18 25 10]
 [ 1 10 13]
 [15 -1 12]
 [24  9 15]
 [-1 13 14]
 [-1 -1 17]
 [30 22 16]
 [27 11 19]
 [ 6 20 18]
 [19  3 21]
 [-1 23 20]
 [17 26 23]
 [21 -1 22]
 [28 14 25]
 [11 27 24]
 [22 30 27]
 [25 18 26]
 [ 9 24 29]
 [31 -1 28]
 [26 17 31]
 [-1 29 30]]

print tri.equations

[[ 0.11219678  0.11219678 -0.98733164 -0.        ]
 [ 0.11219678  0.11219678 -0.98733164 -0.        ]
 [ 0.62007634  0.08858233 -0.77952454 -0.265747  ]
 [ 0.62007634  0.08858233 -0.77952454 -0.265747  ]

  ...
 [ 0.44291167  0.44291167 -0.77952454 -0.265747  ]
 [ 0.44291167  0.44291167 -0.77952454 -0.265747  ]
 [ 0.25777915  0.60148468 -0.75615216 -0.30074234]
 [ 0.25777915  0.60148468 -0.75615216 -0.30074234]
 [ 0.4063027   0.56882377 -0.71509274 -0.36567243]
 [ 0.4063027   0.56882377 -0.71509274 -0.36567243]]

print tri.paraboloid_scale

0.454545454545

print tri.paraboloid_shift

-0.0

print tri.coplanar

[]

# Lineをプロット
plt.triplot(points[:, 0], points[:, 1], tri.simplices.copy())
plt.plot(points[:, 0], points[:, 1], 'o')
plt.xlim([-0.2, 1.2])
plt.ylim([-0.2, 1.2])
# 節点番号の表示
for i, p in enumerate(tri.points):
    plt.text(p[0], p[1], i, ha='right')

# 要素番号の表示
for j, s in enumerate(tri.vertices):
    p = tri.points[s].mean(axis=0)
    plt.text(p[0], p[1], '#%d' % j, ha='center')


plt.show()

png

# メッシュデータを保存
# --------------------
import csv

with open('mesh2.csv', 'w') as f:
    writer = csv.writer(f, lineterminator='n')
    writer.writerow(tri.vertex_to_simplex)