Speed up the find_containing_tetrahedron function by using Data.Vector.

diff --git a/src/Cube.hs b/src/Cube.hs
index 6e31f1911cdef2ed70fe7fcdf836d0bed541c3d8..54fda14a08d078baea20c00f8898e18a2d7a924b 100644 (file)
--- a/src/Cube.hs
+++ b/src/Cube.hs
@@ -1,7 +1,16 @@
 module Cube
 where
 
-import Data.List ( (\\) )
+import Data.Maybe (fromJust)
+import qualified Data.Vector as V (
+  Vector,
+  findIndex,
+  map,
+  minimum,
+  singleton,
+  snoc,
+  unsafeIndex
+  )
 import Test.QuickCheck (Arbitrary(..), Gen, Positive(..), choose)
 
 import Cardinal
@@ -462,39 +471,82 @@ tetrahedron c 23 =
 tetrahedron _ _ = error "asked for a nonexistent tetrahedron"
 
 
+-- Only used in tests, so we don't need the added speed
+-- of Data.Vector.
 tetrahedra :: Cube -> [Tetrahedron]
-tetrahedra c =
-    [ tetrahedron c n | n <- [0..23] ]
-
--- | All completely contained in the front half of the cube.
-front_half_tetrahedra :: Cube -> [Tetrahedron]
-front_half_tetrahedra c =
-  [ tetrahedron c n | n <- [0,1,2,3,6,12,19,21] ]
-
--- | All tetrahedra completely contained in the top half of the cube.
-top_half_tetrahedra :: Cube -> [Tetrahedron]
-top_half_tetrahedra c =
-  [ tetrahedron c n | n <- [4,5,6,7,0,10,16,20] ]
-
--- | All tetrahedra completely contained in the back half of the cube.
-back_half_tetrahedra :: Cube -> [Tetrahedron]
-back_half_tetrahedra c =
-  [ tetrahedron c n | n <- [8,9,10,11,4,14,17,23] ]
-
--- | All tetrahedra completely contained in the down half of the cube.
-down_half_tetrahedra :: Cube -> [Tetrahedron]
-down_half_tetrahedra c =
-  [ tetrahedron c n | n <- [12,13,14,15,2,8,18,22] ]
-
--- | All tetrahedra completely contained in the right half of the cube.
-right_half_tetrahedra :: Cube -> [Tetrahedron]
-right_half_tetrahedra c =
-  [ tetrahedron c n | n <- [16,17,18,19,1,5,9,13] ]
-
--- | All tetrahedra completely contained in the left half of the cube.
-left_half_tetrahedra :: Cube -> [Tetrahedron]
-left_half_tetrahedra c =
-  [ tetrahedron c n | n <- [20,21,22,23,3,7,11,15] ]
+tetrahedra c = [ tetrahedron c n | n <- [0..23] ]
+
+front_left_top_tetrahedra :: Cube -> V.Vector Tetrahedron
+front_left_top_tetrahedra  c =
+  V.singleton (tetrahedron c 0) `V.snoc`
+    (tetrahedron c 3) `V.snoc`
+    (tetrahedron c 6) `V.snoc`
+    (tetrahedron c 7) `V.snoc`
+    (tetrahedron c 20) `V.snoc`
+    (tetrahedron c 21)
+
+front_left_down_tetrahedra :: Cube -> V.Vector Tetrahedron
+front_left_down_tetrahedra  c =
+  V.singleton (tetrahedron c 0) `V.snoc`
+    (tetrahedron c 2) `V.snoc`
+    (tetrahedron c 3) `V.snoc`
+    (tetrahedron c 12) `V.snoc`
+    (tetrahedron c 15) `V.snoc`
+    (tetrahedron c 21)
+
+front_right_top_tetrahedra :: Cube -> V.Vector Tetrahedron
+front_right_top_tetrahedra  c =
+  V.singleton (tetrahedron c 0) `V.snoc`
+    (tetrahedron c 1) `V.snoc`
+    (tetrahedron c 5) `V.snoc`
+    (tetrahedron c 6) `V.snoc`
+    (tetrahedron c 16) `V.snoc`
+    (tetrahedron c 19)
+
+front_right_down_tetrahedra :: Cube -> V.Vector Tetrahedron
+front_right_down_tetrahedra  c =
+  V.singleton (tetrahedron c 1) `V.snoc`
+    (tetrahedron c 2) `V.snoc`
+    (tetrahedron c 12) `V.snoc`
+    (tetrahedron c 13) `V.snoc`
+    (tetrahedron c 18) `V.snoc`
+    (tetrahedron c 19)
+
+back_left_top_tetrahedra :: Cube -> V.Vector Tetrahedron
+back_left_top_tetrahedra  c =
+  V.singleton (tetrahedron c 0) `V.snoc`
+    (tetrahedron c 3) `V.snoc`
+    (tetrahedron c 6) `V.snoc`
+    (tetrahedron c 7) `V.snoc`
+    (tetrahedron c 20) `V.snoc`
+    (tetrahedron c 21)
+
+back_left_down_tetrahedra :: Cube -> V.Vector Tetrahedron
+back_left_down_tetrahedra  c =
+  V.singleton (tetrahedron c 8) `V.snoc`
+    (tetrahedron c 11) `V.snoc`
+    (tetrahedron c 14) `V.snoc`
+    (tetrahedron c 15) `V.snoc`
+    (tetrahedron c 22) `V.snoc`
+    (tetrahedron c 23)
+
+back_right_top_tetrahedra :: Cube -> V.Vector Tetrahedron
+back_right_top_tetrahedra  c =
+  V.singleton (tetrahedron c 4) `V.snoc`
+    (tetrahedron c 5) `V.snoc`
+    (tetrahedron c 9) `V.snoc`
+    (tetrahedron c 10) `V.snoc`
+    (tetrahedron c 16) `V.snoc`
+    (tetrahedron c 17)
+
+back_right_down_tetrahedra :: Cube -> V.Vector Tetrahedron
+back_right_down_tetrahedra  c =
+  V.singleton (tetrahedron c 8) `V.snoc`
+    (tetrahedron c 9) `V.snoc`
+    (tetrahedron c 13) `V.snoc`
+    (tetrahedron c 14) `V.snoc`
+    (tetrahedron c 17) `V.snoc`
+    (tetrahedron c 18)
 
 in_top_half :: Cube -> Point -> Bool
 in_top_half c (_,_,z) =
