1 {-# LANGUAGE BangPatterns #-}
3 -- | The FunctionValues module contains the 'FunctionValues' type and
4 -- the functions used to manipulate it.
6 module FunctionValues (
12 function_values_tests,
13 function_values_properties,
17 import Prelude hiding ( LT )
18 import Test.Tasty ( TestTree, testGroup )
19 import Test.Tasty.HUnit ( Assertion, testCase )
20 import Test.Tasty.QuickCheck ( Arbitrary(..), choose, testProperty )
22 import Assertions ( assertTrue )
23 import Cardinal ( Cardinal(..), cwx, cwy, cwz )
24 import Examples ( trilinear )
25 import Values ( Values3D, dims, idx )
27 -- | The FunctionValues type represents the value of our function f at
28 -- the 27 points surrounding (and including) the center of a
29 -- 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 :: TestTree
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 :: TestTree
431 function_values_properties =
432 testGroup "FunctionValues properties" [
434 "x rotation doesn't affect front"
435 prop_x_rotation_doesnt_affect_front,
437 "x rotation doesn't affect back"
438 prop_x_rotation_doesnt_affect_back,
440 "y rotation doesn't affect left"
441 prop_y_rotation_doesnt_affect_left,
443 "y rotation doesn't affect right"
444 prop_y_rotation_doesnt_affect_right,
446 "z rotation doesn't affect top"
447 prop_z_rotation_doesnt_affect_top,
449 "z rotation doesn't affect down"
450 prop_z_rotation_doesnt_affect_down ]