[yt-svn] commit/yt: 10 new changesets

commits-noreply at bitbucket.org commits-noreply at bitbucket.org
Sat Jun 3 07:23:55 PDT 2017


10 new commits in yt:

https://bitbucket.org/yt_analysis/yt/commits/fde34f266986/
Changeset:   fde34f266986
User:        Alex Lindsay
Date:        2017-05-30 22:11:09+00:00
Summary:     Fix copysign issue and add more documentation
Affected #:  1 file

diff -r c831eed47a1848eb77ea608c2b55671485793851 -r fde34f26698658d2d74a481bf93e0946f170efdb yt/utilities/lib/geometry_utils.pyx
--- a/yt/utilities/lib/geometry_utils.pyx
+++ b/yt/utilities/lib/geometry_utils.pyx
@@ -501,12 +501,12 @@
                              np.ndarray[np.float64_t, ndim=3] triangles):
     cdef np.float64_t p0[3]
     cdef np.float64_t p1[3]
-    cdef np.float64_t p3[3]
+    cdef np.float64_t p2[3]
     cdef np.float64_t E0[3]
     cdef np.float64_t E1[3]
     cdef np.float64_t tri_norm[3]
     cdef np.float64_t plane_norm[3]
-    cdef np.float64_t dp
+    cdef np.float64_t dp, offset
     cdef int i, j, k, count, ntri, nlines
     nlines = 0
     ntri = triangles.shape[0]
@@ -517,28 +517,34 @@
     for i in range(ntri):
         count = 0
 
-        # skip if triangle is close to being parallel to plane
+        # Here p0 holds the triangle's zeroth node coordinates,
+        # p1 holds the first node's coordinates, and
+        # p2 holds the second node's coordinates
         for j in range(3):
             p0[j] = triangles[i, 0, j]
             p1[j] = triangles[i, 1, j]
-            p3[j] = triangles[i, 2, j]
+            p2[j] = triangles[i, 2, j]
             plane_norm[j] = 0.0
         plane_norm[ax] = 1.0
         subtract(p1, p0, E0)
-        subtract(p3, p0, E1)
+        subtract(p2, p0, E1)
         cross(E0, E1, tri_norm)
         dp = dot(tri_norm, plane_norm)
         dp /= L2_norm(tri_norm)
+        # Skip if triangle is close to being parallel to plane.
         if (fabs(dp) > 0.995):
             continue
-        
+
         # Now for each line segment (01, 12, 20) we check to see how many cross
-        # the coordinate.
+        # the coordinate of the slice.
+        # Here, the components of p2 are either +1 or -1 depending on whether the
+        # node's coordinate corresponding to the slice axis is greater than the
+        # coordinate of the slice. p2[0] -> node 0; p2[1] -> node 1; p2[2] -> node2
         for j in range(3):
-            p3[j] = copysign(1.0, triangles[i, j, ax] - coord)
-        if p3[0] * p3[1] < 0: count += 1
-        if p3[1] * p3[2] < 0: count += 1
-        if p3[2] * p3[0] < 0: count += 1
+            p2[j] = copysign(1.0, triangles[i, j, ax] - coord + 0)
+        if p2[0] * p2[1] < 0: count += 1
+        if p2[1] * p2[2] < 0: count += 1
+        if p2[2] * p2[0] < 0: count += 1
         if count == 2:
             nlines += 1
         elif count == 3:
@@ -548,17 +554,22 @@
         points = <PointSet *> malloc(sizeof(PointSet))
         points.count = 0
         points.next = NULL
-        if p3[0] * p3[1] < 0:
+
+        # Here p0 and p1 again hold node coordinates
+        if p2[0] * p2[1] < 0:
+            # intersection of 01 triangle segment with plane
             for j in range(3):
                 p0[j] = triangles[i, 0, j]
                 p1[j] = triangles[i, 1, j]
             get_intersection(p0, p1, ax, coord, points)
-        if p3[1] * p3[2] < 0:
+        if p2[1] * p2[2] < 0:
+            # intersection of 12 triangle segment with plane
             for j in range(3):
                 p0[j] = triangles[i, 1, j]
                 p1[j] = triangles[i, 2, j]
             get_intersection(p0, p1, ax, coord, points)
