src/Cube.hs

   1 module Cube
   2 where
   3
   4 import Face
   5 import FunctionValues
   6 --import Grid
   7 import Point
   8 import ThreeDimensional
   9
  10 data Cube = Cube { h :: Double,
  11                    i :: Int,
  12                    j :: Int,
  13                    k :: Int,
  14                    fv :: FunctionValues }
  15             deriving (Eq)
  16
  17
  18 instance Show Cube where
  19     show c =
  20         "Cube_" ++ (show (i c)) ++ "," ++ (show (j c)) ++ "," ++ (show (k c)) ++
  21         " (Center: " ++ (show (center c)) ++ ")" ++
  22         " (xmin: " ++ (show (xmin c)) ++ ")" ++
  23         " (xmax: " ++ (show (xmax c)) ++ ")" ++
  24         " (ymin: " ++ (show (ymin c)) ++ ")" ++
  25         " (ymax: " ++ (show (ymax c)) ++ ")" ++
  26         " (zmin: " ++ (show (zmin c)) ++ ")" ++
  27         " (zmax: " ++ (show (zmax c)) ++ ")"
  28
  29 empty_cube :: Cube
  30 empty_cube = Cube 0 0 0 0 empty_values
  31
  32
  33 -- | The left-side boundary of the cube. See Sorokina and Zeilfelder,
  34 --   p. 76.
  35 xmin :: Cube -> Double
  36 xmin c = (2*i' - 1)*delta / 2
  37     where
  38       i' = fromIntegral (i c) :: Double
  39       delta = h c
  40
  41 -- | The right-side boundary of the cube. See Sorokina and Zeilfelder,
  42 --   p. 76.
  43 xmax :: Cube -> Double
  44 xmax c = (2*i' + 1)*delta / 2
  45     where
  46       i' = fromIntegral (i c) :: Double
  47       delta = h c
  48
  49 -- | The front boundary of the cube. See Sorokina and Zeilfelder,
  50 --   p. 76.
  51 ymin :: Cube -> Double
  52 ymin c = (2*j' - 1)*delta / 2
  53     where
  54       j' = fromIntegral (j c) :: Double
  55       delta = h c
  56
  57 -- | The back boundary of the cube. See Sorokina and Zeilfelder,
  58 --   p. 76.
  59 ymax :: Cube -> Double
  60 ymax c = (2*j' + 1)*delta / 2
  61     where
  62       j' = fromIntegral (j c) :: Double
  63       delta = h c
  64
  65 -- | The bottom boundary of the cube. See Sorokina and Zeilfelder,
  66 --   p. 76.
  67 zmin :: Cube -> Double
  68 zmin c = (2*k' - 1)*delta / 2
  69     where
  70       k' = fromIntegral (k c) :: Double
  71       delta = h c
  72
  73 -- | The top boundary of the cube. See Sorokina and Zeilfelder,
  74 --   p. 76.
  75 zmax :: Cube -> Double
  76 zmax c = (2*k' + 1)*delta / 2
  77     where
  78       k' = fromIntegral (k c) :: Double
  79       delta = h c
  80
  81 instance ThreeDimensional Cube where
  82     -- | The center of Cube_ijk coincides with v_ijk at
  83     --   (ih, jh, kh). See Sorokina and Zeilfelder, p. 76.
  84     center c = (x, y, z)
  85            where
  86              delta = h c
  87              i' = fromIntegral (i c) :: Double
  88              j' = fromIntegral (j c) :: Double
  89              k' = fromIntegral (k c) :: Double
  90              x = delta * i'
  91              y = delta * j'
  92              z = delta * k'
  93
  94     contains_point c p
  95         | (x_coord p) < (xmin c) = False
  96         | (x_coord p) > (xmax c) = False
  97         | (y_coord p) < (ymin c) = False
  98         | (y_coord p) > (ymax c) = False
  99         | (z_coord p) < (zmin c) = False
 100         | (z_coord p) > (zmax c) = False
 101         | otherwise = True
 102
 103
 104 -- instance Num Cube where
 105 --     (Cube g1 i1 j1 k1 d1) + (Cube _ i2 j2 k2 d2) =
 106 --         Cube g1 (i1 + i2) (j1 + j2) (k1 + k2) (d1 + d2)
 107
 108 --     (Cube g1 i1 j1 k1 d1) - (Cube _ i2 j2 k2 d2) =
 109 --         Cube g1 (i1 - i2) (j1 - j2) (k1 - k2) (d1 - d2)
 110
 111 --     (Cube g1 i1 j1 k1 d1) * (Cube _ i2 j2 k2 d2) =
 112 --         Cube g1 (i1 * i2) (j1 * j2) (k1 * k2) (d1 * d2)
 113
 114 --     abs (Cube g1 i1 j1 k1 d1) =
 115 --         Cube g1 (abs i1) (abs j1) (abs k1) (abs d1)
 116
 117 --     signum (Cube g1 i1 j1 k1 d1) =
 118 --         Cube g1 (signum i1) (signum j1) (signum k1) (signum d1)
 119
 120 --     fromInteger x = empty_cube { datum = (fromIntegral x) }
 121
 122 -- instance Fractional Cube where
 123 --     (Cube g1 i1 j1 k1 d1) / (Cube _ _ _ _ d2) =
 124 --         Cube g1 i1 j1 k1 (d1 / d2)
 125
 126 --     recip (Cube g1 i1 j1 k1 d1) =
 127 --         Cube g1 i1 j1 k1 (recip d1)
 128
 129 --     fromRational q = empty_cube { datum = fromRational q }
 130
 131
 132
 133 -- | Return the cube corresponding to the grid point i,j,k. The list
 134 --   of cubes is stored as [z][y][x] but we'll be requesting it by
 135 --   [x][y][z] so we flip the indices in the last line.
 136 -- cube_at :: Grid -> Int -> Int -> Int -> Cube
 137 -- cube_at g i' j' k'
 138 --         | i' >= length (function_values g) = Cube g i' j' k' 0
 139 --         | i' < 0                          = Cube g i' j' k' 0
 140 --         | j' >= length ((function_values g) !! i') = Cube g i' j' k' 0
 141 --         | j' < 0                                  = Cube g i' j' k' 0
 142 --         | k' >= length (((function_values g) !! i') !! j') = Cube g i' j' k' 0
 143 --         | k' < 0                                          = Cube g i' j' k' 0
 144 --         | otherwise =
 145 --             (((cubes g) !! k') !! j') !! i'
 146
 147
 148
 149
 150
 151
 152 -- Face stuff.
 153
 154 -- | The top (in the direction of z) face of the cube.
 155 top_face :: Cube -> Face
 156 top_face c = Face v0' v1' v2' v3'
 157     where
 158       delta = (1/2)*(h c)
 159       v0' = (center c) + (-delta, delta, delta)
 160       v1' = (center c) + (delta, delta, delta)
 161       v2' = (center c) + (delta, -delta, delta)
 162       v3' = (center c) + (-delta, -delta, delta)
 163
 164
 165
 166 -- | The back (in the direction of x) face of the cube.
 167 back_face :: Cube -> Face
 168 back_face c = Face v0' v1' v2' v3'
 169     where
 170       delta = (1/2)*(h c)
 171       v0' = (center c) + (delta, delta, delta)
 172       v1' = (center c) + (delta, delta, -delta)
 173       v2' = (center c) + (delta, -delta, -delta)
 174       v3' = (center c) + (delta, -delta, delta)
 175
 176
 177 -- The bottom face (in the direction of -z) of the cube.
 178 down_face :: Cube -> Face
 179 down_face c = Face v0' v1' v2' v3'
 180     where
 181       delta = (1/2)*(h c)
 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)
 186
 187
 188
 189 -- | The front (in the direction of -x) face of the cube.
 190 front_face :: Cube -> Face
 191 front_face c = Face v0' v1' v2' v3'
 192     where
 193       delta = (1/2)*(h c)
 194       v0' = (center c) + (-delta, delta, -delta)
 195       v1' = (center c) + (-delta, delta, delta)
 196       v2' = (center c) + (-delta, -delta, delta)
 197       v3' = (center c) + (-delta, -delta, -delta)
 198
 199
 200 -- | The left (in the direction of -y) face of the cube.
 201 left_face :: Cube -> Face
 202 left_face c = Face v0' v1' v2' v3'
 203     where
 204       delta = (1/2)*(h c)
 205       v0' = (center c) + (-delta, -delta, delta)
 206       v1' = (center c) + (delta, -delta, delta)
 207       v2' = (center c) + (delta, -delta, -delta)
 208       v3' = (center c) + (-delta, -delta, -delta)
 209
 210
 211 -- | The right (in the direction of y) face of the cube.
 212 right_face :: Cube -> Face
 213 right_face c = Face v0' v1' v2' v3'
 214     where
 215       delta = (1/2)*(h c)
 216       v0' = (center c) + (-delta, delta, -delta)
 217       v1' = (center c) + (delta, delta, -delta)
 218       v2' = (center c) + (delta, delta, delta)
 219       v3' = (center c) + (-delta, delta, delta)
 220