Granular.jl

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

commit 1d1b5457fd62a9751f2983af9432ade2cfe9b8f5
parent 600cc4b3910ffc0a17661f3061ea3c659c9696aa
Author: Anders Damsgaard <andersd@riseup.net>
Date:   Thu, 27 Apr 2017 15:28:06 -0400

export functions from module namespace, improve docstrings

Diffstat:
Msrc/contact_search.jl | 7+++++++
Msrc/grid.jl | 7+++++++
Msrc/icefloe.jl | 35+++++++++++++++++++++++++++++++++--
Msrc/interaction.jl | 3+++
Msrc/io.jl | 1+
Msrc/simulation.jl | 6++++++
Msrc/temporal.jl | 23+++++++++++++++++++----
Msrc/temporal_integration.jl | 15+++++++++++----
8 files changed, 87 insertions(+), 10 deletions(-)

diff --git a/src/contact_search.jl b/src/contact_search.jl @@ -1,4 +1,5 @@ ## Contact mapping +export findContacts """ Top-level function to perform an inter-ice floe contact search, based on ice floe linear positions and contact radii. @@ -16,6 +17,7 @@ function findContacts!(simulation::Simulation, end end +export interIceFloePositionVector function interIceFloePositionVector(simulation::Simulation, i::Integer, j::Integer) return simulation.ice_floes[j].lin_pos - simulation.ice_floes[i].lin_pos @@ -31,6 +33,11 @@ function findOverlap(simulation::Simulation, i::Integer, j::Integer, simulation.ice_floes[j].contact_radius) end +export findContactsAllToAll +""" +Perform an O(n^2) all-to-all contact search between all ice floes in the +`simulation` object. +""" function findContactsAllToAll!(simulation::Simulation) for i = 1:length(simulation.ice_floes) diff --git a/src/grid.jl b/src/grid.jl @@ -28,6 +28,7 @@ function bilinearInterpolation(field::Array{Float64, 4}, field[i-1, j-1, k, it]*(1. - x_tilde))*(1. - y_tilde) end +export sortIceFloesInOceanGrid """ Find ice-floe positions in ocean grid, based on their center positions. """ @@ -49,6 +50,7 @@ function sortIceFloesInOceanGrid!(simulation::Simulation, verbose=true) end end +export findCellContainingPoint """ Returns the `i`, `j` index of the ocean grid cell containing the `point`. """ @@ -71,6 +73,7 @@ function findCellContainingPoint(ocean::Ocean, point::Array{float, 1}) end +export isPointInCell """ Check if a 2d point is contained inside a cell from the ocean grid. The function uses either an area-based approach (`method = "Area"`), or a @@ -106,6 +109,7 @@ function isPointInCell(ocean::Ocean, i::Int, j::Int, point::Array{float, 1}; end end +export getCellCornerCoordinates """ Returns ocean-grid corner coordinates in the following order (south-west corner, south-east corner, north-east corner, north-west corner). @@ -118,6 +122,7 @@ function getCellCornerCoordinates(ocean::Ocean, i::Int, j::Int) return sw, se, ne, nw end +export areaOfTriangle "Returns the area of an triangle with corner coordinates `a`, `b`, and `c`." function areaOfTriangle(a::Array{float, 1}, b::Array{float, 1}, @@ -129,6 +134,7 @@ function areaOfTriangle(a::Array{float, 1}, ) end +export areaOfQuadrilateral """ Returns the area of a quadrilateral with corner coordinates `a`, `b`, `c`, and `d`. Corners `a` and `c` should be opposite of each other, the same must be @@ -142,6 +148,7 @@ function areaOfQuadrilateral(a::Array{float, 1}, return areaOfTriangle(a, b, c) + areaOfTriangle(c, d, a) end +export conformalQuadrilateralCoordinates """ Returns the non-dimensional coordinates `[x_tilde, y_tilde]` of a point `p` within a quadrilateral with corner coordinates `A`, `B`, `C`, and `D`. diff --git a/src/icefloe.jl b/src/icefloe.jl @@ -1,5 +1,6 @@ ## Manage icefloes in the model +export addIceFloeCylindrical """ Adds a grain to the simulation. Example: @@ -111,30 +112,47 @@ function addIceFloeCylindrical(simulation::Simulation, addIceFloe!(simulation, icefloe, verbose) end +export iceFloeCircumreference +"Returns the circumreference of the ice floe" function iceFloeCircumreference(icefloe::IceFloeCylindrical) return pi*icefloe.areal_radius*2. end +export iceFloeHorizontalSurfaceArea +"Returns the top or bottom (horizontal) surface area of the ice floe" function iceFloeHorizontalSurfaceArea(icefloe::IceFloeCylindrical) return pi*icefloe.areal_radius^2. end +export iceFloeSideSurfaceArea +"Returns the surface area of the ice-floe sides" function iceFloeSideSurfaceArea(icefloe::IceFloeCylindrical) return iceFloeCircumreference(icefloe)*icefloe.thickness end +export iceFloeVolume +"Returns the volume of the ice floe" function iceFloeVolume(icefloe::IceFloeCylindrical) return iceFloeHorizontalSurfaceArea(icefloe)*icefloe.thickness end +export iceFloeMass +"Returns the mass of the ice floe" function iceFloeMass(icefloe::IceFloeCylindrical) return iceFloeVolume(icefloe)*icefloe.density end +export iceFloeMomentOfInertia +"Returns the moment of inertia of the ice floe" function iceFloeMomentOfInertia(icefloe::IceFloeCylindrical) return 0.5*iceFloeMass(icefloe)*icefloe.areal_radius^2. end +export convertIceFloeDataToArrays +""" +Gathers all ice-floe parameters from the `IceFloeCylindrical` type in continuous +arrays in an `IceFloeArrays` structure. +""" function convertIceFloeDataToArrays(simulation::Simulation) ifarr = IceFloeArrays( @@ -220,24 +238,37 @@ function convertIceFloeDataToArrays(simulation::Simulation) return ifarr end +export iceFloeKineticTranslationalEnergy +"Returns the translational kinetic energy of the ice floe" function iceFloeKineticTranslationalEnergy(icefloe::IceFloeCylindrical) return 0.5*icefloe.mass*norm(icefloe.lin_vel)^2. end +export totalIceFloeKineticTranslationalEnergy +""" +Returns the sum of translational kinetic energies of all ice floes in a +simulation +""" function totalIceFloeKineticTranslationalEnergy(simulation::Simulation) - E_sum = 0 + E_sum = 0. for icefloe in simulation.ice_floes E_sum += iceFloeKineticTranslationalEnergy(icefloe) end return E_sum end +export iceFloeKineticRotationalEnergy +"Returns the rotational kinetic energy of the ice floe" function iceFloeKineticRotationalEnergy(icefloe::IceFloeCylindrical) return 0.5*icefloe.moment_of_inertia*norm(icefloe.ang_vel)^2. end +export totalIceFloeKineticRotationalEnergy +""" +Returns the sum of rotational kinetic energies of all ice floes in a simulation +""" function totalIceFloeKineticRotationalEnergy(simulation::Simulation) - E_sum = 0 + E_sum = 0. for icefloe in simulation.ice_floes E_sum += iceFloeKineticRotationalEnergy(icefloe) end diff --git a/src/interaction.jl b/src/interaction.jl @@ -1,5 +1,6 @@ ## Interaction functions +export interact """ Resolve mechanical interaction between all particle pairs. """ @@ -14,6 +15,7 @@ function interact!(simulation::Simulation) end end +export interactIceFloes """ Resolve an grain-to-grain interaction using a prescibed contact law. This function adds the compressive force of the interaction to the ice floe @@ -45,6 +47,7 @@ function interactIceFloes!(simulation::Simulation, norm(force)/simulation.ice_floes[j].side_surface_area; end +export interactNormalLinearElastic """ Resolves linear-elastic interaction between two ice floes in the contact-normal direction. diff --git a/src/io.jl b/src/io.jl @@ -3,6 +3,7 @@ import NetCDF ## IO functions +export writeVTK """ Write a VTK file to disk containing all ice floes in the `simulation` in an unstructured mesh (file type `.vtu`). These files can be read by ParaView and diff --git a/src/simulation.jl b/src/simulation.jl @@ -1,5 +1,6 @@ ## General simulation functions +export createSimulation """ createSimulation([id::String="unnamed", time_iteration::Int=0, @@ -43,6 +44,7 @@ function createSimulation(;id::String="unnamed", ocean) end +export run """ run!(simulation[, verbose::Bool = true, @@ -131,6 +133,7 @@ function run!(simulation::Simulation; end end +export addIceFloe "Add an `icefloe` to the `simulation` object. If `verbose` is true, a short confirmation message will be printed to stdout`." function addIceFloe!(simulation::Simulation, @@ -144,6 +147,7 @@ function addIceFloe!(simulation::Simulation, end end +export removeIceFloe "Remove ice floe with index `i` from the `simulation` object." function removeIceFloe!(simulation::Simulation, i::Integer) if i < 1 @@ -153,6 +157,7 @@ function removeIceFloe!(simulation::Simulation, i::Integer) delete!(simulation.ice_floes, i) end +export zeroForcesAndTorques "Sets the `force` and `torque` values of all ice floes to zero." function zeroForcesAndTorques!(simulation::Simulation) for icefloe in simulation.ice_floes @@ -162,6 +167,7 @@ function zeroForcesAndTorques!(simulation::Simulation) end end +export reportSimulationTimeToStdout "Prints the current simulation time and total time to standard out" function reportSimulationTimeToStdout(simulation::Simulation) print("\r t = ", simulation.time, '/', simulation.time_total, diff --git a/src/temporal.jl b/src/temporal.jl @@ -1,3 +1,4 @@ +export setTotalTime """ setTotalTime!(simulation::Simulation, t::float) @@ -11,6 +12,7 @@ function setTotalTime!(simulation::Simulation, t::float) simulation.time_total = t end +export setCurrentTime """ setCurrentTime!(simulation::Simulation, t::float) @@ -24,6 +26,7 @@ function setCurrentTime!(simulation::Simulation, t::float) simulation.time = t end +export incrementCurrentTime """ incrementCurrentTime!(simulation::Simulation, t::float) @@ -37,19 +40,28 @@ function incrementCurrentTime!(simulation::Simulation, t::float) simulation.time += t end -function setOutputFileInterval!(simulation::Simulation, t::float) - if t <= 0.0 - info("The output file interval should be a positive value (t = $t s)") +export setOutputFileInterval +""" + setOutputFileInterval!(simulation::Simulation, t::float) + +Sets the simulation-time interval between output files are written to disk. If +this value is zero or negative, no output files will be written. +""" +function setOutputFileInterval!(simulation::Simulation, t::float; verbose=true) + if t <= 0.0 && verbose + info("No output files will be written") end simulation.file_time_step = t end +export disableOutputFiles "Disables the write of output files to disk during a simulation." function disableOutputFiles!(simulation::Simulation) simulation.file_time_step = 0.0 end -"Checks if simulation parameters are of reasonable values." +export checkTimeParameters +"Checks if simulation temporal parameters are of reasonable values." function checkTimeParameters(simulation::Simulation) if simulation.time_total <= 0.0 || simulation.time_total <= simulation.time error("Total time should be positive and larger than current time.\n", @@ -61,6 +73,7 @@ function checkTimeParameters(simulation::Simulation) end end +export findSmallestIceFloeMass "Finds the smallest mass of all ice floes in a simulation. Used to determine the optimal time step length." function findSmallestIceFloeMass(simulation::Simulation) @@ -76,6 +89,7 @@ function findSmallestIceFloeMass(simulation::Simulation) return m_min, i_min end +export findLargestIceFloeStiffness "Finds the largest elastic stiffness of all ice floes in a simulation. Used to determine the optimal time step length." function findLargestIceFloeStiffness(simulation::Simulation) @@ -97,6 +111,7 @@ function findLargestIceFloeStiffness(simulation::Simulation) return k_n_max, k_t_max, i_n_max, i_t_max end +export setTimeStep """ Find the computational time step length suitable given the grain radii, contact stiffnesses, and grain density. Uses the scheme by Radjaii et al. 2011. diff --git a/src/temporal_integration.jl b/src/temporal_integration.jl @@ -1,3 +1,4 @@ +export updateIceFloeKinematics """ updateIceFloeKinematics!(simulation::Simulation[, method::String = "Three-term Taylor"]) @@ -31,8 +32,11 @@ function updateIceFloeKinematics!(simulation::Simulation; end end -"Use a two-term Taylor expansion for integrating the kinematic degrees of -freedom for an `icefloe`." +export updateIceFloeKinematicsTwoTermTaylor +""" +Use a two-term Taylor expansion for integrating the kinematic degrees of freedom +for an `icefloe`. +""" function updateIceFloeKinematicsTwoTermTaylor(icefloe::IceFloeCylindrical, simulation::Simulation) icefloe.lin_acc = icefloe.force/icefloe.mass @@ -56,8 +60,11 @@ function updateIceFloeKinematicsTwoTermTaylor(icefloe::IceFloeCylindrical, icefloe.ang_vel += icefloe.ang_acc * simulation.time_step end -"Use a three-term Taylor expansion for integrating the kinematic degrees of -freedom for an `icefloe`." +export updateIceFloeKinematicsThreeTermTaylor +""" +Use a three-term Taylor expansion for integrating the kinematic degrees of +freedom for an `icefloe`. +""" function updateIceFloeKinematicsThreeTermTaylor(icefloe::IceFloeCylindrical, simulation::Simulation)