Speedup mesh_xsection more #6

algrs · 2017-06-26T14:10:48Z

No description provided.

Use a more specific data structure.

maxdumas · 2017-06-26T14:53:59Z

blmath/geometry/primitives/plane.py


        class Graph(object):
-            # A little utility class to build a symmetric graph
+            # A little utility class to build a symmetric graph and calcualate Euler Paths


typo: calcualate -> calculate 😄

What I get for turning off the speel cheeker...

maxdumas · 2017-06-26T15:11:25Z

blmath/geometry/primitives/plane.py

+            e0 = intersection_map.index(f[0], f[1])
+            e1 = intersection_map.index(f[0], f[2])
+            e2 = intersection_map.index(f[1], f[2])
+            if e0 is None:


This is operating on a mathematical certainty that exactly two of the edges of a face that intersects the plane will themselves intersect the plane. Can you make that explicit with a comment? The reasoning for this logic otherwise isn't immediately obvious from just the code.

maxdumas · 2017-06-26T15:14:42Z

blmath/geometry/primitives/plane.py


        # 5: Find the paths for each component
        components = []
        components_closed = []
        while len(G) > 0:
-            # This works because euler_path mutates the graph as it goes
-            path = euler_path(G.d)
+            # This works because pop_euler_path mutates the graph as it goes


This comment doesn't seem necessary anymore now that we have a more appropriate abstraction. 👍

And it's commented up in the function as well. Will remove.

maxdumas

This looks great. The abstractions you added really make this function a lot more readable and less intimidating!

maxdumas · 2017-06-26T15:39:24Z

blmath/geometry/primitives/plane.py

-                    path.append(stack.pop())
-            return path
+        for f in fs:
+            # Since we're dealing with a triangle that intersects the plane, exactly of the edges


exactly two of the edges

algrs

@jbwhite Some thoughts on profiling for you:

algrs · 2017-06-26T15:47:09Z

blmath/geometry/primitives/plane.py

@@ -16,8 +16,8 @@ def __init__(self, point_on_plane, unit_normal):

        unit_normal = vx.normalize(unit_normal)

-        self._r0 = point_on_plane
-        self._n = unit_normal
+        self._r0 = np.asarray(point_on_plane)


For those following along at home, moving the np.array from the getter to the constructor saved about half a million memory allocations. It's a 3 element array, but still those memory allocations added up.

algrs · 2017-06-26T15:49:00Z

blmath/geometry/primitives/plane.py

        if pt is not None:
-            if np.any(pt < np.min((a, b), axis=0)) or np.any(pt > np.max((a, b), axis=0)):
+            if any(np.logical_and(pt > a, pt > b)) or any(np.logical_and(pt < a, pt < b)):


Here's an interesting example of a hotspot that profiling detected. We call this function half a million times, so tiny optimizations help a lot. Since numpy is optimized for large arrays, something like np.min is not as efficient as you'd expect when called many many times on tiny arrays.

@algrs I get that we want things to be readable, but it's definitely more efficient to vectorize the calls to these. It appears to save a few seconds to use:

def line_xsections(self, pts, rays): denoms = np.dot(rays, self.normal) denom_is_zero = denoms == 0 denoms[denom_is_zero] = np.nan p = np.dot(self.reference_point - pts, self.normal) / denoms return np.vstack([p, p, p]).T * rays + pts, ~denom_is_zero def line_segment_xsections(self, a, b): pts, pt_is_valid = self.line_xsections(a, b-a) valid_pts = pts[pt_is_valid] pt_is_out_of_bounds = np.logical_or(np.any(np.logical_and(valid_pts > a, valid_pts > b), axis=1), np.any(np.logical_and(valid_pts < a, valid_pts < b), axis=1)) pt_is_valid[pt_is_valid] = ~pt_is_out_of_bounds pts[~pt_is_valid] = np.nan return pts, pt_is_valid

and then down below instead of the for loop over es you have

pts, pt_is_valid = self.line_segment_xsections(m.v[es[:,0]], m.v[es[:,1]]) valid_pts = pts[pt_is_valid] valid_es = es[pt_is_valid] for val, e in zip(valid_pts, valid_es): if not intersection_map.contains(e[0], e[1]): intersection_map.add(e[0], e[1], val)

And then on top of that I'd try to allow intersection_map to support batch operations

Nothing wrong with vectorizing: #7

algrs · 2017-06-26T15:50:48Z

blmath/geometry/primitives/plane.py

@@ -245,6 +242,11 @@ def mesh_xsection(self, m, neighborhood=None):
            return Polyline(None)
        return Polyline(np.vstack([x.v for x in components]), closed=True)

+    def mesh_intersecting_faces(self, m):


