grid.jl (19258B)
1 #!/usr/bin/env julia 2 using Test 3 import Granular 4 5 # Check the grid interpolation and sorting functions 6 verbose = false 7 8 if Granular.hasNetCDF 9 ocean = Granular.readOceanNetCDF("Baltic/00010101.ocean_month.nc", 10 "Baltic/ocean_hgrid.nc") 11 12 @info "Testing coordinate retrieval functions" 13 sw, se, ne, nw = Granular.getCellCornerCoordinates(ocean.xq, ocean.yq, 1, 1) 14 @test sw ≈ [6., 53.] 15 @test se ≈ [7., 53.] 16 @test ne ≈ [7., 54.] 17 @test nw ≈ [6., 54.] 18 @test Granular.getCellCenterCoordinates(ocean.xh, ocean.yh, 1, 1) ≈ [6.5, 53.5] 19 20 @info "Testing area-determination methods" 21 @test Granular.areaOfTriangle([0., 0.], [1., 0.], [0., 1.]) ≈ .5 22 @test Granular.areaOfTriangle([1., 0.], [0., 1.], [0., 0.]) ≈ .5 23 @test Granular.areaOfQuadrilateral([1., 0.], [0., 1.], [0., 0.], [1., 1.]) ≈ 1. 24 25 @info "Testing area-based cell content determination" 26 @test Granular.isPointInCell(ocean, 1, 1, [6.5, 53.5], sw, se, ne, nw) == true 27 @test Granular.isPointInCell(ocean, 1, 1, [6.5, 53.5]) == true 28 @test Granular.getNonDimensionalCellCoordinates(ocean, 1, 1, [6.5, 53.5]) ≈ 29 [.5, .5] 30 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.5], sw, se, ne, nw) == true 31 @test Granular.getNonDimensionalCellCoordinates(ocean, 1, 1, [6.1, 53.5]) ≈ 32 [.1, .5] 33 @test Granular.isPointInCell(ocean, 1, 1, [6.0, 53.5], sw, se, ne, nw) == true 34 @test Granular.getNonDimensionalCellCoordinates(ocean, 1, 1, [6.0, 53.5]) ≈ 35 [.0, .5] 36 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.7], sw, se, ne, nw) == true 37 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.9], sw, se, ne, nw) == true 38 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.99999], sw, se, ne, nw) == true 39 @test Granular.getNonDimensionalCellCoordinates(ocean, 1, 1, [6.1, 53.99999]) ≈ 40 [.1, .99999] 41 @test Granular.isPointInCell(ocean, 1, 1, [7.5, 53.5], sw, se, ne, nw) == false 42 @test Granular.isPointInCell(ocean, 1, 1, [0.0, 53.5], sw, se, ne, nw) == false 43 x_tilde, _ = Granular.getNonDimensionalCellCoordinates(ocean, 1, 1, [0., 53.5]) 44 @test x_tilde < 0. 45 46 @info "Testing conformal mapping methods" 47 @test Granular.conformalQuadrilateralCoordinates([0., 0.], 48 [5., 0.], 49 [5., 3.], 50 [0., 3.], 51 [2.5, 1.5]) ≈ [0.5, 0.5] 52 @test Granular.conformalQuadrilateralCoordinates([0., 0.], 53 [5., 0.], 54 [5., 3.], 55 [0., 3.], 56 [7.5, 1.5]) ≈ [1.5, 0.5] 57 @test Granular.conformalQuadrilateralCoordinates([0., 0.], 58 [5., 0.], 59 [5., 3.], 60 [0., 3.], 61 [7.5,-1.5]) ≈ [1.5,-0.5] 62 @test_throws ErrorException Granular.conformalQuadrilateralCoordinates([0., 0.], 63 [5., 3.], 64 [0., 3.], 65 [5., 0.], 66 [7.5,-1.5]) 67 68 @info "Checking cell content using conformal mapping methods" 69 @test Granular.isPointInCell(ocean, 1, 1, [6.4, 53.4], sw, se, ne, nw, 70 method="Conformal") == true 71 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.5], sw, se, ne, nw, 72 method="Conformal") == true 73 @test Granular.isPointInCell(ocean, 1, 1, [6.0, 53.5], sw, se, ne, nw, 74 method="Conformal") == true 75 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.7], sw, se, ne, nw, 76 method="Conformal") == true 77 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.9], sw, se, ne, nw, 78 method="Conformal") == true 79 @test Granular.isPointInCell(ocean, 1, 1, [6.1, 53.99999], sw, se, ne, nw, 80 method="Conformal") == true 81 @test Granular.isPointInCell(ocean, 1, 1, [7.5, 53.5], sw, se, ne, nw, 82 method="Conformal") == false 83 @test Granular.isPointInCell(ocean, 1, 1, [0.0, 53.5], sw, se, ne, nw, 84 method="Conformal") == false 85 86 @info "Testing bilinear interpolation scheme on conformal mapping" 87 ocean.u[1, 1, 1, 1] = 1.0 88 ocean.u[2, 1, 1, 1] = 1.0 89 ocean.u[2, 2, 1, 1] = 0.0 90 ocean.u[1, 2, 1, 1] = 0.0 91 val = [NaN, NaN] 92 Granular.bilinearInterpolation!(val, ocean.u[:,:,1,1], ocean.u[:,:,1,1], 93 .5, .5, 1, 1) 94 @time Granular.bilinearInterpolation!(val, ocean.u[:,:,1,1], ocean.u[:,:,1,1], .5, 95 .5, 1, 1) 96 @test val[1] ≈ .5 97 @test val[2] ≈ .5 98 Granular.bilinearInterpolation!(val, ocean.u[:,:,1,1], ocean.u[:,:,1,1], 1., 1., 99 1, 1) 100 @test val[1] ≈ .0 101 @test val[2] ≈ .0 102 Granular.bilinearInterpolation!(val, ocean.u[:,:,1,1], ocean.u[:,:,1,1], 0., 0., 103 1, 1) 104 @test val[1] ≈ 1. 105 @test val[2] ≈ 1. 106 Granular.bilinearInterpolation!(val, ocean.u[:,:,1,1], ocean.u[:,:,1,1], .25, .25, 107 1, 1) 108 @test val[1] ≈ .75 109 @test val[2] ≈ .75 110 Granular.bilinearInterpolation!(val, ocean.u[:,:,1,1], ocean.u[:,:,1,1], .75, .75, 111 1, 1) 112 @test val[1] ≈ .25 113 @test val[2] ≈ .25 114 115 @info "Testing cell binning - Area-based approach" 116 @test Granular.findCellContainingPoint(ocean, [6.2,53.4], method="Area") == (1, 1) 117 @test Granular.findCellContainingPoint(ocean, [7.2,53.4], method="Area") == (2, 1) 118 @test Granular.findCellContainingPoint(ocean, [0.2,53.4], method="Area") == (0, 0) 119 120 @info "Testing cell binning - Conformal mapping" 121 @test Granular.findCellContainingPoint(ocean, [6.2,53.4], method="Conformal") == 122 (1, 1) 123 @test Granular.findCellContainingPoint(ocean, [7.2,53.4], method="Conformal") == 124 (2, 1) 125 @test Granular.findCellContainingPoint(ocean, [0.2, 53.4], method="Conformal") == 126 (0, 0) 127 128 sim = Granular.createSimulation() 129 sim.ocean = Granular.readOceanNetCDF("Baltic/00010101.ocean_month.nc", 130 "Baltic/ocean_hgrid.nc") 131 Granular.addGrainCylindrical!(sim, [6.5, 53.5], 10., 1., verbose=verbose) 132 Granular.addGrainCylindrical!(sim, [6.6, 53.5], 10., 1., verbose=verbose) 133 Granular.addGrainCylindrical!(sim, [7.5, 53.5], 10., 1., verbose=verbose) 134 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 135 @test sim.grains[1].ocean_grid_pos == [1, 1] 136 @test sim.grains[2].ocean_grid_pos == [1, 1] 137 @test sim.grains[3].ocean_grid_pos == [2, 1] 138 @test sim.ocean.grain_list[1, 1] == [1, 2] 139 @test sim.ocean.grain_list[2, 1] == [3] 140 end 141 142 @info "Testing ocean drag" 143 sim = Granular.createSimulation() 144 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 145 sim.ocean.u[:,:,1,1] .= 5. 146 Granular.addGrainCylindrical!(sim, [2.5, 3.5], 1., 1., verbose=verbose) 147 Granular.addGrainCylindrical!(sim, [2.6, 2.5], 1., 1., verbose=verbose) 148 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 149 if !Granular.hasNetCDF 150 ocean = sim.ocean 151 end 152 sim.time = ocean.time[1] 153 Granular.addOceanDrag!(sim) 154 @test sim.grains[1].force[1] > 0. 155 @test sim.grains[1].force[2] ≈ 0. 156 @test sim.grains[2].force[1] > 0. 157 @test sim.grains[2].force[2] ≈ 0. 158 sim.ocean.u[:,:,1,1] .= -5. 159 sim.ocean.v[:,:,1,1] .= 5. 160 Granular.addGrainCylindrical!(sim, [2.5, 3.5], 1., 1., verbose=verbose) 161 Granular.addGrainCylindrical!(sim, [2.6, 2.5], 1., 1., verbose=verbose) 162 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 163 sim.time = ocean.time[1] 164 Granular.addOceanDrag!(sim) 165 @test sim.grains[1].force[1] < 0. 166 @test sim.grains[1].force[2] > 0. 167 @test sim.grains[2].force[1] < 0. 168 @test sim.grains[2].force[2] > 0. 169 170 @info "Testing curl function" 171 ocean.u[1, 1, 1, 1] = 1.0 172 ocean.u[2, 1, 1, 1] = 1.0 173 ocean.u[2, 2, 1, 1] = 0.0 174 ocean.u[1, 2, 1, 1] = 0.0 175 ocean.v[:, :, 1, 1] .= 0.0 176 sw = zeros(2) 177 se = zeros(2) 178 ne = zeros(2) 179 nw = zeros(2) 180 @test Granular.curl(ocean, .5, .5, 1, 1, 1, 1, sw, se, ne, nw) > 0. 181 182 ocean.u[1, 1, 1, 1] = 0.0 183 ocean.u[2, 1, 1, 1] = 0.0 184 ocean.u[2, 2, 1, 1] = 1.0 185 ocean.u[1, 2, 1, 1] = 1.0 186 ocean.v[:, :, 1, 1] .= 0.0 187 @test Granular.curl(ocean, .5, .5, 1, 1, 1, 1, sw, se, ne, nw) < 0. 188 189 @info "Testing atmosphere drag" 190 sim = Granular.createSimulation() 191 sim.atmosphere = Granular.createRegularAtmosphereGrid([4, 4, 2], [4., 4., 2.]) 192 atmosphere = Granular.createRegularAtmosphereGrid([4, 4, 2], [4., 4., 2.]) 193 sim.atmosphere.u[:,:,1,1] .= 5. 194 Granular.addGrainCylindrical!(sim, [2.5, 3.5], 1., 1., verbose=verbose) 195 Granular.addGrainCylindrical!(sim, [2.6, 2.5], 1., 1., verbose=verbose) 196 Granular.sortGrainsInGrid!(sim, sim.atmosphere, verbose=verbose) 197 sim.time = ocean.time[1] 198 Granular.addAtmosphereDrag!(sim) 199 @test sim.grains[1].force[1] > 0. 200 @test sim.grains[1].force[2] ≈ 0. 201 @test sim.grains[2].force[1] > 0. 202 @test sim.grains[2].force[2] ≈ 0. 203 sim.atmosphere.u[:,:,1,1] .= -5. 204 sim.atmosphere.v[:,:,1,1] .= 5. 205 Granular.addGrainCylindrical!(sim, [2.5, 3.5], 1., 1., verbose=verbose) 206 Granular.addGrainCylindrical!(sim, [2.6, 2.5], 1., 1., verbose=verbose) 207 Granular.sortGrainsInGrid!(sim, sim.atmosphere, verbose=verbose) 208 sim.time = ocean.time[1] 209 Granular.addAtmosphereDrag!(sim) 210 @test sim.grains[1].force[1] < 0. 211 @test sim.grains[1].force[2] > 0. 212 @test sim.grains[2].force[1] < 0. 213 @test sim.grains[2].force[2] > 0. 214 215 @info "Testing curl function" 216 atmosphere.u[1, 1, 1, 1] = 1.0 217 atmosphere.u[2, 1, 1, 1] = 1.0 218 atmosphere.u[2, 2, 1, 1] = 0.0 219 atmosphere.u[1, 2, 1, 1] = 0.0 220 atmosphere.v[:, :, 1, 1] .= 0.0 221 @test Granular.curl(atmosphere, .5, .5, 1, 1, 1, 1) > 0. 222 @test Granular.curl(atmosphere, .5, .5, 1, 1, 1, 1, sw, se, ne, nw) > 0. 223 224 atmosphere.u[1, 1, 1, 1] = 0.0 225 atmosphere.u[2, 1, 1, 1] = 0.0 226 atmosphere.u[2, 2, 1, 1] = 1.0 227 atmosphere.u[1, 2, 1, 1] = 1.0 228 atmosphere.v[:, :, 1, 1] .= 0.0 229 @test Granular.curl(atmosphere, .5, .5, 1, 1, 1, 1) < 0. 230 @test Granular.curl(atmosphere, .5, .5, 1, 1, 1, 1, sw, se, ne, nw) < 0. 231 232 233 @info "Testing findEmptyPositionInGridCell" 234 @info "# Insert into empty cell" 235 sim = Granular.createSimulation() 236 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 237 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 238 pos = Granular.findEmptyPositionInGridCell(sim, sim.ocean, 1, 1, 0.5, 239 verbose=true) 240 @test pos != false 241 @test Granular.isPointInCell(sim.ocean, 1, 1, pos) == true 242 243 @info "# Insert into cell with one other ice floe" 244 sim = Granular.createSimulation() 245 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 246 Granular.addGrainCylindrical!(sim, [.25, .25], .25, 1., verbose=verbose) 247 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 248 pos = Granular.findEmptyPositionInGridCell(sim, sim.ocean, 1, 1, .25, 249 verbose=true) 250 @test pos != false 251 @test Granular.isPointInCell(sim.ocean, 1, 1, pos) == true 252 253 @info "# Insert into cell with two other grains" 254 sim = Granular.createSimulation() 255 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 256 Granular.addGrainCylindrical!(sim, [.25, .25], .25, 1., verbose=verbose) 257 Granular.addGrainCylindrical!(sim, [.75, .75], .25, 1., verbose=verbose) 258 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 259 pos = Granular.findEmptyPositionInGridCell(sim, sim.ocean, 1, 1, .25, n_iter=30, 260 verbose=true) 261 @test pos != false 262 @test Granular.isPointInCell(sim.ocean, 1, 1, pos) == true 263 264 @info "# Insert into full cell" 265 sim = Granular.createSimulation() 266 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 267 Granular.addGrainCylindrical!(sim, [.25, .25], 1., 1., verbose=verbose) 268 Granular.addGrainCylindrical!(sim, [.75, .25], 1., 1., verbose=verbose) 269 Granular.addGrainCylindrical!(sim, [.25, .75], 1., 1., verbose=verbose) 270 Granular.addGrainCylindrical!(sim, [.75, .75], 1., 1., verbose=verbose) 271 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 272 pos = Granular.findEmptyPositionInGridCell(sim, sim.ocean, 1, 1, 0.5, 273 verbose=false) 274 @test pos == false 275 276 @info "# Insert into empty cell" 277 sim = Granular.createSimulation() 278 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 279 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 280 pos = Granular.findEmptyPositionInGridCell(sim, sim.ocean, 2, 2, 0.5, 281 verbose=true) 282 @test pos != false 283 @test Granular.isPointInCell(sim.ocean, 2, 2, pos) == true 284 285 @info "# Insert into full cell" 286 sim = Granular.createSimulation() 287 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 288 Granular.addGrainCylindrical!(sim, [1.5, 1.5], 1., 1., verbose=verbose) 289 Granular.addGrainCylindrical!(sim, [1.75, 1.5], 1., 1., verbose=verbose) 290 Granular.addGrainCylindrical!(sim, [1.5, 1.75], 1., 1., verbose=verbose) 291 Granular.addGrainCylindrical!(sim, [1.75, 1.75], 1., 1., verbose=verbose) 292 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 293 pos = Granular.findEmptyPositionInGridCell(sim, sim.ocean, 2, 2, 0.5, 294 verbose=false) 295 @test pos == false 296 297 @info "Test default sorting with ocean/atmosphere grids" 298 sim = Granular.createSimulation() 299 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 300 sim.atmosphere = Granular.createRegularAtmosphereGrid([4, 4, 2], [4., 4.000001, 2.]) 301 Granular.addGrainCylindrical!(sim, [0.5, 0.5], .1, 1., verbose=verbose) 302 Granular.addGrainCylindrical!(sim, [0.7, 0.7], .1, 1., verbose=verbose) 303 Granular.addGrainCylindrical!(sim, [2.6, 2.5], .1, 1., verbose=verbose) 304 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 305 Granular.setTimeStep!(sim) 306 Granular.setTotalTime!(sim, 1.0) 307 Granular.run!(sim, single_step=true, verbose=verbose) 308 @test sim.atmosphere.collocated_with_ocean_grid == false 309 @test sim.grains[1].ocean_grid_pos == [1, 1] 310 @test sim.grains[2].ocean_grid_pos == [1, 1] 311 @test sim.grains[3].ocean_grid_pos == [3, 3] 312 @test sim.ocean.grain_list[1, 1] == [1, 2] 313 @test sim.ocean.grain_list[2, 2] == [] 314 @test sim.ocean.grain_list[3, 3] == [3] 315 @test sim.grains[1].atmosphere_grid_pos == [1, 1] 316 @test sim.grains[2].atmosphere_grid_pos == [1, 1] 317 @test sim.grains[3].atmosphere_grid_pos == [3, 3] 318 @test sim.atmosphere.grain_list[1, 1] == [1, 2] 319 @test sim.atmosphere.grain_list[2, 2] == [] 320 @test sim.atmosphere.grain_list[3, 3] == [3] 321 322 @info "Test optimization when ocean/atmosphere grids are collocated" 323 sim = Granular.createSimulation() 324 sim.ocean = Granular.createRegularOceanGrid([4, 4, 2], [4., 4., 2.]) 325 sim.atmosphere = Granular.createRegularAtmosphereGrid([4, 4, 2], [4., 4., 2.]) 326 Granular.addGrainCylindrical!(sim, [0.5, 0.5], .1, 1., verbose=verbose) 327 Granular.addGrainCylindrical!(sim, [0.7, 0.7], .1, 1., verbose=verbose) 328 Granular.addGrainCylindrical!(sim, [2.6, 2.5], .1, 1., verbose=verbose) 329 Granular.sortGrainsInGrid!(sim, sim.ocean, verbose=verbose) 330 Granular.setTimeStep!(sim) 331 Granular.setTotalTime!(sim, 1.0) 332 Granular.run!(sim, single_step=true, verbose=false) 333 @test sim.atmosphere.collocated_with_ocean_grid == true 334 @test sim.grains[1].ocean_grid_pos == [1, 1] 335 @test sim.grains[2].ocean_grid_pos == [1, 1] 336 @test sim.grains[3].ocean_grid_pos == [3, 3] 337 @test sim.ocean.grain_list[1, 1] == [1, 2] 338 @test sim.ocean.grain_list[2, 2] == [] 339 @test sim.ocean.grain_list[3, 3] == [3] 340 @test sim.grains[1].atmosphere_grid_pos == [1, 1] 341 @test sim.grains[2].atmosphere_grid_pos == [1, 1] 342 @test sim.grains[3].atmosphere_grid_pos == [3, 3] 343 @test sim.atmosphere.grain_list[1, 1] == [1, 2] 344 @test sim.atmosphere.grain_list[2, 2] == [] 345 @test sim.atmosphere.grain_list[3, 3] == [3] 346 347 @info "Testing automatic grid-size adjustment" 348 # ocean grid 349 sim = Granular.createSimulation() 350 @test_throws ErrorException Granular.fitGridToGrains!(sim, sim.ocean) 351 sim = Granular.createSimulation() 352 Granular.addGrainCylindrical!(sim, [0.0, 1.5], .5, 1., verbose=verbose) 353 Granular.addGrainCylindrical!(sim, [2.5, 5.5], 1., 1., verbose=verbose) 354 Granular.fitGridToGrains!(sim, sim.ocean, verbose=true) 355 @test sim.ocean.xq[1,1] ≈ -.5 356 @test sim.ocean.yq[1,1] ≈ 1.0 357 @test sim.ocean.xq[end,end] ≈ 3.5 358 @test sim.ocean.yq[end,end] ≈ 6.5 359 360 sim = Granular.createSimulation() 361 Granular.addGrainCylindrical!(sim, [0.5, 1.5], .5, 1., verbose=verbose) 362 Granular.addGrainCylindrical!(sim, [2.5, 4.5], .5, 1., verbose=verbose) 363 Granular.fitGridToGrains!(sim, sim.ocean, verbose=true) 364 @test sim.ocean.xq[1,1] ≈ 0. 365 @test sim.ocean.yq[1,1] ≈ 1. 366 @test sim.ocean.xq[end,end] ≈ 3. 367 @test sim.ocean.yq[end,end] ≈ 5. 368 369 sim = Granular.createSimulation() 370 Granular.addGrainCylindrical!(sim, [0.5, 0.0], .5, 1., verbose=verbose) 371 Granular.addGrainCylindrical!(sim, [2.0, 4.0], 1., 1., verbose=verbose) 372 Granular.fitGridToGrains!(sim, sim.ocean, padding=.5, verbose=true) 373 @test sim.ocean.xq[1,1] ≈ -.5 374 @test sim.ocean.yq[1,1] ≈ -1. 375 @test sim.ocean.xq[end,end] ≈ 3.5 376 @test sim.ocean.yq[end,end] ≈ 5.5 377 378 # atmosphere grid 379 sim = Granular.createSimulation() 380 @test_throws ErrorException Granular.fitGridToGrains!(sim, sim.atmosphere) 381 sim = Granular.createSimulation() 382 Granular.addGrainCylindrical!(sim, [0.0, 1.5], .5, 1., verbose=verbose) 383 Granular.addGrainCylindrical!(sim, [2.5, 5.5], 1., 1., verbose=verbose) 384 Granular.fitGridToGrains!(sim, sim.atmosphere, verbose=true) 385 @test sim.atmosphere.xq[1,1] ≈ -.5 386 @test sim.atmosphere.yq[1,1] ≈ 1.0 387 @test sim.atmosphere.xq[end,end] ≈ 3.5 388 @test sim.atmosphere.yq[end,end] ≈ 6.5 389 390 sim = Granular.createSimulation() 391 Granular.addGrainCylindrical!(sim, [0.5, 1.5], .5, 1., verbose=verbose) 392 Granular.addGrainCylindrical!(sim, [2.5, 4.5], .5, 1., verbose=verbose) 393 Granular.fitGridToGrains!(sim, sim.atmosphere, verbose=true) 394 @test sim.atmosphere.xq[1,1] ≈ 0. 395 @test sim.atmosphere.yq[1,1] ≈ 1. 396 @test sim.atmosphere.xq[end,end] ≈ 3. 397 @test sim.atmosphere.yq[end,end] ≈ 5. 398 399 sim = Granular.createSimulation() 400 Granular.addGrainCylindrical!(sim, [0.5, 0.0], .5, 1., verbose=verbose) 401 Granular.addGrainCylindrical!(sim, [2.0, 4.0], 1., 1., verbose=verbose) 402 Granular.fitGridToGrains!(sim, sim.atmosphere, padding=.5, verbose=true) 403 @test sim.atmosphere.xq[1,1] ≈ -.5 404 @test sim.atmosphere.yq[1,1] ≈ -1. 405 @test sim.atmosphere.xq[end,end] ≈ 3.5 406 @test sim.atmosphere.yq[end,end] ≈ 5.5 407 408 @info "Testing porosity estimation" 409 sim = Granular.createSimulation() 410 dx = 1.0; dy = 1.0 411 nx = 3; ny = 3 412 sim.ocean = Granular.createRegularOceanGrid([nx, ny, 1], [nx*dx, ny*dy, 1.]) 413 Granular.addGrainCylindrical!(sim, [1.5, 1.5], 0.5*dx, 1.0) 414 A_particle = π*(0.5*dx)^2 415 A_cell = dx^2 416 Granular.findPorosity!(sim, sim.ocean) 417 @test sim.ocean.porosity ≈ [1. 1. 1.; 418 1. (A_cell - A_particle)/A_cell 1.; 419 1. 1. 1]