Granular.jl

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

commit 2e147d502dba2dd18246f4c0d7f92c40a235631d
parent a424c3439fff74481ce41e7d33e29dd7d22d4f81
Author: Anders Damsgaard <andersd@riseup.net>
Date:   Mon,  2 Apr 2018 15:02:49 -0400

Fix srand calls and print_with_color

Diffstat:
Mexamples/double_gyre.jl | 3++-
Mexamples/image.jl | 5+++--
Mexamples/logo.jl | 5+++--
Mexamples/strait.jl | 3++-
Msrc/grain.jl | 71+++++++++++++++++++++++++++++++++++------------------------------------
Msrc/grid.jl | 40+++++++++++++++++++++-------------------
Msrc/io.jl | 27+++++++++++++--------------
Msrc/ocean.jl | 2+-
Msrc/packing.jl | 23++++++++++++-----------
Msrc/simulation.jl | 1-
Mtest/jld.jl | 3++-
Mtest/profiling.jl | 13+++++++------
Mtest/util.jl | 3++-
13 files changed, 103 insertions(+), 96 deletions(-)

diff --git a/examples/double_gyre.jl b/examples/double_gyre.jl @@ -1,6 +1,7 @@ #!/usr/bin/env julia import Granular import Compat +using Compat.Random sim = Granular.createSimulation(id="double_gyre") @@ -59,7 +60,7 @@ Compat.@info "added $(n_walls) fixed ice floes as walls" # Initialize ice floes everywhere floe_padding = .5*r noise_amplitude = .8*floe_padding -Compat.srand(1) +srand(1) for y in (4.*r + noise_amplitude):(2.*r + floe_padding):(L[2] - 4.*r - noise_amplitude) diff --git a/examples/image.jl b/examples/image.jl @@ -4,6 +4,7 @@ import Granular import FileIO import Colors import Compat +using Compat.Random const verbose = true @@ -89,12 +90,12 @@ if forcing == "gyres" end elseif forcing == "down" || forcing == "sandpile" - Compat.srand(1) + srand(1) sim.ocean.u[:, :, 1, 1] = (rand(nx+1, ny+1) - .5)*.1 sim.ocean.v[:, :, 1, 1] = -Ly/5. elseif forcing == "convergent" - Compat.srand(1) + srand(1) sim.ocean.u[:, :, 1, 1] = (rand(nx+1, ny+1) - .5)*.1 for j=1:size(sim.ocean.u, 2) sim.ocean.v[:, j, 1, 1] = -(j/ny - .5)*10. diff --git a/examples/logo.jl b/examples/logo.jl @@ -2,6 +2,7 @@ import Granular import Compat +using Compat.Random const verbose = true @@ -119,12 +120,12 @@ if forcing == "gyres" end elseif forcing == "down" - Compat.srand(1) + srand(1) sim.ocean.u[:, :, 1, 1] = (rand(nx+1, ny+1) - .5)*.1 sim.ocean.v[:, :, 1, 1] = -5. elseif forcing == "convergent" - Compat.srand(1) + srand(1) sim.ocean.u[:, :, 1, 1] = (rand(nx+1, ny+1) - .5)*.1 for j=1:size(sim.ocean.u, 2) sim.ocean.v[:, j, 1, 1] = -(j/ny - .5)*10. diff --git a/examples/strait.jl b/examples/strait.jl @@ -1,6 +1,7 @@ #!/usr/bin/env julia import SeaIce import Compat +using Compat.Random sim = SeaIce.createSimulation(id="strait") n = [10, 10, 2] @@ -87,7 +88,7 @@ dy = sqrt((2.*r_walls)^2. - dx^2.) spacing_to_boundaries = 4.*r floe_padding = .5*r noise_amplitude = floe_padding -Compat.srand(1) +srand(1) for y in (L[2] - r - noise_amplitude):(-(2.*r + floe_padding)):((L[2] - Ly_constriction)/2. + Ly_constriction) for x in (r + noise_amplitude):(2.*r + floe_padding):(Lx - r - diff --git a/src/grain.jl b/src/grain.jl @@ -346,27 +346,27 @@ function convertGrainDataToArrays(simulation::Simulation) ifarr = GrainArrays( # Material properties ## density - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), # Geometrical properties ## thickness - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## contact_radius - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## areal_radius - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## circumreference - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## horizontal_surface_area - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## side_surface_area - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## volume - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## mass - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## moment_of_inertia - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), # Linear kinematic degrees of freedom along horiz plane ## lin_pos @@ -394,58 +394,58 @@ function convertGrainDataToArrays(simulation::Simulation) # Kinematic constraint flags ## fixed - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), ## allow_x_acc - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), ## allow_y_acc - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), ## rotating - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), ## enabled - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), # Rheological parameters ## contact_stiffness_normal - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## contact_stiffness_tangential - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## contact_viscosity_normal - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## contact_viscosity_tangential - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## contact_static_friction - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## contact_dynamic_friction - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## youngs_modulus - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## poissons_ratio - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## tensile_strength - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## shear_strength - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## strength_heal_rate - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## fracture_toughness - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), # Ocean/atmosphere interaction parameters ## ocean_drag_coeff_vert - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## ocean_drag_coeff_horiz - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## atmosphere_drag_coeff_vert - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## atmosphere_drag_coeff_horiz - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), # Interaction ## pressure - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), ## n_contacts - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), ## granular_stress zeros(Float64, 3, length(simulation.grains)), @@ -455,11 +455,11 @@ function convertGrainDataToArrays(simulation::Simulation) zeros(Float64, 3, length(simulation.grains)), ## thermal_energy - Array{Float64}(length(simulation.grains)), + Array{Float64}(undef, length(simulation.grains)), # Visualization parameters ## color - Array{Int}(length(simulation.grains)), + Array{Int}(undef, length(simulation.grains)), ) @@ -607,7 +607,6 @@ function deleteGrainArrays!(ifarr::GrainArrays) ifarr.thermal_energy = f1 ifarr.color = i1 - gc() nothing end diff --git a/src/grid.jl b/src/grid.jl @@ -1,5 +1,6 @@ -import Compat +using Compat using Compat.LinearAlgebra +using Compat.Random """ bilinearInterpolation(field, x_tilde, y_tilde, i, j, k, it) @@ -256,10 +257,10 @@ found the function returns `(0,0)`. not, can be "Conformal" (default) or "Area". """ function findCellContainingPoint(grid::Any, point::Vector{Float64}, - sw = Vector{Float64}(2), - se = Vector{Float64}(2), - ne = Vector{Float64}(2), - nw = Vector{Float64}(2); + sw = Vector{Float64}(undef, 2), + se = Vector{Float64}(undef, 2), + ne = Vector{Float64}(undef, 2), + nw = Vector{Float64}(undef, 2); method::String="Conformal") for i=1:size(grid.xh, 1) for j=1:size(grid.yh, 2) @@ -300,10 +301,10 @@ more robust. This function returns `true` or `false`. """ function isPointInCell(grid::Any, i::Int, j::Int, point::Vector{Float64}, - sw::Vector{Float64} = Vector{Float64}(2), - se::Vector{Float64} = Vector{Float64}(2), - ne::Vector{Float64} = Vector{Float64}(2), - nw::Vector{Float64} = Vector{Float64}(2); + sw::Vector{Float64} = Vector{Float64}(undef, 2), + se::Vector{Float64} = Vector{Float64}(undef, 2), + ne::Vector{Float64} = Vector{Float64}(undef, 2), + nw::Vector{Float64} = Vector{Float64}(undef, 2); method::String="Conformal") if grid.regular_grid @@ -351,10 +352,10 @@ area-based approach (`method = "Area"`), or a conformal mapping approach function returns `true` or `false`. """ function isPointInGrid(grid::Any, point::Vector{Float64}, - sw::Vector{Float64} = Vector{Float64}(2), - se::Vector{Float64} = Vector{Float64}(2), - ne::Vector{Float64} = Vector{Float64}(2), - nw::Vector{Float64} = Vector{Float64}(2); + sw::Vector{Float64} = Vector{Float64}(undef, 2), + se::Vector{Float64} = Vector{Float64}(undef, 2), + ne::Vector{Float64} = Vector{Float64}(undef, 2), + nw::Vector{Float64} = Vector{Float64}(undef, 2); method::String="Conformal") #sw, se, ne, nw = getCellCornerCoordinates(grid.xq, grid.yq, i, j) @@ -585,7 +586,7 @@ function findEmptyPositionInGridCell(simulation::Simulation, for i_iter=1:n_iter overlap_found = false - Compat.srand(i*j*seed*i_iter) + srand(i*j*seed*i_iter) # generate random candidate position x_tilde = rand() y_tilde = rand() @@ -625,7 +626,8 @@ function findEmptyPositionInGridCell(simulation::Simulation, if overlap < 0. if verbose - Compat.@info "overlap with $grain_idx in cell $i,$j" + Compat.@info "overlap with $grain_idx in " * + "cell $i,$j" end overlap_found = true break @@ -1048,10 +1050,10 @@ function findPorosity!(sim::Simulation, grid::Any; verbose::Bool=true) sortGrainsInGrid!(sim, grid, verbose=verbose) end - sw = Vector{Float64}(2) - se = Vector{Float64}(2) - ne = Vector{Float64}(2) - nw = Vector{Float64}(2) + @compat sw = Vector{Float64}(undef, 2) + @compat se = Vector{Float64}(undef, 2) + @compat ne = Vector{Float64}(undef, 2) + @compat nw = Vector{Float64}(undef, 2) cell_area = 0.0 p = zeros(2) r = 0.0 diff --git a/src/io.jl b/src/io.jl @@ -3,7 +3,7 @@ import Compat using Compat.LinearAlgebra hasJLD = false -if typeof(Pkg.installed("JLD")) == VersionNumber +if typeof(Compat.Pkg.installed("JLD")) == VersionNumber import JLD hasJLD = true end @@ -240,31 +240,31 @@ function status(folder::String="."; if write_header println("--------------------------------------" * "--------------------------------------") - print_with_color(:default, "simulation folder \t") - print_with_color(time_color, " time \t") - print_with_color(percentage_color, " completed ") - print_with_color(lastfile_color, "last file \n") + Compat.printstyled("simulation folder \t", color=:default) + Compat.printstyled(" time \t", color=time_color) + Compat.printstyled(" completed ", color=percentage_color) + Compat.printstyled("last file \n", color=lastfile_color) println("--------------------------------------" * "--------------------------------------") end for file in status_files data = readdlm(file) - id = replace(file, ".status.txt", "") - id = replace(id, "./", "") - id = replace(id, r".*/", "") + id = Compat.replace(file, ".status.txt" => "") + id = Compat.replace(id, "./" => "") + id = Compat.replace(id, r".*/" => "") time_s = @sprintf "%6.2fs" data[1] time_h = @sprintf "%5.1fh" data[1]/(60. * 60.) percentage = @sprintf "%3.0f%%" data[2] lastfile = @sprintf "%5d" data[3] if data[2] < 99. - print_with_color(id_color_in_progress, "$id \t") + Compat.printstyled("$id \t", color=id_color_in_progress) else - print_with_color(id_color_complete, "$id \t") + Compat.printstyled("$id \t", color=id_color_complete) end - print_with_color(time_color, "$time_s ($time_h) \t") - print_with_color(percentage_color, "$percentage \t") - print_with_color(lastfile_color, "$lastfile \n") + Compat.printstyled("$time_s ($time_h) \t", color=time_color) + Compat.printstyled("$percentage \t", color=percentage_color) + Compat.printstyled("$lastfile \n", color=lastfile_color) if visualize sim = createSimulation(id) @@ -443,7 +443,6 @@ function writeGrainVTK(simulation::Simulation, deleteGrainArrays!(ifarr) ifarr = 0 - gc() outfiles = WriteVTK.vtk_save(vtkfile) if verbose diff --git a/src/ocean.jl b/src/ocean.jl @@ -3,7 +3,7 @@ using Compat.Test using Compat.LinearAlgebra hasNetCDF = false -if typeof(Pkg.installed("NetCDF")) == VersionNumber +if typeof(Compat.Pkg.installed("NetCDF")) == VersionNumber import NetCDF hasNetCDF = true else diff --git a/src/packing.jl b/src/packing.jl @@ -1,6 +1,7 @@ ## Functions for creating grain packings import Compat using Compat.LinearAlgebra +using Compat.Random export regularPacking! """ @@ -44,7 +45,7 @@ function regularPacking!(simulation::Simulation, r_rand = 0. pos = zeros(2) h = .5 # disc tickness - Compat.srand(seed) + srand(seed) if tiling == "square" dx = r_max * 2. * (1. + padding_factor) # cell size @@ -178,7 +179,7 @@ function irregularPacking!(simulation::Simulation; seed::Integer=1, plot_during_packing::Bool=false, verbose::Bool=true) - Compat.srand(seed) + srand(seed) active_list = Int[] # list of points to originate search from i = 0 @@ -207,9 +208,9 @@ function irregularPacking!(simulation::Simulation; sortGrainsInGrid!(simulation, grid) push!(active_list, 1) else - for i=1:length(simulation.grains) - simulation.grains[i].color=1 - push!(active_list, i) + for idx=1:length(simulation.grains) + simulation.grains[idx].color=1 + push!(active_list, idx) end end @@ -366,7 +367,7 @@ function rasterPacking!(sim::Simulation, h = .5 # disc tickness dx = r_max * 2. * (1. + padding_factor) # cell size dx_padding = r_max * 2. * padding_factor - Compat.srand(seed) + srand(seed) np_init = length(sim.grains) @@ -450,15 +451,15 @@ function rasterMap(sim::Simulation, dx::Real) #i, j = Int.(floor.((grain.lin_pos - origo) ./ dx)) + [1,1] # Find corner indexes for box spanning the grain - min_i, min_j = Int.(floor.((grain.lin_pos - origo - - grain.contact_radius) ./ dx)) + [1,1] - max_i, max_j = Int.(floor.((grain.lin_pos - origo + - grain.contact_radius) ./ dx)) + [1,1] + min_i, min_j = Int.(floor.((grain.lin_pos - origo .- + grain.contact_radius) ./ dx)) .+ [1,1] + max_i, max_j = Int.(floor.((grain.lin_pos - origo .+ + grain.contact_radius) ./ dx)) .+ [1,1] # For each cell in box, check if the grain is contained for i = min_i:max_i for j = min_j:max_j - cell_mid_point = dx .* Vector{Float64}([i,j]) - 0.5 * dx + cell_mid_point = dx .* Vector{Float64}([i,j]) .- 0.5 * dx if (norm(cell_mid_point - grain.lin_pos) - grain.contact_radius < dist) diff --git a/src/simulation.jl b/src/simulation.jl @@ -201,7 +201,6 @@ function run!(simulation::Simulation; reportSimulationTimeToStdout(simulation) println() end - gc() nothing end diff --git a/test/jld.jl b/test/jld.jl @@ -1,9 +1,10 @@ #!/usr/bin/env julia +import Compat Compat.@info "#### $(basename(@__FILE__)) ####" Compat.@info "Determining if JLD is installed" -if typeof(Pkg.installed("JLD")) == VersionNumber +if typeof(Compat.Pkg.installed("JLD")) == VersionNumber Compat.@info "JLD found, proceeding with JLD-specific tests" Compat.@info "Writing simple simulation to JLD file" diff --git a/test/profiling.jl b/test/profiling.jl @@ -48,17 +48,18 @@ function timeSingleStepInDenseSimulation(nx::Int; verbose::Bool=true, fixed=fixed, verbose=false) end end - print_with_color(:green, "number of grains: $(length(sim.grains))\n") + Compat.printstyled("number of grains: $(length(sim.grains))\n", + color=:green) if grid_sorting if include_atmosphere - print_with_color(:green, "using cell-based spatial decomposition " * - " (ocean + atmosphere)\n") + Compat.printstyled("using cell-based spatial decomposition " * + " (ocean + atmosphere)\n", color=:green) else - print_with_color(:green, "using cell-based spatial " * - "decomposition (ocean)\n") + Compat.printstyled("using cell-based spatial " * + "decomposition (ocean)\n", color=:green) end else - print_with_color(:green, "using all-to-all contact search\n") + Compat.printstyled("using all-to-all contact search\n", color=:green) end Granular.setTotalTime!(sim, 1.0) diff --git a/test/util.jl b/test/util.jl @@ -1,5 +1,6 @@ #!/usr/bin/env julia import Compat +using Compat.Random Compat.@info "#### $(basename(@__FILE__)) ####" @@ -14,7 +15,7 @@ Compat.@info "Testing power-law RNG" @test 5 == length(Granular.randpower(5)) @test (5,) == size(Granular.randpower(5)) -Compat.srand(1) +srand(1) for i=1:10^5 @test 0. <= Granular.randpower() <= 1. @test 0. <= Granular.randpower(1, 1., 0., 1.) <= 1.