tcollision-2floes-normal.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-normal.jl (13115B)
---
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 "# One ice floe fixed"
9 sim = Granular.createSimulation(id="test")
10 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
11 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
12 sim.grains[1].lin_vel[1] = 0.1
13 sim.grains[2].fixed = true
14
15 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
16 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
17
18 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
19 # should become more correct
20
21 Granular.setTotalTime!(sim, 10.0)
22 sim_init = deepcopy(sim)
23
24 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
25 Granular.setTimeStep!(sim, epsilon=0.07)
26 tol = 0.2
27 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
28 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
29
30 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
31 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
32 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
33 @test E_kin_rot_init ≈ E_kin_rot_final
34
35
36 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
37 sim = deepcopy(sim_init)
38 Granular.setTimeStep!(sim, epsilon=0.007)
39 tol = 0.02
40 @info "Relative tolerance: $(tol*100.)%"
41 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
42
43 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
44 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
45 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
46 @test E_kin_rot_init ≈ E_kin_rot_final
47
48
49 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
50 sim = deepcopy(sim_init)
51 Granular.setTimeStep!(sim, epsilon=0.07)
52 tol = 0.01
53 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
54 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
55 verbose=verbose)
56
57 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
58 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
59 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
60 @test E_kin_rot_init ≈ E_kin_rot_final
61
62
63 @info "# Ice floes free to move"
64
65 sim = Granular.createSimulation(id="test")
66 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
67 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
68 sim.grains[1].lin_vel[1] = 0.1
69
70 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
71 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
72
73 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
74 # should become more correct
75
76 Granular.setTotalTime!(sim, 10.0)
77 sim_init = deepcopy(sim)
78
79 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
80 Granular.setTimeStep!(sim, epsilon=0.07)
81 tol = 0.2
82 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
83 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
84
85 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
86 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
87 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
88 @test E_kin_rot_init ≈ E_kin_rot_final
89
90
91 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
92 sim = deepcopy(sim_init)
93 Granular.setTimeStep!(sim, epsilon=0.007)
94 tol = 0.02
95 @info "Relative tolerance: $(tol*100.)%"
96 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
97
98 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
99 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
100 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
101 @test E_kin_rot_init ≈ E_kin_rot_final
102
103
104 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
105 sim = deepcopy(sim_init)
106 Granular.setTimeStep!(sim, epsilon=0.07)
107 tol = 0.01
108 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
109 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
110 verbose=verbose)
111
112 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
113 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
114 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
115 @test E_kin_rot_init ≈ E_kin_rot_final
116
117
118 @info "# Adding contact-normal viscosity"
119 @info "# One ice floe fixed"
120 sim = Granular.createSimulation(id="test")
121 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
122 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
123 sim.grains[1].lin_vel[1] = 0.1
124 sim.grains[1].contact_viscosity_normal = 1e4
125 sim.grains[2].contact_viscosity_normal = 1e4
126 sim.grains[2].fixed = true
127
128 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
129 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
130
131 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
132 # should become more correct
133
134 Granular.setTotalTime!(sim, 10.0)
135 sim_init = deepcopy(sim)
136
137
138 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
139 sim = deepcopy(sim_init)
140 Granular.setTimeStep!(sim, epsilon=0.007)
141 tol = 0.02
142 @info "Relative tolerance: $(tol*100.)%"
143 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
144
145 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
146 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
147 @test E_kin_lin_init > E_kin_lin_final
148 @test E_kin_rot_init ≈ E_kin_rot_final
149
150
151 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
152 sim = deepcopy(sim_init)
153 Granular.setTimeStep!(sim, epsilon=0.07)
154 tol = 0.01
155 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
156 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
157 verbose=verbose)
158
159 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
160 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
161 @test E_kin_lin_init > E_kin_lin_final
162 @test E_kin_rot_init ≈ E_kin_rot_final
163
164
165 @info "# Ice floes free to move"
166
167 sim = Granular.createSimulation(id="test")
168 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
169 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
170 sim.grains[1].lin_vel[1] = 0.1
171 sim.grains[1].contact_viscosity_normal = 1e4
172 sim.grains[2].contact_viscosity_normal = 1e4
173
174 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
175 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
176
177 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
178 # should become more correct
179
180 Granular.setTotalTime!(sim, 10.0)
181 sim_init = deepcopy(sim)
182
183 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
184 sim = deepcopy(sim_init)
185 Granular.setTimeStep!(sim, epsilon=0.007)
186 tol = 0.02
187 @info "Relative tolerance: $(tol*100.)%"
188 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
189
190 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
191 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
192 @test E_kin_lin_init > E_kin_lin_final
193 @test E_kin_rot_init ≈ E_kin_rot_final
194
195
196 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
197 sim = deepcopy(sim_init)
198 Granular.setTimeStep!(sim, epsilon=0.07)
199 tol = 0.01
200 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
201 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
202 verbose=verbose)
203
204 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
205 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
206 @test E_kin_lin_init > E_kin_lin_final
207 @test E_kin_rot_init ≈ E_kin_rot_final
208
209
210 @info "# Testing allow_*_acc for fixed grains"
211 sim = Granular.createSimulation(id="test")
212 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
213 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
214 sim.grains[1].lin_vel[1] = 0.1
215 sim.grains[2].fixed = true
216
217 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
218 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
219 grain2_pos_init = sim.grains[2].lin_pos
220
221 Granular.setTotalTime!(sim, 10.0)
222 Granular.setTimeStep!(sim, epsilon=0.07)
223 sim_init = deepcopy(sim)
224 sim.grains[2].allow_y_acc = true # should not influence result
225
226 @info "Two-term Taylor scheme: allow_y_acc"
227 sim = deepcopy(sim_init)
228 sim.grains[2].allow_y_acc = true # should not influence result
229 tol = 0.2
230 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
231
232 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
233 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
234 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
235 @test E_kin_rot_init ≈ E_kin_rot_final
236 @test sim.grains[2].lin_pos ≈ grain2_pos_init
237
238 @info "Two-term Taylor scheme: allow_x_acc"
239 sim = deepcopy(sim_init)
240 sim.grains[2].allow_x_acc = true # should influence result
241 tol = 0.2
242 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
243
244 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
245 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
246 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
247 @test E_kin_rot_init ≈ E_kin_rot_final
248 @test sim.grains[2].lin_pos[1] > grain2_pos_init[1]
249
250 @info "Three-term Taylor scheme: allow_y_acc"
251 sim = deepcopy(sim_init)
252 tol = 0.02
253 sim.grains[2].allow_y_acc = true # should influence result
254 Granular.run!(sim, temporal_integration_method="Three-term Taylor", verbose=verbose)
255
256 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
257 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
258 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
259 @test E_kin_rot_init ≈ E_kin_rot_final
260 @test sim.grains[2].lin_pos ≈ grain2_pos_init
261
262 @info "Three-term Taylor scheme: allow_x_acc"
263 sim = deepcopy(sim_init)
264 tol = 0.02
265 sim.grains[2].allow_x_acc = true # should influence result
266 Granular.run!(sim, temporal_integration_method="Three-term Taylor", verbose=verbose)
267
268 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
269 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
270 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
271 @test E_kin_rot_init ≈ E_kin_rot_final
272 @test sim.grains[2].lin_pos[1] > grain2_pos_init[1]
273
274 #=
275 @info "# Test stability under collision with fixed particles different allow_*_acc"
276 r = 10.
277 i = 1
278 for tensile_strength in [0.0, 200e3]
279 for angle in range(0, 2π, 7)
280 for allow_x_acc in [false, true]
281 for allow_y_acc in [false, true]
282 @info "Test $i"
283 @info "Contact angle: $angle rad"
284 @info "allow_x_acc = $allow_x_acc"
285 @info "allow_y_acc = $allow_y_acc"
286 @info "tensile_strength = $tensile_strength Pa"
287
288 sim = Granular.createSimulation()
289 sim.id = "test-$i-$allow_x_acc-$allow_y_acc-C=$tensile_strength"
290 Granular.addGrainCylindrical!(sim, [0., 0.], r, 1., verbose=verbose)
291 Granular.addGrainCylindrical!(sim, [2.0*r*cos(angle), 2.0*r*sin(angle)],
292 r, 1., verbose=verbose)
293 sim.grains[1].lin_vel = r/10.0 .* [cos(angle), sin(angle)]
294
295 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
296 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
297 grain1_pos_init = sim.grains[1].lin_pos
298 grain2_pos_init = sim.grains[2].lin_pos
299
300 sim.grains[1].fixed = true
301 sim.grains[2].fixed = true
302
303 sim.grains[1].allow_x_acc = allow_x_acc
304 sim.grains[2].allow_x_acc = allow_x_acc
305 sim.grains[1].allow_y_acc = allow_y_acc
306 sim.grains[2].allow_y_acc = allow_y_acc
307
308 sim.grains[1].tensile_strength = tensile_strength
309 sim.grains[2].tensile_strength = tensile_strength
310
311 Granular.setTotalTime!(sim, 20.0)
312 Granular.setTimeStep!(sim, epsilon=0.07)
313 sim_init = deepcopy(sim)
314
315 @info "TY3"
316 sim = deepcopy(sim_init)
317 tol = 0.02
318 Granular.setOutputFileInterval!(sim, 1.0)
319 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
320 verbose=verbose)
321 Granular.render(sim)
322 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
323 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
324 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
325
326 i += 1
327 end
328 end
329 end
330 end
331 =#