Sometimes when profiling in python it's useful to break a large function up into many small ones, so that the profiler can figure out which part of a large function has the slow code. Mostly I removed the subfunctions afterwards, because they didn't add clarity, but some of them turned out to be improvements and got left in.

algrs · 2017-06-26T15:53:33Z

blmath/geometry/primitives/plane.py

+            def __init__(self):
+                self.d = {} # store indicies into self.values here, to make it easier to get inds or values
+                self.values = []
+            def _order(self, u, v):


Here's another bit hotspot optimization. This method gets called a bazillion times, and it's faster to do a quick if statement than the more pythonic sorted((u, v)) or a use of min & max. This optimization saved a second and a half.

@algrs My hunch is that EdgeMap would be a good bit faster if it were built on top of a sparse matrix class. This is particularly true if edges are added in bulk using my proposed line_segment_xsections above.

Yeah, I tried it. Using a sparse matrix for representing a graph is a huge gain when you're doing complex operations on a fixed graph. But here, a dictionary of sets implementation is faster, as you're closer to O(1) for insertion and deletion. There's room to play here, but it's maybe 3/4 of a second of slack; the lower hanging fruit is elsewhere.

algrs · 2017-06-26T15:55:24Z

blmath/geometry/primitives/plane.py

+
+        intersection_map = EdgeMap()
+        for e in es:
+            if not intersection_map.contains(e[0], e[1]):


In the initial implementation, we called line_segment_xsection 700000 times. By not duplicating calls (remember, most intersecting faces will share an edge or two with another intersecting face), we saved 250000 calls and almost 2 seconds. So sometimes your hotspot means "improve the function" and sometimes it means "call the function less often".

algrs · 2017-06-26T15:56:51Z

blmath/geometry/primitives/plane.py

+            # Since we're dealing with a triangle that intersects the plane, exactly two of the edges
+            # will intersect (note that the only other sorts of "intersections" are one edge in
+            # plane or all three edges in plane, which won't be picked up by mesh_intersecting_faces).
+            e0 = intersection_map.index(f[0], f[1])


Wicked non-pythonic code here, but relying on the assumption that there are always three edges and two of them will always intersect lets us save several seconds over a more traditional list comprehension.

jbwhite

Appreciate the comments and explanations!

Tests actually were not running in Python 3. Fixing that revealed a dependency on baiji, which does not run in Python 3. This probably is not an easy dependency to remove, so drop Python 3 support for now.

algrs added 8 commits June 23, 2017 14:44

Speed up graph creation.

d79d57c

speed up plane line segment intersection

af6acdb

speed up plane line segment intersection

cdece83

More speedup of mesh xsection.

86c4b56

Use a more specific data structure.

more speedups

0cf995f

more speedups

e611074

revert not worth it

e201d60

fix lint

fef4a81

algrs requested review from sashaverma and jbwhite June 26, 2017 14:10

maxdumas reviewed Jun 26, 2017

View reviewed changes

fix typo

c35b2e0

maxdumas reviewed Jun 26, 2017

View reviewed changes

updates from PR

129eb9f

maxdumas approved these changes Jun 26, 2017

View reviewed changes

maxdumas reviewed Jun 26, 2017

View reviewed changes

is it a typo or a thinko?

fc644cf

algrs commented Jun 26, 2017

View reviewed changes

jbwhite approved these changes Jun 26, 2017

View reviewed changes

algrs merged commit a295a58 into master Jun 26, 2017

algrs deleted the 128_speedup_xsection branch June 26, 2017 17:07

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Speedup mesh_xsection more #6

Speedup mesh_xsection more #6

algrs commented Jun 26, 2017

maxdumas Jun 26, 2017

algrs Jun 26, 2017

maxdumas Jun 26, 2017 •

edited

Loading

algrs Jun 26, 2017

maxdumas Jun 26, 2017

algrs Jun 26, 2017

maxdumas left a comment

maxdumas Jun 26, 2017

algrs left a comment

algrs Jun 26, 2017

algrs Jun 26, 2017

eerac Jun 29, 2017 •

edited

Loading

algrs Jun 29, 2017

algrs Jun 26, 2017

algrs Jun 26, 2017

eerac Jun 29, 2017 •

edited

Loading

algrs Jun 29, 2017

algrs Jun 26, 2017

algrs Jun 26, 2017

jbwhite left a comment

Speedup mesh_xsection more #6

Speedup mesh_xsection more #6

Conversation

algrs commented Jun 26, 2017

Choose a reason for hiding this comment

Choose a reason for hiding this comment

maxdumas Jun 26, 2017 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

maxdumas left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

algrs left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

eerac Jun 29, 2017 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

eerac Jun 29, 2017 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jbwhite left a comment

Choose a reason for hiding this comment

maxdumas Jun 26, 2017 •

edited

Loading

eerac Jun 29, 2017 •

edited

Loading

eerac Jun 29, 2017 •

edited

Loading