-        if p3[2] * p3[0] < 0:
+        if p2[2] * p2[0] < 0:
+            # intersection of 20 triangle segment with plane
             for j in range(3):
                 p0[j] = triangles[i, 2, j]
                 p1[j] = triangles[i, 0, j]
@@ -568,6 +579,7 @@
         if first == NULL:
             first = points
         last = points
+
     points = first
     cdef np.ndarray[np.float64_t, ndim=3] line_segments
     line_segments = np.empty((nlines, 2, 3), dtype="float64")


https://bitbucket.org/yt_analysis/yt/commits/e5db0be88960/
Changeset:   e5db0be88960
User:        Alex Lindsay
Date:        2017-05-31 22:24:28+00:00
Summary:     Get rid of what I think is unnecessary 3 count.
Affected #:  1 file

diff -r fde34f26698658d2d74a481bf93e0946f170efdb -r e5db0be889609a4c7b54bdea1e9144c85a9c3b36 yt/utilities/lib/geometry_utils.pyx
--- a/yt/utilities/lib/geometry_utils.pyx
+++ b/yt/utilities/lib/geometry_utils.pyx
@@ -478,7 +478,7 @@
 cdef struct PointSet:
     int count
     # First index is point index, second is xyz
-    np.float64_t points[3][3]
+    np.float64_t points[2][3]
     PointSet *next
 
 cdef inline void get_intersection(np.float64_t p0[3], np.float64_t p1[3],
@@ -506,7 +506,7 @@
     cdef np.float64_t E1[3]
     cdef np.float64_t tri_norm[3]
     cdef np.float64_t plane_norm[3]
-    cdef np.float64_t dp, offset
+    cdef np.float64_t dp
     cdef int i, j, k, count, ntri, nlines
     nlines = 0
     ntri = triangles.shape[0]
@@ -548,7 +548,8 @@
         if count == 2:
             nlines += 1
         elif count == 3:
-            nlines += 2
+            raise RuntimeError("I'm impressed you got a plane to intersect all three "
+                               "legs of a triangle")
         else:
             continue
         points = <PointSet *> malloc(sizeof(PointSet))
@@ -589,11 +590,6 @@
             line_segments[k, 0, j] = points.points[0][j]
             line_segments[k, 1, j] = points.points[1][j]
         k += 1
-        if points.count == 3:
-            for j in range(3):
-                line_segments[k, 0, j] = points.points[1][j]
-                line_segments[k, 1, j] = points.points[2][j]
-            k += 1
         last = points
         points = points.next
         free(last)


https://bitbucket.org/yt_analysis/yt/commits/718fa47e0ef5/
Changeset:   718fa47e0ef5
User:        Alex Lindsay
Date:        2017-06-01 02:49:02+00:00
Summary:     Add test
Affected #:  2 files

diff -r e5db0be889609a4c7b54bdea1e9144c85a9c3b36 -r 718fa47e0ef535ac42d5a274e39444cfcca4c38c yt/testing.py
--- a/yt/testing.py
+++ b/yt/testing.py
@@ -303,7 +303,7 @@
     from yt.frontends.stream.api import load_unstructured_mesh
     from yt.frontends.stream.sample_data.tetrahedral_mesh import \
         _connectivity, _coordinates
-    
+
     prng = RandomState(0x4d3d3d3)
 
     # the distance from the origin
@@ -343,6 +343,30 @@
                                 elem_data=elem_data)
     return ds
 
