tcollision-2floes-oblique.jl - Granular.jl - Julia package for granular dynamics simulation
 (HTM) git clone git://src.adamsgaard.dk/Granular.jl
 (DIR) Log
 (DIR) Files
 (DIR) Refs
 (DIR) README
 (DIR) LICENSE
       ---
       tcollision-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