tcollision-5floes-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-5floes-normal.jl (10374B)
---
1 #!/usr/bin/env julia
2 using LinearAlgebra
3
4 # Check for conservation of kinetic energy (=momentum) during a normal collision
5 # between two ice cylindrical grains
6
7 verbose=false
8
9 @info "# One ice floe fixed"
10 sim = Granular.createSimulation(id="test")
11 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
12 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
13 Granular.addGrainCylindrical!(sim, [40.05, 0.], 10., 1., verbose=verbose)
14 Granular.addGrainCylindrical!(sim, [60.05, 0.], 10., 1., verbose=verbose)
15 Granular.addGrainCylindrical!(sim, [80.05, 0.], 10., 1., verbose=verbose)
16 sim.grains[1].lin_vel[1] = 0.1
17 sim.grains[2].fixed = true
18 sim.grains[3].fixed = true
19 sim.grains[4].fixed = true
20 sim.grains[5].fixed = true
21
22 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
23 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
24
25 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
26 # should become more correct
27
28 Granular.setTotalTime!(sim, 10.0)
29 sim_init = deepcopy(sim)
30
31 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
32 Granular.setTimeStep!(sim, epsilon=0.07)
33 tol = 0.2
34 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
35 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
36
37 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
38 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
39 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
40 @test E_kin_rot_init ≈ E_kin_rot_final
41 @test 0. < norm(sim.grains[1].lin_vel)
42 for i=2:5
43 @info "testing ice floe $i"
44 @test 0. ≈ norm(sim.grains[i].lin_vel)
45 end
46
47
48 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
49 sim = deepcopy(sim_init)
50 Granular.setTimeStep!(sim, epsilon=0.007)
51 tol = 0.02
52 @info "Relative tolerance: $(tol*100.)%"
53 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
54
55 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
56 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
57 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
58 @test E_kin_rot_init ≈ E_kin_rot_final
59 @test 0. < norm(sim.grains[1].lin_vel)
60 for i=2:5
61 @info "testing ice floe $i"
62 @test 0. ≈ norm(sim.grains[i].lin_vel)
63 end
64
65
66 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
67 sim = deepcopy(sim_init)
68 Granular.setTimeStep!(sim, epsilon=0.07)
69 tol = 0.01
70 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
71 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
72 verbose=verbose)
73
74 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
75 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
76 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
77 @test E_kin_rot_init ≈ E_kin_rot_final
78 @test 0. < norm(sim.grains[1].lin_vel)
79 for i=2:5
80 @info "testing ice floe $i"
81 @test 0. ≈ norm(sim.grains[i].lin_vel)
82 end
83
84
85 @info "# Ice floes free to move"
86
87 sim = Granular.createSimulation(id="test")
88 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
89 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
90 Granular.addGrainCylindrical!(sim, [40.05, 0.], 10., 1., verbose=verbose)
91 Granular.addGrainCylindrical!(sim, [60.05, 0.], 10., 1., verbose=verbose)
92 Granular.addGrainCylindrical!(sim, [80.05, 0.], 10., 1., verbose=verbose)
93 sim.grains[1].lin_vel[1] = 0.1
94
95 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
96 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
97
98 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
99 # should become more correct
100
101 Granular.setTotalTime!(sim, 40.0)
102 sim_init = deepcopy(sim)
103
104 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
105 Granular.setTimeStep!(sim, epsilon=0.07)
106 tol = 0.2
107 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
108 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
109
110 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
111 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
112 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
113 @test E_kin_rot_init ≈ E_kin_rot_final
114 for i=1:5
115 @info "testing ice floe $i"
116 @test 0. < norm(sim.grains[i].lin_vel)
117 end
118
119
120 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
121 sim = deepcopy(sim_init)
122 Granular.setTimeStep!(sim, epsilon=0.007)
123 tol = 0.02
124 @info "Relative tolerance: $(tol*100.)%"
125 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
126
127 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
128 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
129 @test E_kin_lin_init ≈ E_kin_lin_final atol=E_kin_lin_init*tol
130 @test E_kin_rot_init ≈ E_kin_rot_final
131 for i=1:5
132 @info "testing ice floe $i"
133 @test 0. < norm(sim.grains[i].lin_vel)
134 end
135
136
137 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
138 sim = deepcopy(sim_init)
139 Granular.setTimeStep!(sim, epsilon=0.07)
140 tol = 0.01
141 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
142 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
143 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 atol=E_kin_lin_init*tol
148 @test E_kin_rot_init ≈ E_kin_rot_final
149 for i=1:5
150 @info "testing ice floe $i"
151 @test 0. < norm(sim.grains[i].lin_vel)
152 end
153
154
155 @info "# Adding contact-normal viscosity"
156 @info "# One ice floe fixed"
157 sim = Granular.createSimulation(id="test")
158 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
159 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
160 Granular.addGrainCylindrical!(sim, [40.05, 0.], 10., 1., verbose=verbose)
161 Granular.addGrainCylindrical!(sim, [60.05, 0.], 10., 1., verbose=verbose)
162 Granular.addGrainCylindrical!(sim, [80.05, 0.], 10., 1., verbose=verbose)
163 sim.grains[1].lin_vel[1] = 0.1
164 sim.grains[1].contact_viscosity_normal = 1e4
165 sim.grains[2].contact_viscosity_normal = 1e4
166 sim.grains[3].contact_viscosity_normal = 1e4
167 sim.grains[4].contact_viscosity_normal = 1e4
168 sim.grains[5].contact_viscosity_normal = 1e4
169 sim.grains[2].fixed = true
170 sim.grains[3].fixed = true
171 sim.grains[4].fixed = true
172 sim.grains[5].fixed = true
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
184 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
185 sim = deepcopy(sim_init)
186 Granular.setTimeStep!(sim, epsilon=0.007)
187 tol = 0.02
188 @info "Relative tolerance: $(tol*100.)%"
189 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
190
191 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
192 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
193 @test E_kin_lin_init > E_kin_lin_final
194 @test E_kin_rot_init ≈ E_kin_rot_final
195 @test 0. < norm(sim.grains[1].lin_vel)
196 for i=2:5
197 @info "testing ice floe $i"
198 @test 0. ≈ norm(sim.grains[i].lin_vel)
199 end
200
201
202 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
203 sim = deepcopy(sim_init)
204 Granular.setTimeStep!(sim, epsilon=0.07)
205 tol = 0.01
206 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
207 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
208 verbose=verbose)
209
210 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
211 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
212 @test E_kin_lin_init > E_kin_lin_final
213 @test E_kin_rot_init ≈ E_kin_rot_final
214 @test 0. < norm(sim.grains[1].lin_vel)
215 for i=2:5
216 @info "testing ice floe $i"
217 @test 0. ≈ norm(sim.grains[i].lin_vel)
218 end
219
220
221 @info "# Ice floes free to move"
222
223 sim = Granular.createSimulation(id="test")
224 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., verbose=verbose)
225 Granular.addGrainCylindrical!(sim, [20.05, 0.], 10., 1., verbose=verbose)
226 Granular.addGrainCylindrical!(sim, [40.05, 0.], 10., 1., verbose=verbose)
227 Granular.addGrainCylindrical!(sim, [60.05, 0.], 10., 1., verbose=verbose)
228 Granular.addGrainCylindrical!(sim, [80.05, 0.], 10., 1., verbose=verbose)
229 sim.grains[1].lin_vel[1] = 0.1
230 sim.grains[1].contact_viscosity_normal = 1e4
231 sim.grains[2].contact_viscosity_normal = 1e4
232 sim.grains[3].contact_viscosity_normal = 1e4
233 sim.grains[4].contact_viscosity_normal = 1e4
234 sim.grains[5].contact_viscosity_normal = 1e4
235
236 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim)
237 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim)
238
239 # With decreasing timestep (epsilon towards 0), the explicit integration scheme
240 # should become more correct
241
242 Granular.setTotalTime!(sim, 10.0)
243 sim_init = deepcopy(sim)
244
245 @info "Testing kinetic energy conservation with Two-term Taylor scheme"
246 sim = deepcopy(sim_init)
247 Granular.setTimeStep!(sim, epsilon=0.007)
248 tol = 0.02
249 @info "Relative tolerance: $(tol*100.)%"
250 Granular.run!(sim, temporal_integration_method="Two-term Taylor", verbose=verbose)
251
252 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
253 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
254 @test E_kin_lin_init > E_kin_lin_final
255 @test E_kin_rot_init ≈ E_kin_rot_final
256 for i=1:5
257 @info "testing ice floe $i"
258 @test 0. < norm(sim.grains[i].lin_vel)
259 end
260
261
262 @info "Testing kinetic energy conservation with Three-term Taylor scheme"
263 sim = deepcopy(sim_init)
264 Granular.setTimeStep!(sim, epsilon=0.07)
265 tol = 0.01
266 @info "Relative tolerance: $(tol*100.)% with time step: $(sim.time_step)"
267 Granular.run!(sim, temporal_integration_method="Three-term Taylor",
268 verbose=verbose)
269
270 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim)
271 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim)
272 @test E_kin_lin_init > E_kin_lin_final
273 @test E_kin_rot_init ≈ E_kin_rot_final
274 for i=1:5
275 @info "testing ice floe $i"
276 @test 0. < norm(sim.grains[i].lin_vel)
277 end