import Cardinal
import Comparisons ((~=), (~~=))
import qualified Face (Face(Face, v0, v1, v2, v3))
-import FunctionValues
+import FunctionValues (FunctionValues, eval, rotate)
import Misc (all_equal, disjoint)
import Point
-import Tetrahedron (
- Tetrahedron(..),
- c,
- b0,
- b1,
- b2,
- b3,
- volume
- )
+import Tetrahedron (Tetrahedron(..), c, volume)
import ThreeDimensional
data Cube = Cube { h :: Double,
fv' <- arbitrary :: Gen FunctionValues
(Positive tet_vol) <- arbitrary :: Gen (Positive Double)
return (Cube h' i' j' k' fv' tet_vol)
- where
- coordmin = -268435456 -- -(2^29 / 2)
- coordmax = 268435456 -- +(2^29 / 2)
+ where
+ -- The idea here is that, when cubed in the volume formula,
+ -- these numbers don't overflow 64 bits. This number is not
+ -- magic in any other sense than that it does not cause test
+ -- failures, while 2^23 does.
+ coordmax = 4194304 -- 2^22
+ coordmin = -coordmax
instance Show Cube where
" ymin: " ++ (show (ymin cube)) ++ "\n" ++
" ymax: " ++ (show (ymax cube)) ++ "\n" ++
" zmin: " ++ (show (zmin cube)) ++ "\n" ++
- " zmax: " ++ (show (zmax cube)) ++ "\n" ++
- " fv: " ++ (show (Cube.fv cube)) ++ "\n"
+ " zmax: " ++ (show (zmax cube)) ++ "\n"
where
subscript =
(show (i cube)) ++ "," ++ (show (j cube)) ++ "," ++ (show (k cube))
-- | The left-side boundary of the cube. See Sorokina and Zeilfelder,
-- p. 76.
xmin :: Cube -> Double
-xmin cube = (2*i' - 1)*delta / 2
+xmin cube = (i' - 1/2)*delta
where
i' = fromIntegral (i cube) :: Double
delta = h cube
-- | The right-side boundary of the cube. See Sorokina and Zeilfelder,
-- p. 76.
xmax :: Cube -> Double
-xmax cube = (2*i' + 1)*delta / 2
+xmax cube = (i' + 1/2)*delta
where
i' = fromIntegral (i cube) :: Double
delta = h cube
-- | The front boundary of the cube. See Sorokina and Zeilfelder,
-- p. 76.
ymin :: Cube -> Double
-ymin cube = (2*j' - 1)*delta / 2
+ymin cube = (j' - 1/2)*delta
where
j' = fromIntegral (j cube) :: Double
delta = h cube
-- | The back boundary of the cube. See Sorokina and Zeilfelder,
-- p. 76.
ymax :: Cube -> Double
-ymax cube = (2*j' + 1)*delta / 2
+ymax cube = (j' + 1/2)*delta
where
j' = fromIntegral (j cube) :: Double
delta = h cube
-- | The bottom boundary of the cube. See Sorokina and Zeilfelder,
-- p. 76.
zmin :: Cube -> Double
-zmin cube = (2*k' - 1)*delta / 2
+zmin cube = (k' - 1/2)*delta
where
k' = fromIntegral (k cube) :: Double
delta = h cube
-- | The top boundary of the cube. See Sorokina and Zeilfelder,
-- p. 76.
zmax :: Cube -> Double
-zmax cube = (2*k' + 1)*delta / 2
+zmax cube = (k' + 1/2)*delta
where
k' = fromIntegral (k cube) :: Double
delta = h cube
-- | Since the grid size is necessarily positive, all tetrahedra
--- (which comprise cubes of positive volume) must have positive volume
--- as well.
+-- (which comprise cubes of positive volume) must have positive
+-- volume as well.
prop_all_volumes_positive :: Cube -> Bool
prop_all_volumes_positive cube =
- null nonpositive_volumes
+ all (>= 0) volumes
where
ts = tetrahedra cube
volumes = map volume ts
- nonpositive_volumes = filter (<= 0) volumes
+
-- | In fact, since all of the tetrahedra are identical, we should
-- already know their volumes. There's 24 tetrahedra to a cube, so
-- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Note that the
--- third and fourth indices of c-t1 have been switched. This is
+-- third and fourth indices of c-t3 have been switched. This is
-- because we store the triangles oriented such that their volume is
-- positive. If T and T-tilde share \<v0,v1,v2\> and v3,v3-tilde point
-- in opposite directions, one of them has to have negative volume!
t4 = tetrahedron cube 4
t5 = tetrahedron cube 5
--- -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
--- -- 'prop_c0120_identity1' with tetrahedrons 5 and 6.
+-- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
+-- 'prop_c0120_identity1' with tetrahedrons 5 and 6.
prop_c0120_identity6 :: Cube -> Bool
prop_c0120_identity6 cube =
c t6 0 1 2 0 ~= (c t6 0 0 2 1 + c t5 0 0 1 2) / 2
t6 = tetrahedron cube 6
--- -- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
--- -- 'prop_c0120_identity1' with tetrahedrons 6 and 7.
+-- | Given in Sorokina and Zeilfelder, p. 79, (2.6). Repeats
+-- 'prop_c0120_identity1' with tetrahedrons 6 and 7.
prop_c0120_identity7 :: Cube -> Bool
prop_c0120_identity7 cube =
c t7 0 1 2 0 ~= (c t7 0 0 2 1 + c t6 0 0 1 2) / 2
t6 = tetrahedron cube 6
-
--- | Given in Sorokina and Zeilfelder, p. 78.
-prop_cijk1_identity :: Cube -> Bool
-prop_cijk1_identity cube =
- and [ c t0 i j k 1 ~=
- (c t1 (i+1) j k 0) * ((b0 t0) (v3 t1)) +
- (c t1 i (j+1) k 0) * ((b1 t0) (v3 t1)) +
- (c t1 i j (k+1) 0) * ((b2 t0) (v3 t1)) +
- (c t1 i j k 1) * ((b3 t0) (v3 t1)) | i <- [0..2],
- j <- [0..2],
- k <- [0..2],
- i + j + k == 2]
- where
- t0 = tetrahedron cube 0
- t1 = tetrahedron cube 1
-
-
-- | The function values at the interior should be the same for all
-- tetrahedra.
prop_interior_values_all_identical :: Cube -> Bool
t20 = tetrahedron cube 20
-
-
-
-p78_25_properties :: Test.Framework.Test
-p78_25_properties =
- testGroup "p. 78, Section (2.5) Properties" [
- testProperty "c_ijk1 identity" prop_cijk1_identity ]
-
p79_26_properties :: Test.Framework.Test
p79_26_properties =
testGroup "p. 79, Section (2.6) Properties" [
cube_properties :: Test.Framework.Test
cube_properties =
testGroup "Cube Properties" [
- p78_25_properties,
p79_26_properties,
p79_27_properties,
p79_28_properties,
projects / spline3.git / blobdiff