4 import Test.QuickCheck (Arbitrary(..), Gen, Positive(..), choose)
7 import qualified Face (Face(Face, v0, v1, v2, v3))
10 import Tetrahedron hiding (c)
11 import ThreeDimensional
13 data Cube = Cube { h :: Double,
17 fv :: FunctionValues }
21 instance Arbitrary Cube where
23 (Positive h') <- arbitrary :: Gen (Positive Double)
24 i' <- choose (coordmin, coordmax)
25 j' <- choose (coordmin, coordmax)
26 k' <- choose (coordmin, coordmax)
27 fv' <- arbitrary :: Gen FunctionValues
28 return (Cube h' i' j' k' fv')
30 coordmin = -268435456 -- -(2^29 / 2)
31 coordmax = 268435456 -- +(2^29 / 2)
34 instance Show Cube where
36 "Cube_" ++ subscript ++ "\n" ++
37 " h: " ++ (show (h c)) ++ "\n" ++
38 " Center: " ++ (show (center c)) ++ "\n" ++
39 " xmin: " ++ (show (xmin c)) ++ "\n" ++
40 " xmax: " ++ (show (xmax c)) ++ "\n" ++
41 " ymin: " ++ (show (ymin c)) ++ "\n" ++
42 " ymax: " ++ (show (ymax c)) ++ "\n" ++
43 " zmin: " ++ (show (zmin c)) ++ "\n" ++
44 " zmax: " ++ (show (zmax c)) ++ "\n" ++
45 " fv: " ++ (show (Cube.fv c)) ++ "\n"
48 (show (i c)) ++ "," ++ (show (j c)) ++ "," ++ (show (k c))
51 -- | Returns an empty 'Cube'.
53 empty_cube = Cube 0 0 0 0 empty_values
56 -- | The left-side boundary of the cube. See Sorokina and Zeilfelder,
58 xmin :: Cube -> Double
59 xmin c = (2*i' - 1)*delta / 2
61 i' = fromIntegral (i c) :: Double
64 -- | The right-side boundary of the cube. See Sorokina and Zeilfelder,
66 xmax :: Cube -> Double
67 xmax c = (2*i' + 1)*delta / 2
69 i' = fromIntegral (i c) :: Double
72 -- | The front boundary of the cube. See Sorokina and Zeilfelder,
74 ymin :: Cube -> Double
75 ymin c = (2*j' - 1)*delta / 2
77 j' = fromIntegral (j c) :: Double
80 -- | The back boundary of the cube. See Sorokina and Zeilfelder,
82 ymax :: Cube -> Double
83 ymax c = (2*j' + 1)*delta / 2
85 j' = fromIntegral (j c) :: Double
88 -- | The bottom boundary of the cube. See Sorokina and Zeilfelder,
90 zmin :: Cube -> Double
91 zmin c = (2*k' - 1)*delta / 2
93 k' = fromIntegral (k c) :: Double
96 -- | The top boundary of the cube. See Sorokina and Zeilfelder,
98 zmax :: Cube -> Double
99 zmax c = (2*k' + 1)*delta / 2
101 k' = fromIntegral (k c) :: Double
104 instance ThreeDimensional Cube where
105 -- | The center of Cube_ijk coincides with v_ijk at
106 -- (ih, jh, kh). See Sorokina and Zeilfelder, p. 76.
110 i' = fromIntegral (i c) :: Double
111 j' = fromIntegral (j c) :: Double
112 k' = fromIntegral (k c) :: Double
117 -- | It's easy to tell if a point is within a cube; just make sure
118 -- that it falls on the proper side of each of the cube's faces.
120 | (x_coord p) < (xmin c) = False
121 | (x_coord p) > (xmax c) = False
122 | (y_coord p) < (ymin c) = False
123 | (y_coord p) > (ymax c) = False
124 | (z_coord p) < (zmin c) = False
125 | (z_coord p) > (zmax c) = False
132 -- | The top (in the direction of z) face of the cube.
133 top_face :: Cube -> Face.Face
134 top_face c = Face.Face v0' v1' v2' v3'
137 v0' = (center c) + (delta, -delta, delta)
138 v1' = (center c) + (delta, delta, delta)
139 v2' = (center c) + (-delta, delta, delta)
140 v3' = (center c) + (-delta, -delta, delta)
144 -- | The back (in the direction of x) face of the cube.
145 back_face :: Cube -> Face.Face
146 back_face c = Face.Face v0' v1' v2' v3'
149 v0' = (center c) + (delta, -delta, -delta)
150 v1' = (center c) + (delta, delta, -delta)
151 v2' = (center c) + (delta, delta, delta)
152 v3' = (center c) + (delta, -delta, delta)
155 -- The bottom face (in the direction of -z) of the cube.
156 down_face :: Cube -> Face.Face
157 down_face c = Face.Face v0' v1' v2' v3'
160 v0' = (center c) + (-delta, -delta, -delta)
161 v1' = (center c) + (-delta, delta, -delta)
162 v2' = (center c) + (delta, delta, -delta)
163 v3' = (center c) + (delta, -delta, -delta)
167 -- | The front (in the direction of -x) face of the cube.
168 front_face :: Cube -> Face.Face
169 front_face c = Face.Face v0' v1' v2' v3'
172 v0' = (center c) + (-delta, -delta, delta)
173 v1' = (center c) + (-delta, delta, delta)
174 v2' = (center c) + (-delta, delta, -delta)
175 v3' = (center c) + (-delta, -delta, -delta)
177 -- | The left (in the direction of -y) face of the cube.
178 left_face :: Cube -> Face.Face
179 left_face c = Face.Face v0' v1' v2' v3'
182 v0' = (center c) + (delta, -delta, delta)
183 v1' = (center c) + (-delta, -delta, delta)
184 v2' = (center c) + (-delta, -delta, -delta)
185 v3' = (center c) + (delta, -delta, -delta)
188 -- | The right (in the direction of y) face of the cube.
189 right_face :: Cube -> Face.Face
190 right_face c = Face.Face v0' v1' v2' v3'
193 v0' = (center c) + (-delta, delta, delta)
194 v1' = (center c) + (delta, delta, delta)
195 v2' = (center c) + (delta, delta, -delta)
196 v3' = (center c) + (-delta, delta, -delta)
199 tetrahedron0 :: Cube -> Tetrahedron
201 Tetrahedron (Cube.fv c) v0' v1' v2' v3'
204 v1' = center (front_face c)
205 v2' = Face.v0 (front_face c)
206 v3' = Face.v1 (front_face c)
208 tetrahedron1 :: Cube -> Tetrahedron
210 Tetrahedron fv' v0' v1' v2' v3'
213 v1' = center (front_face c)
214 v2' = Face.v1 (front_face c)
215 v3' = Face.v2 (front_face c)
216 fv' = rotate ccwx (Cube.fv c)
218 tetrahedron2 :: Cube -> Tetrahedron
220 Tetrahedron fv' v0' v1' v2' v3'
223 v1' = center (front_face c)
224 v2' = Face.v2 (front_face c)
225 v3' = Face.v3 (front_face c)
226 fv' = rotate ccwx $ rotate ccwx $ Cube.fv c
228 tetrahedron3 :: Cube -> Tetrahedron
230 Tetrahedron fv' v0' v1' v2' v3'
233 v1' = center (front_face c)
234 v2' = Face.v3 (front_face c)
235 v3' = Face.v0 (front_face c)
236 fv' = rotate cwx (Cube.fv c)
238 tetrahedron4 :: Cube -> Tetrahedron
240 Tetrahedron fv' v0' v1' v2' v3'
243 v1' = center (top_face c)
244 v2' = Face.v0 (top_face c)
245 v3' = Face.v1 (top_face c)
246 fv' = rotate cwy (Cube.fv c)
248 tetrahedron5 :: Cube -> Tetrahedron
250 Tetrahedron fv' v0' v1' v2' v3'
253 v1' = center (top_face c)
254 v2' = Face.v1 (top_face c)
255 v3' = Face.v2 (top_face c)
256 fv' = rotate cwy $ rotate cwz $ Tetrahedron.fv (tetrahedron0 c)
258 tetrahedron6 :: Cube -> Tetrahedron
260 Tetrahedron fv' v0' v1' v2' v3'
263 v1' = center (top_face c)
264 v2' = Face.v2 (top_face c)
265 v3' = Face.v3 (top_face c)
266 fv' = rotate cwy $ rotate cwz $ rotate cwz $ Tetrahedron.fv (tetrahedron0 c)
268 tetrahedron7 :: Cube -> Tetrahedron
270 Tetrahedron fv' v0' v1' v2' v3'
273 v1' = center (top_face c)
274 v2' = Face.v3 (top_face c)
275 v3' = Face.v0 (top_face c)
276 fv' = rotate cwy $ rotate ccwz $ Tetrahedron.fv (tetrahedron0 c)
278 tetrahedron8 :: Cube -> Tetrahedron
280 Tetrahedron fv' v0' v1' v2' v3'
283 v1' = center (back_face c)
284 v2' = Face.v0 (back_face c)
285 v3' = Face.v1 (back_face c)
286 fv' = rotate cwy $ rotate cwy $ (Tetrahedron.fv (tetrahedron0 c))
288 tetrahedron9 :: Cube -> Tetrahedron
290 Tetrahedron fv' v0' v1' v2' v3'
293 v1' = center (back_face c)
294 v2' = Face.v1 (back_face c)
295 v3' = Face.v2 (back_face c)
296 fv' = rotate cwy $ rotate cwy $ rotate cwx $ Tetrahedron.fv (tetrahedron0 c)
298 tetrahedron10 :: Cube -> Tetrahedron
300 Tetrahedron fv' v0' v1' v2' v3'
303 v1' = center (back_face c)
304 v2' = Face.v2 (back_face c)
305 v3' = Face.v3 (back_face c)
306 fv' = rotate cwy $ rotate cwy
309 $ Tetrahedron.fv (tetrahedron0 c)
312 tetrahedron11 :: Cube -> Tetrahedron
314 Tetrahedron fv' v0' v1' v2' v3'
317 v1' = center (back_face c)
318 v2' = Face.v3 (back_face c)
319 v3' = Face.v0 (back_face c)
320 fv' = rotate cwy $ rotate cwy
322 $ Tetrahedron.fv (tetrahedron0 c)
325 tetrahedron12 :: Cube -> Tetrahedron
327 Tetrahedron fv' v0' v1' v2' v3'
330 v1' = center (down_face c)
331 v2' = Face.v0 (down_face c)
332 v3' = Face.v1 (down_face c)
333 fv' = rotate ccwy (Tetrahedron.fv (tetrahedron0 c))
336 tetrahedron13 :: Cube -> Tetrahedron
338 Tetrahedron fv' v0' v1' v2' v3'
341 v1' = center (down_face c)
342 v2' = Face.v1 (down_face c)
343 v3' = Face.v2 (down_face c)
344 fv' = rotate ccwy $ rotate ccwz $ Tetrahedron.fv (tetrahedron0 c)
347 tetrahedron14 :: Cube -> Tetrahedron
349 Tetrahedron fv' v0' v1' v2' v3'
352 v1' = center (down_face c)
353 v2' = Face.v2 (down_face c)
354 v3' = Face.v3 (down_face c)
355 fv' = rotate ccwy $ rotate ccwz
357 $ Tetrahedron.fv (tetrahedron0 c)
360 tetrahedron15 :: Cube -> Tetrahedron
362 Tetrahedron fv' v0' v1' v2' v3'
365 v1' = center (down_face c)
366 v2' = Face.v3 (down_face c)
367 v3' = Face.v0 (down_face c)
368 fv' = rotate ccwy $ rotate cwz $ Tetrahedron.fv (tetrahedron0 c)
371 tetrahedron16 :: Cube -> Tetrahedron
373 Tetrahedron fv' v0' v1' v2' v3'
376 v1' = center (right_face c)
377 v2' = Face.v0 (right_face c)
378 v3' = Face.v1 (right_face c)
379 fv' = rotate ccwz (Tetrahedron.fv (tetrahedron0 c))
382 tetrahedron17 :: Cube -> Tetrahedron
384 Tetrahedron fv' v0' v1' v2' v3'
387 v1' = center (right_face c)
388 v2' = Face.v1 (right_face c)
389 v3' = Face.v2 (right_face c)
390 fv' = rotate ccwz $ rotate cwy $ Tetrahedron.fv (tetrahedron0 c)
393 tetrahedron18 :: Cube -> Tetrahedron
395 Tetrahedron fv' v0' v1' v2' v3'
398 v1' = center (right_face c)
399 v2' = Face.v2 (right_face c)
400 v3' = Face.v3 (right_face c)
401 fv' = rotate ccwz $ rotate cwy
403 $ Tetrahedron.fv (tetrahedron0 c)
406 tetrahedron19 :: Cube -> Tetrahedron
408 Tetrahedron fv' v0' v1' v2' v3'
411 v1' = center (right_face c)
412 v2' = Face.v3 (right_face c)
413 v3' = Face.v0 (right_face c)
414 fv' = rotate ccwz $ rotate ccwy
415 $ Tetrahedron.fv (tetrahedron0 c)
418 tetrahedron20 :: Cube -> Tetrahedron
420 Tetrahedron fv' v0' v1' v2' v3'
423 v1' = center (left_face c)
424 v2' = Face.v0 (left_face c)
425 v3' = Face.v1 (left_face c)
426 fv' = rotate cwz (Tetrahedron.fv (tetrahedron0 c))
429 tetrahedron21 :: Cube -> Tetrahedron
431 Tetrahedron fv' v0' v1' v2' v3'
434 v1' = center (left_face c)
435 v2' = Face.v1 (left_face c)
436 v3' = Face.v2 (left_face c)
437 fv' = rotate cwz $ rotate ccwy $ Tetrahedron.fv (tetrahedron0 c)
440 tetrahedron22 :: Cube -> Tetrahedron
442 Tetrahedron fv' v0' v1' v2' v3'
445 v1' = center (left_face c)
446 v2' = Face.v2 (left_face c)
447 v3' = Face.v3 (left_face c)
448 fv' = rotate cwz $ rotate ccwy
450 $ Tetrahedron.fv (tetrahedron0 c)
453 tetrahedron23 :: Cube -> Tetrahedron
455 Tetrahedron fv' v0' v1' v2' v3'
458 v1' = center (left_face c)
459 v2' = Face.v3 (left_face c)
460 v3' = Face.v0 (left_face c)
461 fv' = rotate cwz $ rotate cwy
462 $ Tetrahedron.fv (tetrahedron0 c)
465 tetrahedra :: Cube -> [Tetrahedron]
493 -- | Takes a 'Cube', and returns all Tetrahedra belonging to it that
494 -- contain the given 'Point'.
495 find_containing_tetrahedra :: Cube -> Point -> [Tetrahedron]
496 find_containing_tetrahedra c p =
497 filter contains_our_point all_tetrahedra
499 contains_our_point = flip contains_point p
500 all_tetrahedra = tetrahedra c