@@ -531,33 +583,41 @@ in_left_half c (_,y,_) =
 --
 find_containing_tetrahedron :: Cube -> Point -> Tetrahedron
 find_containing_tetrahedron c p =
-  head containing_tetrahedra
+  candidates `V.unsafeIndex` (fromJust lucky_idx)
   where
-    candidates = tetrahedra c
-    non_candidates_x =
-        if (in_front_half c p) then
-          back_half_tetrahedra c
+    front_half = in_front_half c p
+    top_half = in_top_half c p
+    left_half = in_left_half c p
+
+    candidates =
+      if front_half then
+
+        if left_half then
+          if top_half then
+            front_left_top_tetrahedra c
+          else
+            front_left_down_tetrahedra c
         else
-          front_half_tetrahedra c
-
-    candidates_x = candidates \\ non_candidates_x
-
-    non_candidates_y =
-      if (in_left_half c p) then
-        right_half_tetrahedra c
-      else
-        left_half_tetrahedra c
-
-    candidates_xy = candidates_x \\ non_candidates_y
-
-    non_candidates_z =
-      if (in_top_half c p) then
-        down_half_tetrahedra c
-      else
-        top_half_tetrahedra c
-
-    candidates_xyz = candidates_xy \\ non_candidates_z
-
-    contains_our_point = flip contains_point p
-    containing_tetrahedra = filter contains_our_point candidates_xyz
-
+          if top_half then
+            front_right_top_tetrahedra c
+          else
+            front_right_down_tetrahedra c
+
+      else -- bottom half
+
+        if left_half then
+          if top_half then
+            back_left_top_tetrahedra c
+          else
+            back_left_down_tetrahedra c
+        else
+          if top_half then
+            back_right_top_tetrahedra c
+          else
+            back_right_down_tetrahedra c
+
+    -- Use the dot product instead of 'distance' here to save a
+    -- sqrt(). So, "distances" below really means "distances squared."
+    distances = V.map ((dot p) . center) candidates
+    shortest_distance = V.minimum distances
+    lucky_idx = V.findIndex (\t -> (center t) `dot` p == shortest_distance)