Why does this numpy construct work and is there a more correct way of doing it

I’m using numpy to extract faces from tetrahdra defined by vertex indices. I have an initial array defining the tehrahedra toplogy

tetrahedra = np.array([[0, 1, 2, 3], [1, 2, 3, 4], [2, 3, 4, 5]])

For each tetrahedra I identify the faces using mask array

face1 = [True, True, True, False]
face2 = [True, True, False, True]
face3 = [False, True, True, True]
face4 = [True, False, True, True]

And I find the following numpy expression yields face defnitions for each tetrahedra

faces = tetrahedra[:,np.reshape(np.r_[tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3 ],tetrahedra[0][face4 ]], (-1,3))]

EDIT : Thanks to @npaulj I now see that this only appears to work because tetrahedra[0] in the indexing notation is actually [0, 1, 2, 3]. This is better expressed by changing the boolean masks to direct index masks as follows,

mask1_ = np.array([0, 1, 2])
mask2_ = np.array([0, 1, 3])
mask3_ = np.array([1, 2, 3])
mask4_ = np.array([0, 2, 3])

and then updating the expression to

faces = tetrahedra[:,np.reshape(np.r_[face1_, face2_, face3_, face4_], (-1,3))]

Now my question is, how is this actually working and is there a preferred/faster way of doing this operation? The output is shown below

THanks in advance for any help with this. I’m tempted to just live with it since it appears to work, but I can’t fugre out why it works which worries me…..

[[[0 1 2]
  [0 1 3]
  [1 2 3]
  [0 2 3]]
 [[1 2 3]
  [1 2 4]
  [2 3 4]
  [1 3 4]]
 [[2 3 4]
  [2 3 5]
  [3 4 5]
  [2 4 5]]

]

edit Cleaned up version is,

face_masks = np.array([[0, 1, 2], [0, 1, 3], [1, 2, 3],[0, 2, 3]])
faces = tetrahedra[:,face_masks]

Answer

All the intermediate steps:

In [77]: tetrahedra = np.array([[0, 1, 2, 3], [1, 2, 3, 4], [2, 3, 4, 5]])
    ...: face1 = [True, True, True, False]
    ...: face2 = [True, True, False, True]
    ...: face3 = [False, True, True, True]
    ...: face4 = [True, False, True, True]

In [78]: tetrahedra[0]
Out[78]: array([0, 1, 2, 3])

In [79]: tetrahedra[0][face1]
Out[79]: array([0, 1, 2])

r_ concatenates these 4 selections into one array:

In [80]: np.r_[tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3 ],tetrahedra[0][face4 ]]
Out[80]: array([0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3])

In [81]: np.reshape(np.r_[tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3
    ...:  ],tetrahedra[0][face4 ]], (-1,3))
Out[81]: 
array([[0, 1, 2],
       [0, 1, 3],
       [1, 2, 3],
       [0, 2, 3]])

And finally just use that to index the columns dimension of tetrahedra. Index a (3,4) with this (4,3) produces a (3,4,3).

vstack could also be used to concatenate the selections, but that’s a minor change:

In [82]: np.vstack([tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3 ],tet
    ...: rahedra[0][face4 ]])
Out[82]: 
array([[0, 1, 2],
       [0, 1, 3],
       [1, 2, 3],
       [0, 2, 3]])

edit

Or if you don’t want to count on tetrahedra[0] being [0,1,2,3], you could just seek the indices of the True elements:

In [106]: np.nonzero([face1,face2,face3,face4])
Out[106]: 
(array([0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3]),
 array([0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3]))

The 2nd element of that tuple [1] is the same indices you with r_.

In [122]: idx = np.nonzero([face1,face2,face3,face4])
In [123]: idx[1]
Out[123]: array([0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3])
In [124]: tetrahedra[np.arange(3)[:,None],idx[1]]
Out[124]: 
array([[0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3],
       [1, 2, 3, 1, 2, 4, 2, 3, 4, 1, 3, 4],
       [2, 3, 4, 2, 3, 5, 3, 4, 5, 2, 4, 5]])
In [125]: tetrahedra[np.arange(3)[:,None],idx[1]].reshape(3,4,3)
Out[125]: 
array([[[0, 1, 2],
        [0, 1, 3],
        [1, 2, 3],
        [0, 2, 3]],

       [[1, 2, 3],
        [1, 2, 4],
        [2, 3, 4],
        [1, 3, 4]],

       [[2, 3, 4],
        [2, 3, 5],
        [3, 4, 5],
        [2, 4, 5]]])

Leave a Reply

Your email address will not be published. Required fields are marked *