+def small_fake_hexahedral_ds():
+    from yt.frontends.stream.api import load_unstructured_mesh
+
+    _coordinates = np.array([[-1., -1., -1.],
+                              [ 0., -1., -1.],
+                              [ -0.,  0., -1.],
+                              [-1.,  -0., -1.],
+                              [-1., -1.,  0.],
+                              [ -0., -1.,  0.],
+                              [ -0.,  0.,  -0.],
+                              [-1.,  0.,  -0.]])
+    _connectivity = np.array([[1, 2, 3, 4, 5, 6, 7, 8]])
+
+    # the distance from the origin
+    node_data = {}
+    dist = np.sum(_coordinates**2, 1)
+    node_data[('connect1', 'test')] = dist[_connectivity-1]
+
+    ds = load_unstructured_mesh(_connectivity-1,
+                                _coordinates,
+                                node_data=node_data)
+    return ds
+
+
 
 def fake_vr_orientation_test_ds(N = 96, scale=1):
     """

diff -r e5db0be889609a4c7b54bdea1e9144c85a9c3b36 -r 718fa47e0ef535ac42d5a274e39444cfcca4c38c yt/visualization/tests/test_mesh_slices.py
--- a/yt/visualization/tests/test_mesh_slices.py
+++ b/yt/visualization/tests/test_mesh_slices.py
@@ -17,11 +17,13 @@
 import numpy as np
 import yt
 from yt.testing import fake_tetrahedral_ds
-from yt.testing import fake_hexahedral_ds
+from yt.testing import fake_hexahedral_ds, small_fake_hexahedral_ds
 from yt.utilities.answer_testing.framework import \
     requires_ds, \
     data_dir_load, \
     GenericImageTest
+from yt.utilities.lib.geometry_utils import triangle_plane_intersect
+from yt.utilities.lib.mesh_triangulation import triangulate_indices
 
 
 def setup():
@@ -70,7 +72,7 @@
     # This test is temporarily disabled because it is seg faulty,
     # see #1394
     return
-    
+
     # Perform I/O in safe place instead of yt main dir
     tmpdir = tempfile.mkdtemp()
     curdir = os.getcwd()
@@ -99,3 +101,25 @@
     os.chdir(curdir)
     # clean up
     shutil.rmtree(tmpdir)
+
+def test_perfect_element_intersection():
+    '''
+    This test tests mesh line annotation where a z=0 slice
+    perfectly intersects the top of a hexahedral element with node
+    z-coordinates containing both -0 and +0. Before
+    https://github.com/yt-project/yt/pull/1437 this test falsely
+    yielded three annotation lines, whereas the correct result is four
+    corresponding to the four edges of the top hex face.
+    '''
+    ds = small_fake_hexahedral_ds()
+    indices = ds.index.meshes[0].connectivity_indices
+    coords = ds.index.meshes[0].connectivity_coords
+    tri_indices = triangulate_indices(indices)
+    tri_coords = coords[tri_indices]
+    lines = triangle_plane_intersect(2, 0, tri_coords)
+    non_zero_lines = 0
+    for i in range(lines.shape[0]):
+        norm = np.linalg.norm(lines[i][0] - lines[i][1])
+        if norm > 1e-8:
+            non_zero_lines += 1
+    np.testing.assert_equal(non_zero_lines, 4)


https://bitbucket.org/yt_analysis/yt/commits/2b8cc3b06bec/
Changeset:   2b8cc3b06bec
User:        Alex Lindsay
Date:        2017-06-01 14:42:10+00:00
Summary:     Change int64 to intp in windows erroring cython function.
Affected #:  1 file

diff -r 718fa47e0ef535ac42d5a274e39444cfcca4c38c -r 2b8cc3b06becd4422f1a088bc27ffbe904851329 yt/utilities/lib/mesh_triangulation.pyx
--- a/yt/utilities/lib/mesh_triangulation.pyx
+++ b/yt/utilities/lib/mesh_triangulation.pyx
@@ -299,7 +299,7 @@
 
 @cython.boundscheck(False)
 @cython.wraparound(False)
-def triangulate_indices(np.int64_t[:, ::1] indices):
+def triangulate_indices(np.intp_t[:, ::1] indices):
     '''
 
     This is like triangulate_mesh, except it only considers the
@@ -309,9 +309,9 @@
     '''
 
     cdef MeshInfoHolder m = MeshInfoHolder(indices)
-    cdef np.int64_t[:, ::1] tri_indices = np.empty((m.num_tri, 3), dtype=np.int64)
+    cdef np.intp_t[:, ::1] tri_indices = np.empty((m.num_tri, 3), dtype=np.intp)
 
-    cdef np.int64_t i, j, k
+    cdef np.intp_t i, j, k
     for i in range(m.num_elem):
         for j in range(m.TPE):
             for k in range(3):


https://bitbucket.org/yt_analysis/yt/commits/c17af541e6ee/
Changeset:   c17af541e6ee
User:        Alex Lindsay
Date:        2017-06-01 16:34:36+00:00
Summary:     Try this again. Numpy, cython, and windows...rage
Affected #:  1 file

diff -r 2b8cc3b06becd4422f1a088bc27ffbe904851329 -r c17af541e6ee4bac5b13155cfd9145192e4f99e1 yt/utilities/lib/mesh_triangulation.pyx
--- a/yt/utilities/lib/mesh_triangulation.pyx
+++ b/yt/utilities/lib/mesh_triangulation.pyx
@@ -299,7 +299,7 @@
 
 @cython.boundscheck(False)
 @cython.wraparound(False)
-def triangulate_indices(np.intp_t[:, ::1] indices):
+def triangulate_indices(np.int_t[:, ::1] indices):
     '''
 
     This is like triangulate_mesh, except it only considers the
@@ -309,9 +309,9 @@
     '''
 
     cdef MeshInfoHolder m = MeshInfoHolder(indices)
-    cdef np.intp_t[:, ::1] tri_indices = np.empty((m.num_tri, 3), dtype=np.intp)
+    cdef np.int_t[:, ::1] tri_indices = np.empty((m.num_tri, 3), dtype=np.int_)
 
-    cdef np.intp_t i, j, k
+    cdef np.int_t i, j, k
     for i in range(m.num_elem):
         for j in range(m.TPE):
             for k in range(3):


https://bitbucket.org/yt_analysis/yt/commits/ab55d6600b3d/
Changeset:   ab55d6600b3d
User:        ngoldbaum
Date:        2017-06-01 17:56:12+00:00
Summary:     fix windows errors
Affected #:  2 files

diff -r c17af541e6ee4bac5b13155cfd9145192e4f99e1 -r ab55d6600b3dc508a955fbce4f0bcd04c0c4ed09 yt/utilities/lib/mesh_triangulation.pyx
--- a/yt/utilities/lib/mesh_triangulation.pyx
+++ b/yt/utilities/lib/mesh_triangulation.pyx
@@ -182,7 +182,7 @@
     cdef np.int64_t TPE  # num tris per element
     cdef int[MAX_NUM_TRI][3] tri_array
 
-    def __cinit__(self, np.int64_t[:, ::1] indices):
+    def __cinit__(self, np.int_t[:, ::1] indices):
         '''
 
         This class is used to store metadata about the type of mesh being used.

diff -r c17af541e6ee4bac5b13155cfd9145192e4f99e1 -r ab55d6600b3dc508a955fbce4f0bcd04c0c4ed09 yt/visualization/plot_modifications.py
--- a/yt/visualization/plot_modifications.py
+++ b/yt/visualization/plot_modifications.py
@@ -1686,7 +1686,7 @@
             elif num_dims == 2 and num_verts == 4:
                 coords, indices = self.promote_2d_to_3d(coords, indices, plot)
 
-            tri_indices = triangulate_indices(indices)
+            tri_indices = triangulate_indices(indices.astype(np.int_))
             points = coords[tri_indices]
         
             tfc = TriangleFacetsCallback(points, plot_args=self.plot_args)


https://bitbucket.org/yt_analysis/yt/commits/3d0bdeb10534/
Changeset:   3d0bdeb10534
User:        Alex Lindsay
Date:        2017-06-01 18:29:16+00:00
Summary:     Re-enable test_mesh_slices (fixes #1394)
Affected #:  1 file

diff -r ab55d6600b3dc508a955fbce4f0bcd04c0c4ed09 -r 3d0bdeb1053459eeda468675a6e7c8869f4f5163 yt/visualization/tests/test_mesh_slices.py
--- a/yt/visualization/tests/test_mesh_slices.py
+++ b/yt/visualization/tests/test_mesh_slices.py
@@ -69,10 +69,6 @@
         yield compare(ds, field, "answers_multi_region_%s_%s" % (field[0], field[1]))
 
 def test_mesh_slices():
-    # This test is temporarily disabled because it is seg faulty,
-    # see #1394
-    return
-
     # Perform I/O in safe place instead of yt main dir
     tmpdir = tempfile.mkdtemp()
     curdir = os.getcwd()


https://bitbucket.org/yt_analysis/yt/commits/0c468e904325/
Changeset:   0c468e904325
User:        Alex Lindsay
Date:        2017-06-02 19:51:32+00:00
Summary:     Replace old hexahedral ds with new one.
Affected #:  1 file

This diff is so big that we needed to truncate the remainder.

https://bitbucket.org/yt_analysis/yt/commits/e925d996de11/
Changeset:   e925d996de11
User:        Alex Lindsay
Date:        2017-06-03 01:43:08+00:00
Summary:     Address nathan review
Affected #:  2 files

diff -r 0c468e904325e891482f5003c43f5808adc32375 -r e925d996de1183abc903afc98ef6fe3a2d68e376 yt/utilities/lib/geometry_utils.pyx
--- a/yt/utilities/lib/geometry_utils.pyx
+++ b/yt/utilities/lib/geometry_utils.pyx
@@ -541,6 +541,8 @@
         # node's coordinate corresponding to the slice axis is greater than the
         # coordinate of the slice. p2[0] -> node 0; p2[1] -> node 1; p2[2] -> node2
         for j in range(3):
+            # Add 0 so that any -0s become +0s. Necessary for consistent determination
+            # of plane intersection
             p2[j] = copysign(1.0, triangles[i, j, ax] - coord + 0)
         if p2[0] * p2[1] < 0: count += 1
         if p2[1] * p2[2] < 0: count += 1
@@ -548,8 +550,9 @@
         if count == 2:
             nlines += 1
         elif count == 3:
-            raise RuntimeError("I'm impressed you got a plane to intersect all three "
-                               "legs of a triangle")
+            raise RuntimeError("It should be geometrically impossible for a plane to"
+                               "to intersect all three legs of a triangle. Please contact"
+                               "yt developers with your mesh")
         else:
             continue
         points = <PointSet *> malloc(sizeof(PointSet))

diff -r 0c468e904325e891482f5003c43f5808adc32375 -r e925d996de1183abc903afc98ef6fe3a2d68e376 yt/visualization/tests/test_mesh_slices.py
--- a/yt/visualization/tests/test_mesh_slices.py
+++ b/yt/visualization/tests/test_mesh_slices.py
@@ -99,14 +99,13 @@
     shutil.rmtree(tmpdir)
 
 def test_perfect_element_intersection():
-    '''
-    This test tests mesh line annotation where a z=0 slice
-    perfectly intersects the top of a hexahedral element with node
-    z-coordinates containing both -0 and +0. Before
-    https://github.com/yt-project/yt/pull/1437 this test falsely
-    yielded three annotation lines, whereas the correct result is four
-    corresponding to the four edges of the top hex face.
-    '''
+    # This test tests mesh line annotation where a z=0 slice
+    # perfectly intersects the top of a hexahedral element with node
+    # z-coordinates containing both -0 and +0. Before
+    # https://github.com/yt-project/yt/pull/1437 this test falsely
+    # yielded three annotation lines, whereas the correct result is four
+    # corresponding to the four edges of the top hex face.
+
     ds = small_fake_hexahedral_ds()
     indices = ds.index.meshes[0].connectivity_indices
     coords = ds.index.meshes[0].connectivity_coords


https://bitbucket.org/yt_analysis/yt/commits/026d5035a546/
Changeset:   026d5035a546
User:        ngoldbaum
Date:        2017-06-03 14:23:37+00:00
Summary:     Merge pull request #1437 from lindsayad/mesh_annotations_1273

Fix mesh line annotation for triangulations with a side in the plane
Affected #:  6 files

This diff is so big that we needed to truncate the remainder.

Repository URL: https://bitbucket.org/yt_analysis/yt/

--

This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.


More information about the yt-svn mailing list