4 find_containing_tetrahedron,
10 import Data.Maybe (fromJust)
11 import qualified Data.Vector as V (
20 import Prelude hiding (LT)
21 import Test.Framework (Test, testGroup)
22 import Test.Framework.Providers.QuickCheck2 (testProperty)
23 import Test.QuickCheck (Arbitrary(..), Gen, Positive(..), choose)
26 import Comparisons ((~=), (~~=))
27 import qualified Face (Face(Face, v0, v1, v2, v3))
29 import Misc (all_equal, disjoint)
31 import Tetrahedron (Tetrahedron(..), c, volume)
32 import ThreeDimensional
34 data Cube = Cube { h :: Double,
39 tetrahedra_volume :: Double }
43 instance Arbitrary Cube where
45 (Positive h') <- arbitrary :: Gen (Positive Double)
46 i' <- choose (coordmin, coordmax)
47 j' <- choose (coordmin, coordmax)
48 k' <- choose (coordmin, coordmax)
49 fv' <- arbitrary :: Gen FunctionValues
50 (Positive tet_vol) <- arbitrary :: Gen (Positive Double)
51 return (Cube h' i' j' k' fv' tet_vol)
53 coordmin = -268435456 -- -(2^29 / 2)
54 coordmax = 268435456 -- +(2^29 / 2)
57 instance Show Cube where
59 "Cube_" ++ subscript ++ "\n" ++
60 " h: " ++ (show (h cube)) ++ "\n" ++
61 " Center: " ++ (show (center cube)) ++ "\n" ++
62 " xmin: " ++ (show (xmin cube)) ++ "\n" ++
63 " xmax: " ++ (show (xmax cube)) ++ "\n" ++
64 " ymin: " ++ (show (ymin cube)) ++ "\n" ++
65 " ymax: " ++ (show (ymax cube)) ++ "\n" ++
66 " zmin: " ++ (show (zmin cube)) ++ "\n" ++
67 " zmax: " ++ (show (zmax cube)) ++ "\n"
70 (show (i cube)) ++ "," ++ (show (j cube)) ++ "," ++ (show (k cube))
73 -- | The left-side boundary of the cube. See Sorokina and Zeilfelder,
75 xmin :: Cube -> Double
76 xmin cube = (i' - 1/2)*delta
78 i' = fromIntegral (i cube) :: Double
81 -- | The right-side boundary of the cube. See Sorokina and Zeilfelder,
83 xmax :: Cube -> Double
84 xmax cube = (i' + 1/2)*delta
86 i' = fromIntegral (i cube) :: Double
89 -- | The front boundary of the cube. See Sorokina and Zeilfelder,
91 ymin :: Cube -> Double
92 ymin cube = (j' - 1/2)*delta
94 j' = fromIntegral (j cube) :: Double
97 -- | The back boundary of the cube. See Sorokina and Zeilfelder,
99 ymax :: Cube -> Double
100 ymax cube = (j' + 1/2)*delta
102 j' = fromIntegral (j cube) :: Double
105 -- | The bottom boundary of the cube. See Sorokina and Zeilfelder,
107 zmin :: Cube -> Double
108 zmin cube = (k' - 1/2)*delta
110 k' = fromIntegral (k cube) :: Double
113 -- | The top boundary of the cube. See Sorokina and Zeilfelder,
115 zmax :: Cube -> Double
116 zmax cube = (k' + 1/2)*delta
118 k' = fromIntegral (k cube) :: Double
121 instance ThreeDimensional Cube where
122 -- | The center of Cube_ijk coincides with v_ijk at
123 -- (ih, jh, kh). See Sorokina and Zeilfelder, p. 76.
124 center cube = (x, y, z)
127 i' = fromIntegral (i cube) :: Double
128 j' = fromIntegral (j cube) :: Double
129 k' = fromIntegral (k cube) :: Double
134 -- | It's easy to tell if a point is within a cube; just make sure
135 -- that it falls on the proper side of each of the cube's faces.
136 contains_point cube (x, y, z)
137 | x < (xmin cube) = False
138 | x > (xmax cube) = False
139 | y < (ymin cube) = False
140 | y > (ymax cube) = False
141 | z < (zmin cube) = False
142 | z > (zmax cube) = False
149 -- | The top (in the direction of z) face of the cube.
150 top_face :: Cube -> Face.Face
151 top_face cube = Face.Face v0' v1' v2' v3'
153 delta = (1/2)*(h cube)
154 v0' = (center cube) + (delta, -delta, delta)
155 v1' = (center cube) + (delta, delta, delta)
156 v2' = (center cube) + (-delta, delta, delta)
157 v3' = (center cube) + (-delta, -delta, delta)
161 -- | The back (in the direction of x) face of the cube.
162 back_face :: Cube -> Face.Face
163 back_face cube = Face.Face v0' v1' v2' v3'
165 delta = (1/2)*(h cube)
166 v0' = (center cube) + (delta, -delta, -delta)
167 v1' = (center cube) + (delta, delta, -delta)
168 v2' = (center cube) + (delta, delta, delta)
169 v3' = (center cube) + (delta, -delta, delta)
172 -- The bottom face (in the direction of -z) of the cube.
173 down_face :: Cube -> Face.Face
174 down_face cube = Face.Face v0' v1' v2' v3'
176 delta = (1/2)*(h cube)
177 v0' = (center cube) + (-delta, -delta, -delta)
178 v1' = (center cube) + (-delta, delta, -delta)
179 v2' = (center cube) + (delta, delta, -delta)
180 v3' = (center cube) + (delta, -delta, -delta)
184 -- | The front (in the direction of -x) face of the cube.
185 front_face :: Cube -> Face.Face
186 front_face cube = Face.Face v0' v1' v2' v3'
188 delta = (1/2)*(h cube)
189 v0' = (center cube) + (-delta, -delta, delta)
190 v1' = (center cube) + (-delta, delta, delta)
191 v2' = (center cube) + (-delta, delta, -delta)
192 v3' = (center cube) + (-delta, -delta, -delta)
194 -- | The left (in the direction of -y) face of the cube.
195 left_face :: Cube -> Face.Face
196 left_face cube = Face.Face v0' v1' v2' v3'
198 delta = (1/2)*(h cube)
199 v0' = (center cube) + (delta, -delta, delta)
200 v1' = (center cube) + (-delta, -delta, delta)
201 v2' = (center cube) + (-delta, -delta, -delta)
202 v3' = (center cube) + (delta, -delta, -delta)
205 -- | The right (in the direction of y) face of the cube.
206 right_face :: Cube -> Face.Face
207 right_face cube = Face.Face v0' v1' v2' v3'
209 delta = (1/2)*(h cube)
210 v0' = (center cube) + (-delta, delta, delta)
211 v1' = (center cube) + (delta, delta, delta)
212 v2' = (center cube) + (delta, delta, -delta)
213 v3' = (center cube) + (-delta, delta, -delta)
216 tetrahedron :: Cube -> Int -> Tetrahedron
219 Tetrahedron (fv cube) v0' v1' v2' v3' vol
222 v1' = center (front_face cube)
223 v2' = Face.v0 (front_face cube)
224 v3' = Face.v1 (front_face cube)
225 vol = tetrahedra_volume cube
228 Tetrahedron fv' v0' v1' v2' v3' vol
231 v1' = center (front_face cube)
232 v2' = Face.v1 (front_face cube)
233 v3' = Face.v2 (front_face cube)
234 fv' = rotate ccwx (fv cube)
235 vol = tetrahedra_volume cube
238 Tetrahedron fv' v0' v1' v2' v3' vol
241 v1' = center (front_face cube)
242 v2' = Face.v2 (front_face cube)
243 v3' = Face.v3 (front_face cube)
244 fv' = rotate ccwx $ rotate ccwx $ fv cube
245 vol = tetrahedra_volume cube
248 Tetrahedron fv' v0' v1' v2' v3' vol
251 v1' = center (front_face cube)
252 v2' = Face.v3 (front_face cube)
253 v3' = Face.v0 (front_face cube)
254 fv' = rotate cwx (fv cube)
255 vol = tetrahedra_volume cube
258 Tetrahedron fv' v0' v1' v2' v3' vol
261 v1' = center (top_face cube)
262 v2' = Face.v0 (top_face cube)
263 v3' = Face.v1 (top_face cube)
264 fv' = rotate cwy (fv cube)
265 vol = tetrahedra_volume cube
268 Tetrahedron fv' v0' v1' v2' v3' vol
271 v1' = center (top_face cube)
272 v2' = Face.v1 (top_face cube)
273 v3' = Face.v2 (top_face cube)
274 fv' = rotate cwy $ rotate cwz $ fv cube
275 vol = tetrahedra_volume cube
278 Tetrahedron fv' v0' v1' v2' v3' vol
281 v1' = center (top_face cube)
282 v2' = Face.v2 (top_face cube)
283 v3' = Face.v3 (top_face cube)
284 fv' = rotate cwy $ rotate cwz
287 vol = tetrahedra_volume cube
290 Tetrahedron fv' v0' v1' v2' v3' vol
293 v1' = center (top_face cube)
294 v2' = Face.v3 (top_face cube)
295 v3' = Face.v0 (top_face cube)
296 fv' = rotate cwy $ rotate ccwz $ fv cube
297 vol = tetrahedra_volume cube
300 Tetrahedron fv' v0' v1' v2' v3' vol
303 v1' = center (back_face cube)
304 v2' = Face.v0 (back_face cube)
305 v3' = Face.v1 (back_face cube)
306 fv' = rotate cwy $ rotate cwy $ fv cube
307 vol = tetrahedra_volume cube
310 Tetrahedron fv' v0' v1' v2' v3' vol
313 v1' = center (back_face cube)
314 v2' = Face.v1 (back_face cube)
315 v3' = Face.v2 (back_face cube)
316 fv' = rotate cwy $ rotate cwy
319 vol = tetrahedra_volume cube
321 tetrahedron cube 10 =
322 Tetrahedron fv' v0' v1' v2' v3' vol
325 v1' = center (back_face cube)
326 v2' = Face.v2 (back_face cube)
327 v3' = Face.v3 (back_face cube)
328 fv' = rotate cwy $ rotate cwy
333 vol = tetrahedra_volume cube
335 tetrahedron cube 11 =
336 Tetrahedron fv' v0' v1' v2' v3' vol
339 v1' = center (back_face cube)
340 v2' = Face.v3 (back_face cube)
341 v3' = Face.v0 (back_face cube)
342 fv' = rotate cwy $ rotate cwy
345 vol = tetrahedra_volume cube
347 tetrahedron cube 12 =
348 Tetrahedron fv' v0' v1' v2' v3' vol
351 v1' = center (down_face cube)
352 v2' = Face.v0 (down_face cube)
353 v3' = Face.v1 (down_face cube)
354 fv' = rotate ccwy $ fv cube
355 vol = tetrahedra_volume cube
357 tetrahedron cube 13 =
358 Tetrahedron fv' v0' v1' v2' v3' vol
361 v1' = center (down_face cube)
362 v2' = Face.v1 (down_face cube)
363 v3' = Face.v2 (down_face cube)
364 fv' = rotate ccwy $ rotate ccwz $ fv cube
365 vol = tetrahedra_volume cube
367 tetrahedron cube 14 =
368 Tetrahedron fv' v0' v1' v2' v3' vol
371 v1' = center (down_face cube)
372 v2' = Face.v2 (down_face cube)
373 v3' = Face.v3 (down_face cube)
374 fv' = rotate ccwy $ rotate ccwz
377 vol = tetrahedra_volume cube
379 tetrahedron cube 15 =
380 Tetrahedron fv' v0' v1' v2' v3' vol
383 v1' = center (down_face cube)
384 v2' = Face.v3 (down_face cube)
385 v3' = Face.v0 (down_face cube)
386 fv' = rotate ccwy $ rotate cwz $ fv cube
387 vol = tetrahedra_volume cube
389 tetrahedron cube 16 =
390 Tetrahedron fv' v0' v1' v2' v3' vol
393 v1' = center (right_face cube)
394 v2' = Face.v0 (right_face cube)
395 v3' = Face.v1 (right_face cube)
396 fv' = rotate ccwz $ fv cube
397 vol = tetrahedra_volume cube
399 tetrahedron cube 17 =
400 Tetrahedron fv' v0' v1' v2' v3' vol
403 v1' = center (right_face cube)
404 v2' = Face.v1 (right_face cube)
405 v3' = Face.v2 (right_face cube)
406 fv' = rotate ccwz $ rotate cwy $ fv cube
407 vol = tetrahedra_volume cube
409 tetrahedron cube 18 =
410 Tetrahedron fv' v0' v1' v2' v3' vol
413 v1' = center (right_face cube)
414 v2' = Face.v2 (right_face cube)
415 v3' = Face.v3 (right_face cube)
416 fv' = rotate ccwz $ rotate cwy
419 vol = tetrahedra_volume cube
421 tetrahedron cube 19 =
422 Tetrahedron fv' v0' v1' v2' v3' vol
425 v1' = center (right_face cube)
426 v2' = Face.v3 (right_face cube)
427 v3' = Face.v0 (right_face cube)
428 fv' = rotate ccwz $ rotate ccwy
430 vol = tetrahedra_volume cube
432 tetrahedron cube 20 =
433 Tetrahedron fv' v0' v1' v2' v3' vol
436 v1' = center (left_face cube)
437 v2' = Face.v0 (left_face cube)
438 v3' = Face.v1 (left_face cube)
439 fv' = rotate cwz $ fv cube
440 vol = tetrahedra_volume cube
442 tetrahedron cube 21 =
443 Tetrahedron fv' v0' v1' v2' v3' vol
446 v1' = center (left_face cube)
447 v2' = Face.v1 (left_face cube)
448 v3' = Face.v2 (left_face cube)
449 fv' = rotate cwz $ rotate ccwy $ fv cube
450 vol = tetrahedra_volume cube
452 tetrahedron cube 22 =
453 Tetrahedron fv' v0' v1' v2' v3' vol
456 v1' = center (left_face cube)
457 v2' = Face.v2 (left_face cube)
458 v3' = Face.v3 (left_face cube)
459 fv' = rotate cwz $ rotate ccwy
462 vol = tetrahedra_volume cube
464 tetrahedron cube 23 =
465 Tetrahedron fv' v0' v1' v2' v3' vol
468 v1' = center (left_face cube)
469 v2' = Face.v3 (left_face cube)
470 v3' = Face.v0 (left_face cube)
471 fv' = rotate cwz $ rotate cwy
473 vol = tetrahedra_volume cube
475 -- Feels dirty, but whatever.
476 tetrahedron _ _ = error "asked for a nonexistent tetrahedron"
479 -- Only used in tests, so we don't need the added speed
481 tetrahedra :: Cube -> [Tetrahedron]
482 tetrahedra cube = [ tetrahedron cube n | n <- [0..23] ]
484 front_left_top_tetrahedra :: Cube -> V.Vector Tetrahedron
485 front_left_top_tetrahedra cube =
486 V.singleton (tetrahedron cube 0) `V.snoc`
487 (tetrahedron cube 3) `V.snoc`
488 (tetrahedron cube 6) `V.snoc`
489 (tetrahedron cube 7) `V.snoc`
490 (tetrahedron cube 20) `V.snoc`
491 (tetrahedron cube 21)
493 front_left_down_tetrahedra :: Cube -> V.Vector Tetrahedron
494 front_left_down_tetrahedra cube =
495 V.singleton (tetrahedron cube 0) `V.snoc`
496 (tetrahedron cube 2) `V.snoc`
497 (tetrahedron cube 3) `V.snoc`
498 (tetrahedron cube 12) `V.snoc`
499 (tetrahedron cube 15) `V.snoc`
500 (tetrahedron cube 21)
502 front_right_top_tetrahedra :: Cube -> V.Vector Tetrahedron
503 front_right_top_tetrahedra cube =
504 V.singleton (tetrahedron cube 0) `V.snoc`
505 (tetrahedron cube 1) `V.snoc`
506 (tetrahedron cube 5) `V.snoc`
507 (tetrahedron cube 6) `V.snoc`
508 (tetrahedron cube 16) `V.snoc`
509 (tetrahedron cube 19)
511 front_right_down_tetrahedra :: Cube -> V.Vector Tetrahedron
512 front_right_down_tetrahedra cube =
513 V.singleton (tetrahedron cube 1) `V.snoc`
514 (tetrahedron cube 2) `V.snoc`
515 (tetrahedron cube 12) `V.snoc`
516 (tetrahedron cube 13) `V.snoc`
517 (tetrahedron cube 18) `V.snoc`
518 (tetrahedron cube 19)
520 back_left_top_tetrahedra :: Cube -> V.Vector Tetrahedron
521 back_left_top_tetrahedra cube =
522 V.singleton (tetrahedron cube 0) `V.snoc`
523 (tetrahedron cube 3) `V.snoc`
524 (tetrahedron cube 6) `V.snoc`
525 (tetrahedron cube 7) `V.snoc`
526 (tetrahedron cube 20) `V.snoc`
527 (tetrahedron cube 21)
529 back_left_down_tetrahedra :: Cube -> V.Vector Tetrahedron
530 back_left_down_tetrahedra cube =
531 V.singleton (tetrahedron cube 8) `V.snoc`
532 (tetrahedron cube 11) `V.snoc`
533 (tetrahedron cube 14) `V.snoc`
534 (tetrahedron cube 15) `V.snoc`
535 (tetrahedron cube 22) `V.snoc`
536 (tetrahedron cube 23)
538 back_right_top_tetrahedra :: Cube -> V.Vector Tetrahedron
539 back_right_top_tetrahedra cube =
540 V.singleton (tetrahedron cube 4) `V.snoc`
541 (tetrahedron cube 5) `V.snoc`
542 (tetrahedron cube 9) `V.snoc`
543 (tetrahedron cube 10) `V.snoc`
544 (tetrahedron cube 16) `V.snoc`
545 (tetrahedron cube 17)
547 back_right_down_tetrahedra :: Cube -> V.Vector Tetrahedron
548 back_right_down_tetrahedra cube =
549 V.singleton (tetrahedron cube 8) `V.snoc`
550 (tetrahedron cube 9) `V.snoc`
551 (tetrahedron cube 13) `V.snoc`
552 (tetrahedron cube 14) `V.snoc`
553 (tetrahedron cube 17) `V.snoc`
554 (tetrahedron cube 18)
556 in_top_half :: Cube -> Point -> Bool
557 in_top_half cube (_,_,z) =
558 distance_from_top <= distance_from_bottom
560 distance_from_top = abs $ (zmax cube) - z
561 distance_from_bottom = abs $ (zmin cube) - z
563 in_front_half :: Cube -> Point -> Bool
564 in_front_half cube (x,_,_) =
565 distance_from_front <= distance_from_back
567 distance_from_front = abs $ (xmin cube) - x
568 distance_from_back = abs $ (xmax cube) - x
571 in_left_half :: Cube -> Point -> Bool
572 in_left_half cube (_,y,_) =
573 distance_from_left <= distance_from_right
575 distance_from_left = abs $ (ymin cube) - y
576 distance_from_right = abs $ (ymax cube) - y
579 -- | Takes a 'Cube', and returns the Tetrahedra belonging to it that
580 -- contain the given 'Point'. This should be faster than checking
581 -- every tetrahedron individually, since we determine which half
582 -- (hemisphere?) of the cube the point lies in three times: once in
583 -- each dimension. This allows us to eliminate non-candidates
586 -- This can throw an exception, but the use of 'head' might
587 -- save us some unnecessary computations.
589 find_containing_tetrahedron :: Cube -> Point -> Tetrahedron
590 find_containing_tetrahedron cube p =
591 candidates `V.unsafeIndex` (fromJust lucky_idx)
593 front_half = in_front_half cube p
594 top_half = in_top_half cube p
595 left_half = in_left_half cube p
602 front_left_top_tetrahedra cube
604 front_left_down_tetrahedra cube
607 front_right_top_tetrahedra cube
609 front_right_down_tetrahedra cube
615 back_left_top_tetrahedra cube
617 back_left_down_tetrahedra cube
620 back_right_top_tetrahedra cube
622 back_right_down_tetrahedra cube
624 -- Use the dot product instead of 'distance' here to save a
625 -- sqrt(). So, "distances" below really means "distances squared."
626 distances = V.map ((dot p) . center) candidates
627 shortest_distance = V.minimum distances
628 lucky_idx = V.findIndex
629 (\t -> (center t) `dot` p == shortest_distance)
641 prop_opposite_octant_tetrahedra_disjoint1 :: Cube -> Bool
642 prop_opposite_octant_tetrahedra_disjoint1 cube =
643 disjoint (front_left_top_tetrahedra cube) (front_right_down_tetrahedra cube)
645 prop_opposite_octant_tetrahedra_disjoint2 :: Cube -> Bool
646 prop_opposite_octant_tetrahedra_disjoint2 cube =
647 disjoint (back_left_top_tetrahedra cube) (back_right_down_tetrahedra cube)
649 prop_opposite_octant_tetrahedra_disjoint3 :: Cube -> Bool
650 prop_opposite_octant_tetrahedra_disjoint3 cube =
651 disjoint (front_left_top_tetrahedra cube) (back_right_top_tetrahedra cube)
653 prop_opposite_octant_tetrahedra_disjoint4 :: Cube -> Bool
654 prop_opposite_octant_tetrahedra_disjoint4 cube =
655 disjoint (front_left_down_tetrahedra cube) (back_right_down_tetrahedra cube)
657 prop_opposite_octant_tetrahedra_disjoint5 :: Cube -> Bool
658 prop_opposite_octant_tetrahedra_disjoint5 cube =
659 disjoint (front_left_top_tetrahedra cube) (back_left_down_tetrahedra cube)
661 prop_opposite_octant_tetrahedra_disjoint6 :: Cube -> Bool
662 prop_opposite_octant_tetrahedra_disjoint6 cube =
663 disjoint (front_right_top_tetrahedra cube) (back_right_down_tetrahedra cube)
666 -- | Since the grid size is necessarily positive, all tetrahedra
667 -- (which comprise cubes of positive volume) must have positive
669 prop_all_volumes_positive :: Cube -> Bool
670 prop_all_volumes_positive cube =
674 volumes = map volume ts
677 -- | In fact, since all of the tetrahedra are identical, we should
678 -- already know their volumes. There's 24 tetrahedra to a cube, so
679 -- we'd expect the volume of each one to be (1/24)*h^3.
680 prop_all_volumes_exact :: Cube -> Bool
681 prop_all_volumes_exact cube =
682 and [volume t ~~= (1/24)*(delta^(3::Int)) | t <- tetrahedra cube]
686 -- | All tetrahedron should have their v0 located at the center of the cube.
687 prop_v0_all_equal :: Cube -> Bool
688 prop_v0_all_equal cube = (v0 t0) == (v0 t1)
690 t0 = head (tetrahedra cube) -- Doesn't matter which two we choose.
691 t1 = head $ tail (tetrahedra cube)
694 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Note that the
695 -- third and fourth indices of c-t3 have been switched. This is
696 -- because we store the triangles oriented such that their volume is
697 -- positive. If T and T-tilde share \<v0,v1,v2\> and v3,v3-tilde point
698 -- in opposite directions, one of them has to have negative volume!
699 prop_c0120_identity1 :: Cube -> Bool
700 prop_c0120_identity1 cube =
701 c t0 0 1 2 0 ~= (c t0 0 0 2 1 + c t3 0 0 1 2) / 2
703 t0 = tetrahedron cube 0
704 t3 = tetrahedron cube 3
707 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
708 -- 'prop_c0120_identity1' with tetrahedrons 1 and 2.
709 prop_c0120_identity2 :: Cube -> Bool
710 prop_c0120_identity2 cube =
711 c t1 0 1 2 0 ~= (c t1 0 0 2 1 + c t0 0 0 1 2) / 2
713 t0 = tetrahedron cube 0
714 t1 = tetrahedron cube 1
716 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
717 -- 'prop_c0120_identity1' with tetrahedrons 1 and 2.
718 prop_c0120_identity3 :: Cube -> Bool
719 prop_c0120_identity3 cube =
720 c t2 0 1 2 0 ~= (c t2 0 0 2 1 + c t1 0 0 1 2) / 2
722 t1 = tetrahedron cube 1
723 t2 = tetrahedron cube 2
725 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
726 -- 'prop_c0120_identity1' with tetrahedrons 2 and 3.
727 prop_c0120_identity4 :: Cube -> Bool
728 prop_c0120_identity4 cube =
729 c t3 0 1 2 0 ~= (c t3 0 0 2 1 + c t2 0 0 1 2) / 2
731 t2 = tetrahedron cube 2
732 t3 = tetrahedron cube 3
735 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
736 -- 'prop_c0120_identity1' with tetrahedrons 4 and 5.
737 prop_c0120_identity5 :: Cube -> Bool
738 prop_c0120_identity5 cube =
739 c t5 0 1 2 0 ~= (c t5 0 0 2 1 + c t4 0 0 1 2) / 2
741 t4 = tetrahedron cube 4
742 t5 = tetrahedron cube 5
744 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
745 -- 'prop_c0120_identity1' with tetrahedrons 5 and 6.
746 prop_c0120_identity6 :: Cube -> Bool
747 prop_c0120_identity6 cube =
748 c t6 0 1 2 0 ~= (c t6 0 0 2 1 + c t5 0 0 1 2) / 2
750 t5 = tetrahedron cube 5
751 t6 = tetrahedron cube 6
754 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
755 -- 'prop_c0120_identity1' with tetrahedrons 6 and 7.
756 prop_c0120_identity7 :: Cube -> Bool
757 prop_c0120_identity7 cube =
758 c t7 0 1 2 0 ~= (c t7 0 0 2 1 + c t6 0 0 1 2) / 2
760 t6 = tetrahedron cube 6
761 t7 = tetrahedron cube 7
764 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). See
765 -- 'prop_c0120_identity1'.
766 prop_c0210_identity1 :: Cube -> Bool
767 prop_c0210_identity1 cube =
768 c t0 0 2 1 0 ~= (c t0 0 1 1 1 + c t3 0 1 1 1) / 2
770 t0 = tetrahedron cube 0
771 t3 = tetrahedron cube 3
774 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). See
775 -- 'prop_c0120_identity1'.
776 prop_c0300_identity1 :: Cube -> Bool
777 prop_c0300_identity1 cube =
778 c t0 0 3 0 0 ~= (c t0 0 2 0 1 + c t3 0 2 1 0) / 2
780 t0 = tetrahedron cube 0
781 t3 = tetrahedron cube 3
784 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). See
785 -- 'prop_c0120_identity1'.
786 prop_c1110_identity :: Cube -> Bool
787 prop_c1110_identity cube =
788 c t0 1 1 1 0 ~= (c t0 1 0 1 1 + c t3 1 0 1 1) / 2
790 t0 = tetrahedron cube 0
791 t3 = tetrahedron cube 3
794 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). See
795 -- 'prop_c0120_identity1'.
796 prop_c1200_identity1 :: Cube -> Bool
797 prop_c1200_identity1 cube =
798 c t0 1 2 0 0 ~= (c t0 1 1 0 1 + c t3 1 1 1 0) / 2
800 t0 = tetrahedron cube 0
801 t3 = tetrahedron cube 3
804 -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). See
805 -- 'prop_c0120_identity1'.
806 prop_c2100_identity1 :: Cube -> Bool
807 prop_c2100_identity1 cube =
808 c t0 2 1 0 0 ~= (c t0 2 0 0 1 + c t3 2 0 1 0) / 2
810 t0 = tetrahedron cube 0
811 t3 = tetrahedron cube 3
815 -- | Given in Sorokina and Zeilfelder, p. 79, (2.7). Note that the
816 -- third and fourth indices of c-t3 have been switched. This is
817 -- because we store the triangles oriented such that their volume is
818 -- positive. If T and T-tilde share \<v0,v1,v2\> and v3,v3-tilde
819 -- point in opposite directions, one of them has to have negative
821 prop_c0102_identity1 :: Cube -> Bool
822 prop_c0102_identity1 cube =
823 c t0 0 1 0 2 ~= (c t0 0 0 1 2 + c t1 0 0 2 1) / 2
825 t0 = tetrahedron cube 0
826 t1 = tetrahedron cube 1
829 -- | Given in Sorokina and Zeilfelder, p. 79, (2.7). See
830 -- 'prop_c0102_identity1'.
831 prop_c0201_identity1 :: Cube -> Bool
832 prop_c0201_identity1 cube =
833 c t0 0 2 0 1 ~= (c t0 0 1 1 1 + c t1 0 1 1 1) / 2
835 t0 = tetrahedron cube 0
836 t1 = tetrahedron cube 1
839 -- | Given in Sorokina and Zeilfelder, p. 79, (2.7). See
840 -- 'prop_c0102_identity1'.
841 prop_c0300_identity2 :: Cube -> Bool
842 prop_c0300_identity2 cube =
843 c t0 0 3 0 0 ~= (c t0 0 2 1 0 + c t1 0 2 0 1) / 2
845 t0 = tetrahedron cube 0
846 t1 = tetrahedron cube 1
849 -- | Given in Sorokina and Zeilfelder, p. 79, (2.7). See
850 -- 'prop_c0102_identity1'.
851 prop_c1101_identity :: Cube -> Bool
852 prop_c1101_identity cube =
853 c t0 1 1 0 1 ~= (c t0 1 0 1 1 + c t1 1 0 1 1) / 2
855 t0 = tetrahedron cube 0
856 t1 = tetrahedron cube 1
859 -- | Given in Sorokina and Zeilfelder, p. 79, (2.7). See
860 -- 'prop_c0102_identity1'.
861 prop_c1200_identity2 :: Cube -> Bool
862 prop_c1200_identity2 cube =
863 c t0 1 2 0 0 ~= (c t0 1 1 1 0 + c t1 1 1 0 1) / 2
865 t0 = tetrahedron cube 0
866 t1 = tetrahedron cube 1
869 -- | Given in Sorokina and Zeilfelder, p. 79, (2.7). See
870 -- 'prop_c0102_identity1'.
871 prop_c2100_identity2 :: Cube -> Bool
872 prop_c2100_identity2 cube =
873 c t0 2 1 0 0 ~= (c t0 2 0 1 0 + c t1 2 0 0 1) / 2
875 t0 = tetrahedron cube 0
876 t1 = tetrahedron cube 1
879 -- | Given in Sorokina and Zeilfelder, p. 79, (2.8). The third and
880 -- fourth indices of c-t6 have been switched. This is because we
881 -- store the triangles oriented such that their volume is
882 -- positive. If T and T-tilde share \<v0,v1,v2\> and v3,v3-tilde
883 -- point in opposite directions, one of them has to have negative
885 prop_c3000_identity :: Cube -> Bool
886 prop_c3000_identity cube =
887 c t0 3 0 0 0 ~= c t0 2 1 0 0 + c t6 2 1 0 0
888 - ((c t0 2 0 1 0 + c t0 2 0 0 1)/ 2)
890 t0 = tetrahedron cube 0
891 t6 = tetrahedron cube 6
894 -- | Given in Sorokina and Zeilfelder, p. 79, (2.8). See
895 -- 'prop_c3000_identity'.
896 prop_c2010_identity :: Cube -> Bool
897 prop_c2010_identity cube =
898 c t0 2 0 1 0 ~= c t0 1 1 1 0 + c t6 1 1 0 1
899 - ((c t0 1 0 2 0 + c t0 1 0 1 1)/ 2)
901 t0 = tetrahedron cube 0
902 t6 = tetrahedron cube 6
905 -- | Given in Sorokina and Zeilfelder, p. 79, (2.8). See
906 -- 'prop_c3000_identity'.
907 prop_c2001_identity :: Cube -> Bool
908 prop_c2001_identity cube =
909 c t0 2 0 0 1 ~= c t0 1 1 0 1 + c t6 1 1 1 0
910 - ((c t0 1 0 0 2 + c t0 1 0 1 1)/ 2)
912 t0 = tetrahedron cube 0
913 t6 = tetrahedron cube 6
916 -- | Given in Sorokina and Zeilfelder, p. 79, (2.8). See
917 -- 'prop_c3000_identity'.
918 prop_c1020_identity :: Cube -> Bool
919 prop_c1020_identity cube =
920 c t0 1 0 2 0 ~= c t0 0 1 2 0 + c t6 0 1 0 2
921 - ((c t0 0 0 3 0 + c t0 0 0 2 1)/ 2)
923 t0 = tetrahedron cube 0
924 t6 = tetrahedron cube 6
927 -- | Given in Sorokina and Zeilfelder, p. 79, (2.8). See
928 -- 'prop_c3000_identity'.
929 prop_c1002_identity :: Cube -> Bool
930 prop_c1002_identity cube =
931 c t0 1 0 0 2 ~= c t0 0 1 0 2 + c t6 0 1 2 0
932 - ((c t0 0 0 0 3 + c t0 0 0 1 2)/ 2)
934 t0 = tetrahedron cube 0
935 t6 = tetrahedron cube 6
938 -- | Given in Sorokina and Zeilfelder, p. 79, (2.8). See
939 -- 'prop_c3000_identity'.
940 prop_c1011_identity :: Cube -> Bool
941 prop_c1011_identity cube =
942 c t0 1 0 1 1 ~= c t0 0 1 1 1 + c t6 0 1 1 1 -
943 ((c t0 0 0 1 2 + c t0 0 0 2 1)/ 2)
945 t0 = tetrahedron cube 0
946 t6 = tetrahedron cube 6
949 -- | The function values at the interior should be the same for all
951 prop_interior_values_all_identical :: Cube -> Bool
952 prop_interior_values_all_identical cube =
953 all_equal [ eval (function_values tet) I | tet <- tetrahedra cube ]
956 -- | We know what (c t6 2 1 0 0) should be from Sorokina and Zeilfelder, p. 87.
957 -- This test checks the rotation works as expected.
958 prop_c_tilde_2100_rotation_correct :: Cube -> Bool
959 prop_c_tilde_2100_rotation_correct cube =
962 t0 = tetrahedron cube 0
963 t6 = tetrahedron cube 6
965 -- What gets computed for c2100 of t6.
966 expr1 = eval (function_values t6) $
968 (1/12)*(T + R + L + D) +
969 (1/64)*(FT + FR + FL + FD) +
972 (1/96)*(RT + LD + LT + RD) +
973 (1/192)*(BT + BR + BL + BD)
975 -- What should be computed for c2100 of t6.
976 expr2 = eval (function_values t0) $
978 (1/12)*(F + R + L + B) +
979 (1/64)*(FT + RT + LT + BT) +
982 (1/96)*(FR + FL + BR + BL) +
983 (1/192)*(FD + RD + LD + BD)
986 -- | We know what (c t6 2 1 0 0) should be from Sorokina and
987 -- Zeilfelder, p. 87. This test checks the actual value based on
988 -- the FunctionValues of the cube.
990 -- If 'prop_c_tilde_2100_rotation_correct' passes, then this test is
992 prop_c_tilde_2100_correct :: Cube -> Bool
993 prop_c_tilde_2100_correct cube =
994 c t6 2 1 0 0 == expected
996 t0 = tetrahedron cube 0
997 t6 = tetrahedron cube 6
998 fvs = function_values t0
999 expected = eval fvs $
1001 (1/12)*(F + R + L + B) +
1002 (1/64)*(FT + RT + LT + BT) +
1005 (1/96)*(FR + FL + BR + BL) +
1006 (1/192)*(FD + RD + LD + BD)
1009 -- Tests to check that the correct edges are incidental.
1010 prop_t0_shares_edge_with_t1 :: Cube -> Bool
1011 prop_t0_shares_edge_with_t1 cube =
1012 (v1 t0) == (v1 t1) && (v3 t0) == (v2 t1)
1014 t0 = tetrahedron cube 0
1015 t1 = tetrahedron cube 1
1017 prop_t0_shares_edge_with_t3 :: Cube -> Bool
1018 prop_t0_shares_edge_with_t3 cube =
1019 (v1 t0) == (v1 t3) && (v2 t0) == (v3 t3)
1021 t0 = tetrahedron cube 0
1022 t3 = tetrahedron cube 3
1024 prop_t0_shares_edge_with_t6 :: Cube -> Bool
1025 prop_t0_shares_edge_with_t6 cube =
1026 (v2 t0) == (v3 t6) && (v3 t0) == (v2 t6)
1028 t0 = tetrahedron cube 0
1029 t6 = tetrahedron cube 6
1031 prop_t1_shares_edge_with_t2 :: Cube -> Bool
1032 prop_t1_shares_edge_with_t2 cube =
1033 (v1 t1) == (v1 t2) && (v3 t1) == (v2 t2)
1035 t1 = tetrahedron cube 1
1036 t2 = tetrahedron cube 2
1038 prop_t1_shares_edge_with_t19 :: Cube -> Bool
1039 prop_t1_shares_edge_with_t19 cube =
1040 (v2 t1) == (v3 t19) && (v3 t1) == (v2 t19)
1042 t1 = tetrahedron cube 1
1043 t19 = tetrahedron cube 19
1045 prop_t2_shares_edge_with_t3 :: Cube -> Bool
1046 prop_t2_shares_edge_with_t3 cube =
1047 (v1 t1) == (v1 t2) && (v3 t1) == (v2 t2)
1049 t1 = tetrahedron cube 1
1050 t2 = tetrahedron cube 2
1052 prop_t2_shares_edge_with_t12 :: Cube -> Bool
1053 prop_t2_shares_edge_with_t12 cube =
1054 (v2 t2) == (v3 t12) && (v3 t2) == (v2 t12)
1056 t2 = tetrahedron cube 2
1057 t12 = tetrahedron cube 12
1059 prop_t3_shares_edge_with_t21 :: Cube -> Bool
1060 prop_t3_shares_edge_with_t21 cube =
1061 (v2 t3) == (v3 t21) && (v3 t3) == (v2 t21)
1063 t3 = tetrahedron cube 3
1064 t21 = tetrahedron cube 21
1066 prop_t4_shares_edge_with_t5 :: Cube -> Bool
1067 prop_t4_shares_edge_with_t5 cube =
1068 (v1 t4) == (v1 t5) && (v3 t4) == (v2 t5)
1070 t4 = tetrahedron cube 4
1071 t5 = tetrahedron cube 5
1073 prop_t4_shares_edge_with_t7 :: Cube -> Bool
1074 prop_t4_shares_edge_with_t7 cube =
1075 (v1 t4) == (v1 t7) && (v2 t4) == (v3 t7)
1077 t4 = tetrahedron cube 4
1078 t7 = tetrahedron cube 7
1080 prop_t4_shares_edge_with_t10 :: Cube -> Bool
1081 prop_t4_shares_edge_with_t10 cube =
1082 (v2 t4) == (v3 t10) && (v3 t4) == (v2 t10)
1084 t4 = tetrahedron cube 4
1085 t10 = tetrahedron cube 10
1087 prop_t5_shares_edge_with_t6 :: Cube -> Bool
1088 prop_t5_shares_edge_with_t6 cube =
1089 (v1 t5) == (v1 t6) && (v3 t5) == (v2 t6)
1091 t5 = tetrahedron cube 5
1092 t6 = tetrahedron cube 6
1094 prop_t5_shares_edge_with_t16 :: Cube -> Bool
1095 prop_t5_shares_edge_with_t16 cube =
1096 (v2 t5) == (v3 t16) && (v3 t5) == (v2 t16)
1098 t5 = tetrahedron cube 5
1099 t16 = tetrahedron cube 16
1101 prop_t6_shares_edge_with_t7 :: Cube -> Bool
1102 prop_t6_shares_edge_with_t7 cube =
1103 (v1 t6) == (v1 t7) && (v3 t6) == (v2 t7)
1105 t6 = tetrahedron cube 6
1106 t7 = tetrahedron cube 7
1108 prop_t7_shares_edge_with_t20 :: Cube -> Bool
1109 prop_t7_shares_edge_with_t20 cube =
1110 (v2 t7) == (v3 t20) && (v2 t7) == (v3 t20)
1112 t7 = tetrahedron cube 7
1113 t20 = tetrahedron cube 20
1116 p79_26_properties :: Test.Framework.Test
1118 testGroup "p. 79, Section (2.6) Properties" [
1119 testProperty "c0120 identity1" prop_c0120_identity1,
1120 testProperty "c0120 identity2" prop_c0120_identity2,
1121 testProperty "c0120 identity3" prop_c0120_identity3,
1122 testProperty "c0120 identity4" prop_c0120_identity4,
1123 testProperty "c0120 identity5" prop_c0120_identity5,
1124 testProperty "c0120 identity6" prop_c0120_identity6,
1125 testProperty "c0120 identity7" prop_c0120_identity7,
1126 testProperty "c0210 identity1" prop_c0210_identity1,
1127 testProperty "c0300 identity1" prop_c0300_identity1,
1128 testProperty "c1110 identity" prop_c1110_identity,
1129 testProperty "c1200 identity1" prop_c1200_identity1,
1130 testProperty "c2100 identity1" prop_c2100_identity1]
1132 p79_27_properties :: Test.Framework.Test
1134 testGroup "p. 79, Section (2.7) Properties" [
1135 testProperty "c0102 identity1" prop_c0102_identity1,
1136 testProperty "c0201 identity1" prop_c0201_identity1,
1137 testProperty "c0300 identity2" prop_c0300_identity2,
1138 testProperty "c1101 identity" prop_c1101_identity,
1139 testProperty "c1200 identity2" prop_c1200_identity2,
1140 testProperty "c2100 identity2" prop_c2100_identity2 ]
1143 p79_28_properties :: Test.Framework.Test
1145 testGroup "p. 79, Section (2.8) Properties" [
1146 testProperty "c3000 identity" prop_c3000_identity,
1147 testProperty "c2010 identity" prop_c2010_identity,
1148 testProperty "c2001 identity" prop_c2001_identity,
1149 testProperty "c1020 identity" prop_c1020_identity,
1150 testProperty "c1002 identity" prop_c1002_identity,
1151 testProperty "c1011 identity" prop_c1011_identity ]
1154 edge_incidence_tests :: Test.Framework.Test
1155 edge_incidence_tests =
1156 testGroup "Edge Incidence Tests" [
1157 testProperty "t0 shares edge with t6" prop_t0_shares_edge_with_t6,
1158 testProperty "t0 shares edge with t1" prop_t0_shares_edge_with_t1,
1159 testProperty "t0 shares edge with t3" prop_t0_shares_edge_with_t3,
1160 testProperty "t1 shares edge with t2" prop_t1_shares_edge_with_t2,
1161 testProperty "t1 shares edge with t19" prop_t1_shares_edge_with_t19,
1162 testProperty "t2 shares edge with t3" prop_t2_shares_edge_with_t3,
1163 testProperty "t2 shares edge with t12" prop_t2_shares_edge_with_t12,
1164 testProperty "t3 shares edge with t21" prop_t3_shares_edge_with_t21,
1165 testProperty "t4 shares edge with t5" prop_t4_shares_edge_with_t5,
1166 testProperty "t4 shares edge with t7" prop_t4_shares_edge_with_t7,
1167 testProperty "t4 shares edge with t10" prop_t4_shares_edge_with_t10,
1168 testProperty "t5 shares edge with t6" prop_t5_shares_edge_with_t6,
1169 testProperty "t5 shares edge with t16" prop_t5_shares_edge_with_t16,
1170 testProperty "t6 shares edge with t7" prop_t6_shares_edge_with_t7,
1171 testProperty "t7 shares edge with t20" prop_t7_shares_edge_with_t20 ]
1173 cube_properties :: Test.Framework.Test
1175 testGroup "Cube Properties" [
1179 edge_incidence_tests,
1180 testProperty "opposite octant tetrahedra are disjoint (1)"
1181 prop_opposite_octant_tetrahedra_disjoint1,
1182 testProperty "opposite octant tetrahedra are disjoint (2)"
1183 prop_opposite_octant_tetrahedra_disjoint2,
1184 testProperty "opposite octant tetrahedra are disjoint (3)"
1185 prop_opposite_octant_tetrahedra_disjoint3,
1186 testProperty "opposite octant tetrahedra are disjoint (4)"
1187 prop_opposite_octant_tetrahedra_disjoint4,
1188 testProperty "opposite octant tetrahedra are disjoint (5)"
1189 prop_opposite_octant_tetrahedra_disjoint5,
1190 testProperty "opposite octant tetrahedra are disjoint (6)"
1191 prop_opposite_octant_tetrahedra_disjoint6,
1192 testProperty "all volumes positive" prop_all_volumes_positive,
1193 testProperty "all volumes exact" prop_all_volumes_exact,
1194 testProperty "v0 all equal" prop_v0_all_equal,
1195 testProperty "interior values all identical"
1196 prop_interior_values_all_identical,
1197 testProperty "c-tilde_2100 rotation correct"
1198 prop_c_tilde_2100_rotation_correct,
1199 testProperty "c-tilde_2100 correct"
1200 prop_c_tilde_2100_correct ]