1 {-# LANGUAGE BangPatterns #-}
2 -- | The FunctionValues module contains the 'FunctionValues' type and
3 -- the functions used to manipulate it.
4 module FunctionValues (
10 function_values_tests,
11 function_values_properties,
16 import Prelude hiding (LT)
17 import Test.HUnit (Assertion)
18 import Test.Framework (Test, testGroup)
19 import Test.Framework.Providers.HUnit (testCase)
20 import Test.Framework.Providers.QuickCheck2 (testProperty)
21 import Test.QuickCheck (Arbitrary(..), choose)
23 import Assertions (assertTrue)
24 import Cardinal ( Cardinal(..), cwx, cwy, cwz )
25 import Examples (trilinear)
26 import Values (Values3D, dims, idx)
28 -- | The FunctionValues type represents the value of our function f at
29 -- the 27 points surrounding (and including) the center of a
30 -- cube. Each value of f can be accessed by the name of its
33 FunctionValues { front :: !Double,
39 front_left :: !Double,
40 front_right :: !Double,
41 front_down :: !Double,
44 back_right :: !Double,
49 right_down :: !Double,
51 front_left_down :: !Double,
52 front_left_top :: !Double,
53 front_right_down :: !Double,
54 front_right_top :: !Double,
55 back_left_down :: !Double,
56 back_left_top :: !Double,
57 back_right_down :: !Double,
58 back_right_top :: !Double,
63 instance Arbitrary FunctionValues where
65 front' <- choose (min_double, max_double)
66 back' <- choose (min_double, max_double)
67 left' <- choose (min_double, max_double)
68 right' <- choose (min_double, max_double)
69 top' <- choose (min_double, max_double)
70 down' <- choose (min_double, max_double)
71 front_left' <- choose (min_double, max_double)
72 front_right' <- choose (min_double, max_double)
73 front_top' <- choose (min_double, max_double)
74 front_down' <- choose (min_double, max_double)
75 back_left' <- choose (min_double, max_double)
76 back_right' <- choose (min_double, max_double)
77 back_top' <- choose (min_double, max_double)
78 back_down' <- choose (min_double, max_double)
79 left_top' <- choose (min_double, max_double)
80 left_down' <- choose (min_double, max_double)
81 right_top' <- choose (min_double, max_double)
82 right_down' <- choose (min_double, max_double)
83 front_left_top' <- choose (min_double, max_double)
84 front_left_down' <- choose (min_double, max_double)
85 front_right_top' <- choose (min_double, max_double)
86 front_right_down' <- choose (min_double, max_double)
87 back_left_top' <- choose (min_double, max_double)
88 back_left_down' <- choose (min_double, max_double)
89 back_right_top' <- choose (min_double, max_double)
90 back_right_down' <- choose (min_double, max_double)
91 interior' <- choose (min_double, max_double)
93 return empty_values { front = front',
99 front_left = front_left',
100 front_right = front_right',
101 front_top = front_top',
102 front_down = front_down',
103 back_left = back_left',
104 back_right = back_right',
105 back_top = back_top',
106 back_down = back_down',
107 left_top = left_top',
108 left_down = left_down',
109 right_top = right_top',
110 right_down = right_down',
111 front_left_top = front_left_top',
112 front_left_down = front_left_down',
113 front_right_top = front_right_top',
114 front_right_down = front_right_down',
115 back_left_top = back_left_top',
116 back_left_down = back_left_down',
117 back_right_top = back_right_top',
118 back_right_down = back_right_down',
119 interior = interior' }
121 -- | We perform addition with the function values contained in a
122 -- FunctionValues object. If we choose random doubles near the machine
123 -- min/max, we risk overflowing or underflowing the 'Double'. This
124 -- places a reasonably safe limit on the maximum size of our generated
129 -- | See 'max_double'.
131 min_double = (-1) * max_double
134 -- | Return a 'FunctionValues' with a bunch of zeros for data points.
135 empty_values :: FunctionValues
137 FunctionValues 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
140 -- | The eval function is where the magic happens for the
141 -- FunctionValues type. Given a 'Cardinal' direction and a
142 -- 'FunctionValues' object, eval will return the value of the
143 -- function f in that 'Cardinal' direction. Note that 'Cardinal' can
144 -- be a composite type; eval is what performs the \"arithmetic\" on
145 -- 'Cardinal' directions.
146 eval :: FunctionValues -> Cardinal -> Double
153 eval f FL = front_left f
154 eval f FR = front_right f
155 eval f FD = front_down f
156 eval f FT = front_top f
157 eval f BL = back_left f
158 eval f BR = back_right f
159 eval f BD = back_down f
160 eval f BT = back_top f
161 eval f LD = left_down f
162 eval f LT = left_top f
163 eval f RD = right_down f
164 eval f RT = right_top f
165 eval f FLD = front_left_down f
166 eval f FLT = front_left_top f
167 eval f FRD = front_right_down f
168 eval f FRT = front_right_top f
169 eval f BLD = back_left_down f
170 eval f BLT = back_left_top f
171 eval f BRD = back_right_down f
172 eval f BRT = back_right_top f
173 eval f I = interior f
174 eval _ (Scalar x) = x
175 eval f (Sum x y) = (eval f x) + (eval f y)
176 eval f (Difference x y) = (eval f x) - (eval f y)
177 eval f (Product x y) = (eval f x) * (eval f y)
178 eval f (Quotient x y) = (eval f x) / (eval f y)
181 -- | Takes a three-dimensional list of 'Double' and a set of 3D
182 -- coordinates (i,j,k), and returns the value at (i,j,k) in the
183 -- supplied list. If there is no such value, we calculate one
184 -- according to Sorokina and Zeilfelder, remark 7.3, p. 99.
186 -- We specifically do not consider values more than one unit away
191 -- >>> value_at Examples.trilinear 0 0 0
194 -- >>> value_at Examples.trilinear (-1) 0 0
197 -- >>> value_at Examples.trilinear 0 0 4
200 -- >>> value_at Examples.trilinear 1 3 0
203 value_at :: Values3D -> Int -> Int -> Int -> Double
204 value_at v3d !i !j !k
205 -- Put the most common case first!
206 | (valid_i i) && (valid_j j) && (valid_k k) =
209 -- The next three are from the first line in (7.3). Analogous cases
210 -- have been added where the indices are one-too-big. These are the
211 -- "one index is bad" cases.
215 -- We're one-dimensional in our first coordinate, so just
216 -- return the data point that we do have. If we try to use
217 -- the formula from remark 7.3, we go into an infinite loop.
222 2*(value_at v3d 0 j k) - (value_at v3d 1 j k)
224 2*(value_at v3d (i-1) j k) - (value_at v3d (i-2) j k)
229 -- We're one-dimensional in our second coordinate, so just
230 -- return the data point that we do have. If we try to use
231 -- the formula from remark 7.3, we go into an infinite loop.
236 2*(value_at v3d i 0 k) - (value_at v3d i 1 k)
238 2*(value_at v3d i (j-1) k) - (value_at v3d i (j-2) k)
243 -- We're one-dimensional in our third coordinate, so just
244 -- return the data point that we do have. If we try to use
245 -- the formula from remark 7.3, we go into an infinite loop.
250 2*(value_at v3d i j 0) - (value_at v3d i j 1)
252 2*(value_at v3d i j (k-1)) - (value_at v3d i j (k-2))
254 (dim_i, dim_j, dim_k) = dims v3d
256 valid_i :: Int -> Bool
257 valid_i i' = (i' >= 0) && (i' < dim_i)
259 valid_j :: Int -> Bool
260 valid_j j' = (j' >= 0) && (j' < dim_j)
262 valid_k :: Int -> Bool
263 valid_k k' = (k' >= 0) && (k' < dim_k)
267 -- | Given a three-dimensional list of 'Double' and a set of 3D
268 -- coordinates (i,j,k), constructs and returns the 'FunctionValues'
269 -- object centered at (i,j,k)
270 make_values :: Values3D -> Int -> Int -> Int -> FunctionValues
271 make_values values !i !j !k =
272 empty_values { front = value_at values (i-1) j k,
273 back = value_at values (i+1) j k,
274 left = value_at values i (j-1) k,
275 right = value_at values i (j+1) k,
276 down = value_at values i j (k-1),
277 top = value_at values i j (k+1),
278 front_left = value_at values (i-1) (j-1) k,
279 front_right = value_at values (i-1) (j+1) k,
280 front_down =value_at values (i-1) j (k-1),
281 front_top = value_at values (i-1) j (k+1),
282 back_left = value_at values (i+1) (j-1) k,
283 back_right = value_at values (i+1) (j+1) k,
284 back_down = value_at values (i+1) j (k-1),
285 back_top = value_at values (i+1) j (k+1),
286 left_down = value_at values i (j-1) (k-1),
287 left_top = value_at values i (j-1) (k+1),
288 right_down = value_at values i (j+1) (k-1),
289 right_top = value_at values i (j+1) (k+1),
290 front_left_down = value_at values (i-1) (j-1) (k-1),
291 front_left_top = value_at values (i-1) (j-1) (k+1),
292 front_right_down = value_at values (i-1) (j+1) (k-1),
293 front_right_top = value_at values (i-1) (j+1) (k+1),
294 back_left_down = value_at values (i+1) (j-1) (k-1),
295 back_left_top = value_at values (i+1) (j-1) (k+1),
296 back_right_down = value_at values (i+1) (j+1) (k-1),
297 back_right_top = value_at values (i+1) (j+1) (k+1),
298 interior = value_at values i j k }
300 -- | Takes a 'FunctionValues' and a function that transforms one
301 -- 'Cardinal' to another (called a rotation). Then it applies the
302 -- rotation to each element of the 'FunctionValues' object, and
303 -- returns the result.
304 rotate :: (Cardinal -> Cardinal) -> FunctionValues -> FunctionValues
306 FunctionValues { front = eval fv (rotation F),
307 back = eval fv (rotation B),
308 left = eval fv (rotation L),
309 right = eval fv (rotation R),
310 down = eval fv (rotation D),
311 top = eval fv (rotation T),
312 front_left = eval fv (rotation FL),
313 front_right = eval fv (rotation FR),
314 front_down = eval fv (rotation FD),
315 front_top = eval fv (rotation FT),
316 back_left = eval fv (rotation BL),
317 back_right = eval fv (rotation BR),
318 back_down = eval fv (rotation BD),
319 back_top = eval fv (rotation BT),
320 left_down = eval fv (rotation LD),
321 left_top = eval fv (rotation LT),
322 right_down = eval fv (rotation RD),
323 right_top = eval fv (rotation RT),
324 front_left_down = eval fv (rotation FLD),
325 front_left_top = eval fv (rotation FLT),
326 front_right_down = eval fv (rotation FRD),
327 front_right_top = eval fv (rotation FRT),
328 back_left_down = eval fv (rotation BLD),
329 back_left_top = eval fv (rotation BLT),
330 back_right_down = eval fv (rotation BRD),
331 back_right_top = eval fv (rotation BRT),
332 interior = interior fv }
336 -- | Ensure that the trilinear values wind up where we think they
338 test_directions :: Assertion
340 assertTrue "all direction functions work" (and equalities)
342 fvs = make_values trilinear 1 1 1
343 equalities = [ interior fvs == 4,
351 front_right fvs == 1,
355 back_right fvs == 11,
362 front_left_down fvs == 1,
363 front_left_top fvs == 1,
364 front_right_down fvs == 1,
365 front_right_top fvs == 1,
366 back_left_down fvs == 3,
367 back_left_top fvs == 3,
368 back_right_down fvs == 7,
369 back_right_top fvs == 15]
372 function_values_tests :: Test.Framework.Test
373 function_values_tests =
374 testGroup "FunctionValues Tests"
375 [ testCase "test directions" test_directions ]
378 prop_x_rotation_doesnt_affect_front :: FunctionValues -> Bool
379 prop_x_rotation_doesnt_affect_front fv0 =
386 prop_x_rotation_doesnt_affect_back :: FunctionValues -> Bool
387 prop_x_rotation_doesnt_affect_back fv0 =
395 prop_y_rotation_doesnt_affect_left :: FunctionValues -> Bool
396 prop_y_rotation_doesnt_affect_left fv0 =
403 prop_y_rotation_doesnt_affect_right :: FunctionValues -> Bool
404 prop_y_rotation_doesnt_affect_right fv0 =
412 prop_z_rotation_doesnt_affect_down :: FunctionValues -> Bool
413 prop_z_rotation_doesnt_affect_down fv0 =
421 prop_z_rotation_doesnt_affect_top :: FunctionValues -> Bool
422 prop_z_rotation_doesnt_affect_top fv0 =
430 function_values_properties :: Test.Framework.Test
431 function_values_properties =
432 let tp = testProperty
434 testGroup "FunctionValues Properties" [
435 tp "x rotation doesn't affect front" prop_x_rotation_doesnt_affect_front,
436 tp "x rotation doesn't affect back" prop_x_rotation_doesnt_affect_back,
437 tp "y rotation doesn't affect left" prop_y_rotation_doesnt_affect_left,
438 tp "y rotation doesn't affect right" prop_y_rotation_doesnt_affect_right,
439 tp "z rotation doesn't affect top" prop_z_rotation_doesnt_affect_top,
440 tp "z rotation doesn't affect down" prop_z_rotation_doesnt_affect_down ]