Granular.jl

Julia package for granular dynamics simulation
git clone git://src.adamsgaard.dk/Granular.jl # fast
git clone https://src.adamsgaard.dk/Granular.jl.git # slow
Log | Files | Refs | README | LICENSE Back to index

collision-2floes-oblique.jl (25824B)


      1 #!/usr/bin/env julia
      2 
      3 # Check for conservation of kinetic energy (=momentum) during a normal collision 
      4 # between two ice cylindrical grains 
      5 
      6 verbose=false
      7 
      8 @info "## Contact-normal elasticity only"
      9 @info "# One ice floe fixed"
     10 sim = Granular.createSimulation(id="test")
     11 Granular.addGrainCylindrical!(sim, [0., 10.], 10., 1., verbose=verbose)
     12 Granular.addGrainCylindrical!(sim, [19., 0.], 10., 1., verbose=verbose)
     13 sim.grains[1].lin_vel[1] = 0.1
     14 sim.grains[1].contact_dynamic_friction = 0.
     15 sim.grains[2].contact_dynamic_friction = 0.
     16 sim.grains[2].fixed = true
     17 
     18 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
     19 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
     20 
     21 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
     22 # should become more correct
     23 
     24 Granular.setTotalTime!(sim, 30.0)
     25 #sim.file_time_step = 1.
     26 sim_init = deepcopy(sim)
     27 
     28 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
     29 Granular.setTimeStep!(sim, epsilon=0.07)
     30 tol = 0.1
     31 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
     32 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
     33 
     34 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
     35 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
     36 E_thermal_final = Granular.totalGrainThermalEnergy(sim)
     37 @test E_kin_lin_init ≈ E_kin_lin_final+E_thermal_final atol=E_kin_lin_init*tol
     38 @test E_kin_rot_init ≈ E_kin_rot_final
     39 
     40 
     41 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
     42 sim = deepcopy(sim_init)
     43 Granular.setTimeStep!(sim, epsilon=0.007)
     44 tol = 0.01
     45 @info "Relative tolerance: $(tol*100.)%"
     46 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
     47 
     48 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
     49 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
     50 E_thermal_final = Granular.totalGrainThermalEnergy(sim)
     51 @test E_kin_lin_init ≈ E_kin_lin_final+E_thermal_final atol=E_kin_lin_init*tol
     52 @test E_kin_rot_init ≈ E_kin_rot_final
     53 
     54 
     55 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
     56 sim = deepcopy(sim_init)
     57 Granular.setTimeStep!(sim, epsilon=0.07)
     58 tol = 0.01
     59 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
     60 Granular.run!(sim, temporal_integration_method="Three-term Taylor", verbose=verbose)
     61 
     62 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
     63 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
     64 E_thermal_final = Granular.totalGrainThermalEnergy(sim)
     65 @test E_kin_lin_init ≈ E_kin_lin_final+E_thermal_final atol=E_kin_lin_init*tol
     66 @test E_kin_rot_init ≈ E_kin_rot_final
     67 
     68 @info "# Ice floes free to move"
     69 
     70 sim = Granular.createSimulation(id="test")
     71 Granular.addGrainCylindrical!(sim, [0., 10.], 10., 1., verbose=verbose)
     72 Granular.addGrainCylindrical!(sim, [19.0, 0.], 10., 1., verbose=verbose)
     73 sim.grains[1].lin_vel[1] = 0.1
     74 sim.grains[1].contact_dynamic_friction = 0.
     75 sim.grains[2].contact_dynamic_friction = 0.
     76 
     77 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
     78 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
     79 
     80 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
     81 # should become more correct
     82 
     83 Granular.setTotalTime!(sim, 30.0)
     84 sim_init = deepcopy(sim)
     85 
     86 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
     87 Granular.setTimeStep!(sim, epsilon=0.07)
     88 tol = 0.1
     89 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
     90 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
     91 
     92 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
     93 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
     94 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
     95 @test E_kin_rot_init ≈ E_kin_rot_final
     96 
     97 
     98 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
     99 sim = deepcopy(sim_init)
    100 Granular.setTimeStep!(sim, epsilon=0.007)
    101 tol = 0.01
    102 @info "Relative tolerance: $(tol*100.)%"
    103 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
    104 
    105 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    106 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    107 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
    108 @test E_kin_rot_init ≈ E_kin_rot_final
    109 
    110 
    111 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    112 sim = deepcopy(sim_init)
    113 Granular.setTimeStep!(sim, epsilon=0.07)
    114 tol = 0.01
    115 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    116 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    117     verbose=verbose)
    118 
    119 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    120 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    121 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
    122 @test E_kin_rot_init ≈ E_kin_rot_final
    123 
    124 
    125 @info "## Contact-normal elasticity and tangential viscosity and friction"
    126 Granular.setTotalTime!(sim, 30.0)
    127 sim_init.grains[1].contact_viscosity_tangential = 1e6
    128 sim_init.grains[2].contact_viscosity_tangential = 1e6
    129 sim_init.grains[1].contact_dynamic_friction = 1e2  # no Coulomb sliding
    130 sim_init.grains[2].contact_dynamic_friction = 1e2  # no Coulomb sliding
    131 sim_init.grains[2].fixed = true
    132 sim = deepcopy(sim_init)
    133 
    134 
    135 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    136 Granular.setTimeStep!(sim, epsilon=0.07)
    137 tol = 0.1
    138 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    139 Granular.setOutputFileInterval!(sim, 1.0)
    140 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    141             verbose=verbose)
    142 
    143 @test sim.grains[1].ang_pos[3] < 0.
    144 @test sim.grains[1].ang_vel[3] < 0.
    145 @test sim.grains[2].ang_pos[3] ≈ 0.
    146 @test sim.grains[2].ang_vel[3] ≈ 0.
    147 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    148 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    149 println(E_kin_lin_init)
    150 println(E_kin_lin_final)
    151 println(E_kin_rot_init)
    152 println(E_kin_rot_final)
    153 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol
    154 
    155 @info "mu_d = 0."
    156 sim = deepcopy(sim_init)
    157 sim.grains[1].contact_dynamic_friction = 0.
    158 Granular.setTimeStep!(sim, epsilon=0.07)
    159 tol = 0.01
    160 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    161 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    162 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    163 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    164             verbose=verbose)
    165 
    166 @test sim.grains[1].ang_pos[3] ≈ 0.
    167 @test sim.grains[1].ang_vel[3] ≈ 0.
    168 @test sim.grains[2].ang_pos[3] ≈ 0.
    169 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    170 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    171 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
    172 @test E_kin_rot_init ≈ E_kin_rot_final
    173 
    174 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    175 sim = deepcopy(sim_init)
    176 Granular.setTimeStep!(sim, epsilon=0.007)
    177 tol = 0.1
    178 @info "Relative tolerance: $(tol*100.)%"
    179 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    180             verbose=verbose)
    181 
    182 @test sim.grains[1].ang_pos[3] < 0.
    183 @test sim.grains[1].ang_vel[3] < 0.
    184 @test sim.grains[2].ang_pos[3] ≈ 0.
    185 @test sim.grains[2].ang_vel[3] ≈ 0.
    186 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    187 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    188 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol
    189 
    190 
    191 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    192 sim = deepcopy(sim_init)
    193 Granular.setTimeStep!(sim, epsilon=0.07)
    194 tol = 0.09
    195 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    196 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    197             verbose=verbose)
    198 
    199 @test sim.grains[1].ang_pos[3] < 0.
    200 @test sim.grains[1].ang_vel[3] < 0.
    201 @test sim.grains[2].ang_pos[3] ≈ 0.
    202 @test sim.grains[2].ang_vel[3] ≈ 0.
    203 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    204 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    205 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol
    206 
    207 @info "# Ice floes free to move"
    208 
    209 sim = Granular.createSimulation(id="test")
    210 Granular.addGrainCylindrical!(sim, [0., 10.], 10., 1., verbose=verbose)
    211 Granular.addGrainCylindrical!(sim, [19.0, 0.], 10., 1., verbose=verbose)
    212 sim.grains[1].lin_vel[1] = 0.1
    213 sim.grains[1].contact_viscosity_tangential = 1e4
    214 sim.grains[2].contact_viscosity_tangential = 1e4
    215 
    216 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    217 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    218 
    219 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    220 # should become more correct
    221 
    222 Granular.setTotalTime!(sim, 30.0)
    223 sim_init = deepcopy(sim)
    224 
    225 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    226 Granular.setTimeStep!(sim, epsilon=0.07)
    227 tol = 0.1
    228 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    229 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    230             verbose=verbose)
    231 
    232 @test sim.grains[1].ang_pos[3] < 0.
    233 @test sim.grains[1].ang_vel[3] < 0.
    234 @test sim.grains[2].ang_pos[3] < 0.
    235 @test sim.grains[2].ang_vel[3] < 0.
    236 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    237 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    238 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    239 
    240 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    241 sim = deepcopy(sim_init)
    242 Granular.setTimeStep!(sim, epsilon=0.007)
    243 tol = 0.04
    244 @info "Relative tolerance: $(tol*100.)%"
    245 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    246             verbose=verbose)
    247 
    248 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    249 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    250 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    251 
    252 
    253 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    254 sim = deepcopy(sim_init)
    255 Granular.setTimeStep!(sim, epsilon=0.07)
    256 tol = 0.04
    257 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    258 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    259             verbose=verbose)
    260 
    261 @test sim.grains[1].ang_pos[3] < 0.
    262 @test sim.grains[1].ang_vel[3] < 0.
    263 @test sim.grains[2].ang_pos[3] < 0.
    264 @test sim.grains[2].ang_vel[3] < 0.
    265 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    266 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    267 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    268 
    269 
    270 @info "# Ice floes free to move, mirrored"
    271 
    272 sim = Granular.createSimulation(id="test")
    273 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
    274 Granular.addGrainCylindrical!(sim, [19.0, 10.], 10., 1., verbose=verbose)
    275 sim.grains[2].lin_vel[1] = -0.1
    276 sim.grains[1].contact_viscosity_tangential = 1e4
    277 sim.grains[2].contact_viscosity_tangential = 1e4
    278 
    279 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    280 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    281 
    282 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    283 # should become more correct
    284 
    285 Granular.setTotalTime!(sim, 30.0)
    286 sim_init = deepcopy(sim)
    287 
    288 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    289 Granular.setTimeStep!(sim, epsilon=0.07)
    290 tol = 0.1
    291 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    292 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    293             verbose=verbose)
    294 
    295 @test sim.grains[1].ang_pos[3] > 0.
    296 @test sim.grains[1].ang_vel[3] > 0.
    297 @test sim.grains[2].ang_pos[3] > 0.
    298 @test sim.grains[2].ang_vel[3] > 0.
    299 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    300 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    301 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    302 
    303 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    304 sim = deepcopy(sim_init)
    305 Granular.setTimeStep!(sim, epsilon=0.007)
    306 tol = 0.04
    307 @info "Relative tolerance: $(tol*100.)%"
    308 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    309             verbose=verbose)
    310 
    311 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    312 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    313 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    314 
    315 
    316 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    317 sim = deepcopy(sim_init)
    318 Granular.setTimeStep!(sim, epsilon=0.07)
    319 tol = 0.04
    320 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    321 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    322             verbose=verbose)
    323 
    324 @test sim.grains[1].ang_pos[3] > 0.
    325 @test sim.grains[1].ang_vel[3] > 0.
    326 @test sim.grains[2].ang_pos[3] > 0.
    327 @test sim.grains[2].ang_vel[3] > 0.
    328 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    329 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    330 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    331 
    332 
    333 @info "# Ice floes free to move, mirrored #2"
    334 
    335 sim = Granular.createSimulation(id="test")
    336 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
    337 Granular.addGrainCylindrical!(sim, [19.0, -10.], 10., 1., verbose=verbose)
    338 sim.grains[2].lin_vel[1] = -0.1
    339 
    340 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    341 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    342 
    343 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    344 # should become more correct
    345 
    346 Granular.setTotalTime!(sim, 30.0)
    347 sim_init = deepcopy(sim)
    348 
    349 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    350 Granular.setTimeStep!(sim, epsilon=0.07)
    351 tol = 0.1
    352 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    353 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    354             verbose=verbose)
    355 
    356 @test sim.grains[1].ang_pos[3] < 0.
    357 @test sim.grains[1].ang_vel[3] < 0.
    358 @test sim.grains[2].ang_pos[3] < 0.
    359 @test sim.grains[2].ang_vel[3] < 0.
    360 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    361 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    362 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    363 
    364 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    365 sim = deepcopy(sim_init)
    366 Granular.setTimeStep!(sim, epsilon=0.007)
    367 tol = 0.04
    368 @info "Relative tolerance: $(tol*100.)%"
    369 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    370             verbose=verbose)
    371 
    372 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    373 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    374 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    375 
    376 
    377 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    378 sim = deepcopy(sim_init)
    379 Granular.setTimeStep!(sim, epsilon=0.07)
    380 tol = 0.04
    381 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    382 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    383             verbose=verbose)
    384 
    385 @test sim.grains[1].ang_pos[3] < 0.
    386 @test sim.grains[1].ang_vel[3] < 0.
    387 @test sim.grains[2].ang_pos[3] < 0.
    388 @test sim.grains[2].ang_vel[3] < 0.
    389 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    390 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    391 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    392 
    393 
    394 @info "# Tangential elasticity, no tangential viscosity, no Coulomb slip"
    395 
    396 sim = Granular.createSimulation(id="test")
    397 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
    398 Granular.addGrainCylindrical!(sim, [19.0, -10.], 10., 1., verbose=verbose)
    399 sim.grains[2].lin_vel[1] = -0.1
    400 sim.grains[1].contact_dynamic_friction = 1e3  # disable Coulomb slip
    401 sim.grains[2].contact_dynamic_friction = 1e3  # disable Coulomb slip
    402 sim.grains[1].contact_viscosity_tangential = 0.  # disable tan. viscosity
    403 sim.grains[2].contact_viscosity_tangential = 0.  # disable tan. viscosity
    404 sim.grains[1].contact_stiffness_tangential = 
    405     sim.grains[1].contact_stiffness_normal  # enable tangential elasticity
    406 sim.grains[2].contact_stiffness_tangential = 
    407     sim.grains[2].contact_stiffness_normal  # enable tangential elasticity
    408 
    409 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    410 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    411 
    412 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    413 # should become more correct
    414 
    415 Granular.setTotalTime!(sim, 30.0)
    416 sim_init = deepcopy(sim)
    417 
    418 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    419 Granular.setTimeStep!(sim, epsilon=0.07)
    420 tol = 0.1
    421 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    422 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    423             verbose=verbose)
    424 
    425 @test sim.grains[1].ang_pos[3] < 0.
    426 @test sim.grains[1].ang_vel[3] < 0.
    427 @test sim.grains[2].ang_pos[3] < 0.
    428 @test sim.grains[2].ang_vel[3] < 0.
    429 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    430 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    431 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    432 
    433 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    434 sim = deepcopy(sim_init)
    435 Granular.setTimeStep!(sim, epsilon=0.007)
    436 tol = 0.04
    437 @info "Relative tolerance: $(tol*100.)%"
    438 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    439             verbose=verbose)
    440 
    441 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    442 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    443 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    444 
    445 
    446 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    447 sim = deepcopy(sim_init)
    448 Granular.setTimeStep!(sim, epsilon=0.07)
    449 tol = 0.04
    450 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    451 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    452             verbose=verbose)
    453 
    454 @test sim.grains[1].ang_pos[3] < 0.
    455 @test sim.grains[1].ang_vel[3] < 0.
    456 @test sim.grains[2].ang_pos[3] < 0.
    457 @test sim.grains[2].ang_vel[3] < 0.
    458 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    459 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    460 @test E_kin_lin_init+E_kin_rot_init ≈ E_kin_lin_final+E_kin_rot_final atol=E_kin_lin_init*tol 
    461 
    462 
    463 @info "# Tangential elasticity, no tangential viscosity, Coulomb slip"
    464 
    465 sim = Granular.createSimulation(id="test")
    466 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
    467 Granular.addGrainCylindrical!(sim, [19.0, -10.], 10., 1., verbose=verbose)
    468 sim.grains[2].lin_vel[1] = -0.1
    469 sim.grains[1].contact_dynamic_friction = 0.1  # enable Coulomb slip
    470 sim.grains[2].contact_dynamic_friction = 0.1  # enable Coulomb slip
    471 sim.grains[1].contact_viscosity_tangential = 0.  # disable tan. viscosity
    472 sim.grains[2].contact_viscosity_tangential = 0.  # disable tan. viscosity
    473 sim.grains[1].contact_stiffness_tangential = 
    474     sim.grains[1].contact_stiffness_normal  # enable tangential elasticity
    475 sim.grains[2].contact_stiffness_tangential = 
    476     sim.grains[2].contact_stiffness_normal  # enable tangential elasticity
    477 
    478 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    479 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    480 
    481 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    482 # should become more correct
    483 
    484 Granular.setTotalTime!(sim, 30.0)
    485 sim_init = deepcopy(sim)
    486 
    487 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    488 sim = deepcopy(sim_init)
    489 Granular.setTimeStep!(sim, epsilon=0.007)
    490 tol = 0.02
    491 @info "Relative tolerance: $(tol*100.)%"
    492 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    493             verbose=verbose)
    494 
    495 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    496 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    497 @test E_kin_lin_init+E_kin_rot_init > E_kin_lin_final+E_kin_rot_final
    498 
    499 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    500 sim = deepcopy(sim_init)
    501 Granular.setTimeStep!(sim, epsilon=0.07)
    502 tol = 0.03
    503 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    504 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    505             verbose=verbose)
    506 
    507 @test sim.grains[1].ang_pos[3] < 0.
    508 @test sim.grains[1].ang_vel[3] < 0.
    509 @test sim.grains[2].ang_pos[3] < 0.
    510 @test sim.grains[2].ang_vel[3] < 0.
    511 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    512 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    513 @test E_kin_lin_init+E_kin_rot_init > E_kin_lin_final+E_kin_rot_final
    514 
    515 
    516 @info "# Tangential elasticity, tangential viscosity, no Coulomb slip"
    517 
    518 sim = Granular.createSimulation(id="test")
    519 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
    520 Granular.addGrainCylindrical!(sim, [19.0, -10.], 10., 1., verbose=verbose)
    521 sim.grains[2].lin_vel[1] = -0.1
    522 sim.grains[1].contact_dynamic_friction = 1e3  # disable Coulomb slip
    523 sim.grains[2].contact_dynamic_friction = 1e3  # disable Coulomb slip
    524 sim.grains[1].contact_viscosity_tangential = 1e4  # enable tan. viscosity
    525 sim.grains[2].contact_viscosity_tangential = 1e4  # enable tan. viscosity
    526 sim.grains[1].contact_stiffness_tangential = 
    527     sim.grains[1].contact_stiffness_normal  # enable tangential elasticity
    528 sim.grains[2].contact_stiffness_tangential = 
    529     sim.grains[2].contact_stiffness_normal  # enable tangential elasticity
    530 
    531 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    532 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    533 
    534 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    535 # should become more correct
    536 
    537 Granular.setTotalTime!(sim, 30.0)
    538 sim_init = deepcopy(sim)
    539 
    540 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    541 sim = deepcopy(sim_init)
    542 Granular.setTimeStep!(sim, epsilon=0.007)
    543 tol = 0.02
    544 @info "Relative tolerance: $(tol*100.)%"
    545 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    546             verbose=verbose)
    547 
    548 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    549 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    550 @test E_kin_lin_init+E_kin_rot_init > E_kin_lin_final+E_kin_rot_final
    551 
    552 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    553 sim = deepcopy(sim_init)
    554 Granular.setTimeStep!(sim, epsilon=0.07)
    555 tol = 0.03
    556 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    557 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    558             verbose=verbose)
    559 
    560 @test sim.grains[1].ang_pos[3] < 0.
    561 @test sim.grains[1].ang_vel[3] < 0.
    562 @test sim.grains[2].ang_pos[3] < 0.
    563 @test sim.grains[2].ang_vel[3] < 0.
    564 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    565 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    566 @test E_kin_lin_init+E_kin_rot_init > E_kin_lin_final+E_kin_rot_final
    567 
    568 
    569 @info "# Tangential elasticity, tangential viscosity, Coulomb slip"
    570 
    571 sim = Granular.createSimulation(id="test")
    572 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
    573 Granular.addGrainCylindrical!(sim, [19.0, -10.], 10., 1., verbose=verbose)
    574 sim.grains[2].lin_vel[1] = -0.1
    575 sim.grains[1].contact_dynamic_friction = 0.1  # enable Coulomb slip
    576 sim.grains[2].contact_dynamic_friction = 0.1  # enable Coulomb slip
    577 sim.grains[1].contact_viscosity_tangential = 1e4  # enable tan. viscosity
    578 sim.grains[2].contact_viscosity_tangential = 1e4  # enable tan. viscosity
    579 sim.grains[1].contact_stiffness_tangential = 
    580     sim.grains[1].contact_stiffness_normal  # enable tangential elasticity
    581 sim.grains[2].contact_stiffness_tangential = 
    582     sim.grains[2].contact_stiffness_normal  # enable tangential elasticity
    583 
    584 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
    585 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
    586 
    587 # With decreasing timestep (epsilon towards 0), the explicit integration scheme 
    588 # should become more correct
    589 
    590 Granular.setTotalTime!(sim, 30.0)
    591 sim_init = deepcopy(sim)
    592 
    593 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
    594 sim = deepcopy(sim_init)
    595 Granular.setTimeStep!(sim, epsilon=0.007)
    596 tol = 0.02
    597 @info "Relative tolerance: $(tol*100.)%"
    598 Granular.run!(sim, temporal_integration_method="Two-term Taylor",
    599             verbose=verbose)
    600 
    601 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    602 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    603 @test E_kin_lin_init+E_kin_rot_init > E_kin_lin_final+E_kin_rot_final
    604 
    605 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
    606 sim = deepcopy(sim_init)
    607 Granular.setTimeStep!(sim, epsilon=0.07)
    608 tol = 0.03
    609 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
    610 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
    611             verbose=verbose)
    612 
    613 @test sim.grains[1].ang_pos[3] < 0.
    614 @test sim.grains[1].ang_vel[3] < 0.
    615 @test sim.grains[2].ang_pos[3] < 0.
    616 @test sim.grains[2].ang_vel[3] < 0.
    617 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
    618 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
    619 @test E_kin_lin_init+E_kin_rot_init > E_kin_lin_final+E_kin_rot_final