Granular.jl

Julia package for granular dynamics simulation
git clone git://src.adamsgaard.dk/Granular.jl
Log | Files | Refs | README | LICENSE

commit 4947424934abc67010b9a1e184b5464373f1f94a
parent 0c238ad9a3d7510cd2b9b285d482bf0607083dba
Author: Anders Damsgaard <andersd@riseup.net>
Date:   Fri,  3 Nov 2017 14:13:40 -0400

add function to set BCs and add tests of BC functionality

Diffstat:
Msrc/atmosphere.jl | 4++--
Msrc/datatypes.jl | 36++++++++++++++++++++----------------
Msrc/grid.jl | 115+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Msrc/ocean.jl | 6+++---
Atest/periodic-boundaries.jl | 104+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mtest/runtests.jl | 1+
6 files changed, 245 insertions(+), 21 deletions(-)

diff --git a/src/atmosphere.jl b/src/atmosphere.jl @@ -17,7 +17,7 @@ function createEmptyAtmosphere() Array{Vector{Int}}(1, 1), - 0, 0, 0, 0, + 1, 1, 1, 1, false) end @@ -129,7 +129,7 @@ function createRegularAtmosphereGrid(n::Vector{Int}, zl, u, v, Array{Array{Int, 1}}(size(xh, 1), size(xh, 2)), - 0, 0, 0, 0, + 1, 1, 1, 1, false) end diff --git a/src/datatypes.jl b/src/datatypes.jl @@ -166,17 +166,17 @@ h-points. During read, the velocities are interpolated to the cell corners * `grain_list::Array{Float64, Int}`: indexes of grains contained in the ocean grid cells. * `bc_west::Integer`: Boundary condition type for the west edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic * `bc_south::Integer`: Boundary condition type for the south edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic * `bc_east::Integer`: Boundary condition type for the east edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic * `bc_north::Integer`: Boundary condition type for the north edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic =# mutable struct Ocean input_file::Any @@ -243,17 +243,17 @@ cell corners (q-points). * `grain_list::Array{Float64, Int}`: interface height relative to mean sea level [m], dimensions correspond to placement in `[xh, yh, zi, time]`. * `bc_west::Integer`: Boundary condition type for the west edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic * `bc_south::Integer`: Boundary condition type for the south edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic * `bc_east::Integer`: Boundary condition type for the east edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic * `bc_north::Integer`: Boundary condition type for the north edge of the grid. - 0: inactive, - 1: periodic + 1: inactive, + 2: periodic =# mutable struct Atmosphere input_file::Any @@ -306,3 +306,7 @@ mutable struct Simulation Nc_max::Int end + +# Mappings between boundary condition keys (Integers) and strings +const grid_bc_strings = ["inactive", "periodic"] +const grid_bc_flags = Dict(zip(grid_bc_strings, 1:length(grid_bc_strings))) diff --git a/src/grid.jl b/src/grid.jl @@ -617,3 +617,118 @@ function copyGridSortingInfo!(ocean::Ocean, atmosphere::Atmosphere, atmosphere.grain_list = deepcopy(ocean.grain_list) nothing end + +export setGridBoundaryConditions! +""" + setGridBoundaryConditions!(grid, grid_face, mode) + +Set boundary conditions for the granular phase at the edges of `Ocean` or +`Atmosphere` grids. The target boundary can be selected through the `grid_face` +argument, or the same boundary condition can be applied to all grid boundaries +at once. + +When the center coordinate of grains crosses an inactive boundary (`mode = +"inactive"`), the grain is disabled (`GrainCylindrical.enabled = false`). This +keeps the grain in memory, but stops it from moving or interacting with other +grains. *By default, all boundaries are inactive*. + +If the center coordinate of a grain crosses a periodic boundary (`mode = +periodic`), the grain is repositioned to the opposite side of the model domain. +Grains can interact mechanically across the periodic boundary. + +# Arguments +* `grid::Any`: `Ocean` or `Atmosphere` grid to apply the boundary condition to. +* `grid_face::String`: Grid face to apply the boundary condition to. Valid + values are any combination and sequence of `"west"` (-x), `"south"` (-y), + `"east"` (+x), `"north"` (+y). The values may be delimited in any way. + Also, and by default, all boundaries can be selected with `"all"` (-x, -y, + +x, +y), which overrides any other face selection. +* `mode::String`: Boundary behavior, accepted values are `"inactive"` and + `"periodic"`. You cannot specify more than one mode at a time, so if + several modes are desired as boundary conditions for the grid, several calls + to this function should be made. +* `verbose::Bool`: Confirm boundary conditions by reporting values to console. + +# Examples +Set all boundaries for the ocean grid to be periodic: + + setGridBoundaryConditions!(ocean, "periodic", "all") + +Set the south-north boundaries to be inactive, but the west-east boundaries to +be periodic: + + setGridBoundaryConditions!(ocean, "inactive", "south north") + setGridBoundaryConditions!(ocean, "periodic", "west east") + +""" +function setGridBoundaryConditions!(grid::Any, + mode::String, + grid_face::String = "all"; + verbose::Bool=true) + + something_changed = false + + if length(mode) <= 1 + error("The mode string is required ('$mode')") + end + + if !(mode in grid_bc_strings) + error("Mode '$mode' not recognized as a valid boundary condition type") + end + + if contains(grid_face, "west") + grid.bc_west = grid_bc_flags[mode] + something_changed = true + end + + if contains(grid_face, "south") + grid.bc_south = grid_bc_flags[mode] + something_changed = true + end + + if contains(grid_face, "east") + grid.bc_east = grid_bc_flags[mode] + something_changed = true + end + + if contains(grid_face, "north") + grid.bc_north = grid_bc_flags[mode] + something_changed = true + end + + if grid_face == "all" + grid.bc_west = grid_bc_flags[mode] + grid.bc_south = grid_bc_flags[mode] + grid.bc_east = grid_bc_flags[mode] + grid.bc_north = grid_bc_flags[mode] + something_changed = true + end + + if !something_changed + error("grid_face string '$grid_face' not understood, must be east, " * + "west, north, and/or south.") + end + + if verbose + reportGridBoundaryConditions(grid) + end + nothing +end + +export reportGridBoundaryConditions +""" + reportGridBoundaryConditions(grid) + +Report the boundary conditions for the grid to the console. +""" +function reportGridBoundaryConditions(grid::Any) + println("West (-x): " * grid_bc_strings[grid.bc_west] * + "\t($(grid.bc_west))") + println("East (+x): " * grid_bc_strings[grid.bc_east] * + "\t($(grid.bc_east))") + println("South (-y): " * grid_bc_strings[grid.bc_south] * + "\t($(grid.bc_south))") + println("North (+y): " * grid_bc_strings[grid.bc_north] * + "\t($(grid.bc_north))") + nothing +end diff --git a/src/ocean.jl b/src/ocean.jl @@ -18,7 +18,7 @@ function createEmptyOcean() zeros(1,1,1,1), zeros(1,1,1,1), Array{Array{Int, 1}}(1, 1), - 0, 0, 0, 0) + 1, 1, 1, 1) end export readOceanNetCDF @@ -61,7 +61,7 @@ function readOceanNetCDF(velocity_file::String, grid_file::String) h, e, Array{Array{Int, 1}}(size(xh, 1), size(xh, 2)), - 0, 0, 0, 0 + 1, 1, 1, 1 ) return ocean end @@ -242,7 +242,7 @@ function createRegularOceanGrid(n::Array{Int, 1}, zl, zi, u, v, h, e, Array{Array{Int, 1}}(size(xh, 1), size(xh, 2)), - 0, 0, 0, 0) + 1, 1, 1, 1) end export addOceanDrag! diff --git a/test/periodic-boundaries.jl b/test/periodic-boundaries.jl @@ -0,0 +1,104 @@ +#!/usr/bin/env julia + +info("#### $(basename(@__FILE__)) ####") + +verbose=false +i = 0 + +info("Testing assignment and reporting of grid boundary conditions") +ocean = Granular.createEmptyOcean() + +Test.@test ocean.bc_west == 1 +Test.@test ocean.bc_east == 1 +Test.@test ocean.bc_north == 1 +Test.@test ocean.bc_south == 1 + +const originalSTDOUT = STDOUT +(out_r, out_w) = redirect_stdout() +Granular.reportGridBoundaryConditions(ocean) +close(out_w) +redirect_stdout(originalSTDOUT) +output = convert(String, readavailable(out_r)) +Test.@test output == """West (-x): inactive\t(1) +East (+x): inactive\t(1) +South (-y): inactive\t(1) +North (+y): inactive\t(1) +""" + +(out_r, out_w) = redirect_stdout() +Granular.setGridBoundaryConditions!(ocean, "periodic", "south, west", verbose=true) +close(out_w) +redirect_stdout(originalSTDOUT) +output = convert(String, readavailable(out_r)) +Test.@test output == """West (-x): periodic\t(2) +East (+x): inactive\t(1) +South (-y): periodic\t(2) +North (+y): inactive\t(1) +""" +Test.@test ocean.bc_west == 2 +Test.@test ocean.bc_east == 1 +Test.@test ocean.bc_north == 1 +Test.@test ocean.bc_south == 2 + +(out_r, out_w) = redirect_stdout() +Granular.setGridBoundaryConditions!(ocean, "inactive", "all", verbose=false) +close(out_w) +redirect_stdout(originalSTDOUT) +output = convert(String, readavailable(out_r)) +Test.@test output == "" +Test.@test ocean.bc_west == 1 +Test.@test ocean.bc_east == 1 +Test.@test ocean.bc_north == 1 +Test.@test ocean.bc_south == 1 + +(out_r, out_w) = redirect_stdout() +Granular.setGridBoundaryConditions!(ocean, "periodic", "all") +close(out_w) +output = convert(String, readavailable(out_r)) +redirect_stdout(originalSTDOUT) +Test.@test output == """West (-x): periodic\t(2) +East (+x): periodic\t(2) +South (-y): periodic\t(2) +North (+y): periodic\t(2) +""" +Test.@test ocean.bc_west == 2 +Test.@test ocean.bc_east == 2 +Test.@test ocean.bc_north == 2 +Test.@test ocean.bc_south == 2 + +(out_r, out_w) = redirect_stdout() +Granular.setGridBoundaryConditions!(ocean, "inactive") +close(out_w) +output = convert(String, readavailable(out_r)) +redirect_stdout(originalSTDOUT) +Test.@test output == """West (-x): inactive\t(1) +East (+x): inactive\t(1) +South (-y): inactive\t(1) +North (+y): inactive\t(1) +""" +Test.@test ocean.bc_west == 1 +Test.@test ocean.bc_east == 1 +Test.@test ocean.bc_north == 1 +Test.@test ocean.bc_south == 1 + +(out_r, out_w) = redirect_stdout() +Granular.setGridBoundaryConditions!(ocean, "periodic") +close(out_w) +output = convert(String, readavailable(out_r)) +redirect_stdout(originalSTDOUT) +Test.@test output == """West (-x): periodic\t(2) +East (+x): periodic\t(2) +South (-y): periodic\t(2) +North (+y): periodic\t(2) +""" +Test.@test ocean.bc_west == 2 +Test.@test ocean.bc_east == 2 +Test.@test ocean.bc_north == 2 +Test.@test ocean.bc_south == 2 + +Test.@test_throws ErrorException Granular.setGridBoundaryConditions!(ocean, + "periodic", + "asdf") + +Test.@test_throws ErrorException Granular.setGridBoundaryConditions!(ocean, + "asdf") diff --git a/test/runtests.jl b/test/runtests.jl @@ -13,6 +13,7 @@ include("netcdf.jl") include("vtk.jl") include("jld.jl") include("grid.jl") +include("periodic-boundaries.jl") include("ocean.jl") include("atmosphere.jl") include("memory-management.jl")