import qualified Face (Face(Face, v0, v1, v2, v3))
import FunctionValues (FunctionValues, eval, rotate)
import Misc (all_equal, disjoint)
-import Point
+import Point (Point(..), dot)
import Tetrahedron (Tetrahedron(..), c, volume)
import ThreeDimensional
instance ThreeDimensional Cube where
-- | The center of Cube_ijk coincides with v_ijk at
-- (ih, jh, kh). See Sorokina and Zeilfelder, p. 76.
- center cube = (x, y, z)
+ center cube = Point x y z
where
delta = h cube
i' = fromIntegral (i cube) :: Double
-- | It's easy to tell if a point is within a cube; just make sure
-- that it falls on the proper side of each of the cube's faces.
- contains_point cube (x, y, z)
+ contains_point cube (Point x y z)
| x < (xmin cube) = False
| x > (xmax cube) = False
| y < (ymin cube) = False
where
delta = (1/2)*(h cube)
cc = center cube
- v0' = cc + (delta, -delta, delta)
- v1' = cc + (delta, delta, delta)
- v2' = cc + (-delta, delta, delta)
- v3' = cc + (-delta, -delta, delta)
+ v0' = cc + ( Point delta (-delta) delta )
+ v1' = cc + ( Point delta delta delta )
+ v2' = cc + ( Point (-delta) delta delta )
+ v3' = cc + ( Point (-delta) (-delta) delta )
where
delta = (1/2)*(h cube)
cc = center cube
- v0' = cc + (delta, -delta, -delta)
- v1' = cc + (delta, delta, -delta)
- v2' = cc + (delta, delta, delta)
- v3' = cc + (delta, -delta, delta)
+ v0' = cc + ( Point delta (-delta) (-delta) )
+ v1' = cc + ( Point delta delta (-delta) )
+ v2' = cc + ( Point delta delta delta )
+ v3' = cc + ( Point delta (-delta) delta )
-- The bottom face (in the direction of -z) of the cube.
where
delta = (1/2)*(h cube)
cc = center cube
- v0' = cc + (-delta, -delta, -delta)
- v1' = cc + (-delta, delta, -delta)
- v2' = cc + (delta, delta, -delta)
- v3' = cc + (delta, -delta, -delta)
+ v0' = cc + ( Point (-delta) (-delta) (-delta) )
+ v1' = cc + ( Point (-delta) delta (-delta) )
+ v2' = cc + ( Point delta delta (-delta) )
+ v3' = cc + ( Point delta (-delta) (-delta) )
where
delta = (1/2)*(h cube)
cc = center cube
- v0' = cc + (-delta, -delta, delta)
- v1' = cc + (-delta, delta, delta)
- v2' = cc + (-delta, delta, -delta)
- v3' = cc + (-delta, -delta, -delta)
+ v0' = cc + ( Point (-delta) (-delta) delta )
+ v1' = cc + ( Point (-delta) delta delta )
+ v2' = cc + ( Point (-delta) delta (-delta) )
+ v3' = cc + ( Point (-delta) (-delta) (-delta) )
-- | The left (in the direction of -y) face of the cube.
left_face :: Cube -> Face.Face
where
delta = (1/2)*(h cube)
cc = center cube
- v0' = cc + (delta, -delta, delta)
- v1' = cc + (-delta, -delta, delta)
- v2' = cc + (-delta, -delta, -delta)
- v3' = cc + (delta, -delta, -delta)
+ v0' = cc + ( Point delta (-delta) delta )
+ v1' = cc + ( Point (-delta) (-delta) delta )
+ v2' = cc + ( Point (-delta) (-delta) (-delta) )
+ v3' = cc + ( Point delta (-delta) (-delta) )
-- | The right (in the direction of y) face of the cube.
where
delta = (1/2)*(h cube)
cc = center cube
- v0' = cc + (-delta, delta, delta)
- v1' = cc + (delta, delta, delta)
- v2' = cc + (delta, delta, -delta)
- v3' = cc + (-delta, delta, -delta)
+ v0' = cc + ( Point (-delta) delta delta)
+ v1' = cc + ( Point delta delta delta )
+ v2' = cc + ( Point delta delta (-delta) )
+ v3' = cc + ( Point (-delta) delta (-delta) )
tetrahedron :: Cube -> Int -> Tetrahedron
(tetrahedron cube 18)
in_top_half :: Cube -> Point -> Bool
-in_top_half cube (_,_,z) =
+in_top_half cube (Point _ _ z) =
distance_from_top <= distance_from_bottom
where
distance_from_top = abs $ (zmax cube) - z
distance_from_bottom = abs $ (zmin cube) - z
in_front_half :: Cube -> Point -> Bool
-in_front_half cube (x,_,_) =
+in_front_half cube (Point x _ _) =
distance_from_front <= distance_from_back
where
distance_from_front = abs $ (xmin cube) - x
in_left_half :: Cube -> Point -> Bool
-in_left_half cube (_,y,_) =
+in_left_half cube (Point _ y _) =
distance_from_left <= distance_from_right
where
distance_from_left = abs $ (ymin cube) - y
-- This can throw an exception, but the use of 'head' might
-- save us some unnecessary computations.
--
+{-# INLINE find_containing_tetrahedron #-}
find_containing_tetrahedron :: Cube -> Point -> Tetrahedron
find_containing_tetrahedron cube p =
candidates `V.unsafeIndex` (fromJust lucky_idx)
else
back_right_down_tetrahedra cube
- -- Use the dot product instead of 'distance' here to save a
- -- sqrt(). So, "distances" below really means "distances squared."
+ -- Use the dot product instead of Euclidean 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