 |
tcohesion.jl - Granular.jl - Julia package for granular dynamics simulation |
 |
git clone git://src.adamsgaard.dk/Granular.jl |
 |
Log |
|
Files |
|
Refs |
|
README |
|
LICENSE |
|
--- |
|
tcohesion.jl (8100B) |
|
--- |
|
1 #!/usr/bin/env julia |
|
2 using Test |
|
3 import Granular |
|
4 |
|
5 # Check for conservation of kinetic energy (=momentum) during a normal c… |
|
6 # between two ice cylindrical grains |
|
7 |
|
8 verbose=false |
|
9 |
|
10 sim_init = Granular.createSimulation() |
|
11 Granular.addGrainCylindrical!(sim_init, [0., 0.], 10., 1.) |
|
12 Granular.addGrainCylindrical!(sim_init, [18., 0.], 10., 1.) |
|
13 sim_init.grains[1].youngs_modulus = 1e-5 # repulsion is negligible |
|
14 sim_init.grains[2].youngs_modulus = 1e-5 # repulsion is negligible |
|
15 Granular.setTimeStep!(sim_init, verbose=verbose) |
|
16 |
|
17 @info "# Check contact age scheme" |
|
18 sim = deepcopy(sim_init) |
|
19 Granular.setTotalTime!(sim, 10.) |
|
20 sim.time_step = 1. |
|
21 Granular.run!(sim, verbose=verbose) |
|
22 Granular.removeSimulationFiles(sim) |
|
23 @test sim.grains[1].contact_age[1] ≈ sim.time |
|
24 |
|
25 @info "# Check if bonds add tensile strength" |
|
26 sim = Granular.createSimulation(id="cohesion") |
|
27 Granular.addGrainCylindrical!(sim, [0., 0.], 10., 1., tensile_strength=5… |
|
28 Granular.addGrainCylindrical!(sim, [20.1, 0.], 10., 1., tensile_strength… |
|
29 sim.grains[1].lin_vel[1] = 0.1 |
|
30 Granular.setTimeStep!(sim) |
|
31 Granular.setTotalTime!(sim, 10.) |
|
32 Granular.run!(sim, verbose=verbose) |
|
33 Granular.removeSimulationFiles(sim) |
|
34 @test sim.grains[1].lin_vel[1] > 0. |
|
35 @test sim.grains[1].lin_vel[2] ≈ 0. |
|
36 @test sim.grains[2].lin_vel[1] > 0. |
|
37 @test sim.grains[2].lin_vel[2] ≈ 0. |
|
38 @test sim.grains[1].ang_vel ≈ zeros(3) |
|
39 @test sim.grains[2].ang_vel ≈ zeros(3) |
|
40 |
|
41 @info "# Add shear strength and test bending resistance (one grain rotat… |
|
42 sim = Granular.createSimulation(id="cohesion") |
|
43 Granular.addGrainCylindrical!(sim, [0., 0.], 10.1, 1., tensile_strength=… |
|
44 shear_strength=500e3) |
|
45 Granular.addGrainCylindrical!(sim, [20., 0.], 10., 1., tensile_strength=… |
|
46 shear_strength=500e3) |
|
47 sim.grains[1].ang_vel[3] = 0.1 |
|
48 Granular.findContacts!(sim) # make sure contact is registered |
|
49 sim.grains[1].contact_radius=10.0 # decrease radius so there isn't compr… |
|
50 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
51 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim) |
|
52 Granular.setTimeStep!(sim) |
|
53 Granular.setTotalTime!(sim, 5.) |
|
54 #Granular.setOutputFileInterval!(sim, 0.1) |
|
55 Granular.run!(sim, verbose=verbose) |
|
56 Granular.removeSimulationFiles(sim) |
|
57 @test sim.grains[1].lin_vel[1] ≈ 0. |
|
58 @test sim.grains[1].lin_vel[2] ≈ 0. |
|
59 @test sim.grains[2].lin_vel[1] ≈ 0. |
|
60 @test sim.grains[2].lin_vel[2] ≈ 0. |
|
61 @test sim.grains[1].ang_vel[3] != 0. |
|
62 @test sim.grains[2].ang_vel[3] != 0. |
|
63 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
64 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim) |
|
65 E_therm_final = Granular.totalGrainThermalEnergy(sim) |
|
66 @test E_kin_lin_init ≈ E_kin_lin_final |
|
67 @test E_kin_rot_init > E_kin_rot_final + E_therm_final |
|
68 |
|
69 @info "# Add shear strength and test bending resistance (one grain rotat… |
|
70 sim = Granular.createSimulation(id="cohesion") |
|
71 Granular.addGrainCylindrical!(sim, [0., 0.], 10.1, 1., tensile_strength=… |
|
72 shear_strength=500e3) |
|
73 Granular.addGrainCylindrical!(sim, [20., 0.], 10., 1., tensile_strength=… |
|
74 shear_strength=500e3) |
|
75 sim.grains[2].ang_vel[3] = 0.1 |
|
76 Granular.findContacts!(sim) # make sure contact is registered |
|
77 sim.grains[1].contact_radius=10.0 # decrease radius so there isn't compr… |
|
78 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
79 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim) |
|
80 Granular.setTimeStep!(sim) |
|
81 Granular.setTotalTime!(sim, 5.) |
|
82 #Granular.setOutputFileInterval!(sim, 0.1) |
|
83 Granular.run!(sim, verbose=verbose) |
|
84 Granular.removeSimulationFiles(sim) |
|
85 @test sim.grains[1].lin_vel[1] ≈ 0. |
|
86 @test sim.grains[1].lin_vel[2] ≈ 0. |
|
87 @test sim.grains[2].lin_vel[1] ≈ 0. |
|
88 @test sim.grains[2].lin_vel[2] ≈ 0. |
|
89 @test sim.grains[1].ang_vel[3] != 0. |
|
90 @test sim.grains[2].ang_vel[3] != 0. |
|
91 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
92 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim) |
|
93 E_therm_final = Granular.totalGrainThermalEnergy(sim) |
|
94 @test E_kin_lin_init ≈ E_kin_lin_final |
|
95 @test E_kin_rot_init > E_kin_rot_final + E_therm_final |
|
96 |
|
97 @info "# Add shear strength and test bending resistance (both grains rot… |
|
98 sim = Granular.createSimulation(id="cohesion") |
|
99 Granular.addGrainCylindrical!(sim, [0., 0.], 10.0000001, 1., tensile_str… |
|
100 shear_strength=500e3) |
|
101 Granular.addGrainCylindrical!(sim, [20., 0.], 10., 1., tensile_strength=… |
|
102 shear_strength=500e3) |
|
103 sim.grains[1].ang_vel[3] = 0.1 |
|
104 sim.grains[2].ang_vel[3] = -0.1 |
|
105 Granular.findContacts!(sim) # make sure contact is registered |
|
106 sim.grains[1].contact_radius=10.0 # decrease radius so there isn't compr… |
|
107 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
108 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim) |
|
109 Granular.setTimeStep!(sim) |
|
110 Granular.setTotalTime!(sim, 5.) |
|
111 #Granular.setOutputFileInterval!(sim, 0.1) |
|
112 Granular.run!(sim, verbose=verbose) |
|
113 Granular.removeSimulationFiles(sim) |
|
114 @test sim.grains[1].lin_vel[1] ≈ 0. |
|
115 @test sim.grains[1].lin_vel[2] ≈ 0. |
|
116 @test sim.grains[2].lin_vel[1] ≈ 0. |
|
117 @test sim.grains[2].lin_vel[2] ≈ 0. |
|
118 @test sim.grains[1].ang_vel[3] != 0. |
|
119 @test sim.grains[2].ang_vel[3] != 0. |
|
120 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
121 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim) |
|
122 E_therm_final = Granular.totalGrainThermalEnergy(sim) |
|
123 @test E_kin_lin_init ≈ E_kin_lin_final |
|
124 @test E_kin_rot_init > E_kin_rot_final + E_therm_final |
|
125 |
|
126 @info "# Break bond through bending I" |
|
127 sim = Granular.createSimulation(id="cohesion") |
|
128 Granular.addGrainCylindrical!(sim, [0., 0.], 10.0000001, 1., tensile_str… |
|
129 shear_strength=500e3) |
|
130 Granular.addGrainCylindrical!(sim, [20., 0.], 10., 1., tensile_strength=… |
|
131 shear_strength=500e3) |
|
132 sim.grains[1].ang_vel[3] = 100 |
|
133 sim.grains[2].ang_vel[3] = -100 |
|
134 Granular.findContacts!(sim) # make sure contact is registered |
|
135 sim.grains[1].contact_radius=10.0 # decrease radius so there isn't compr… |
|
136 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
137 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim) |
|
138 Granular.setTimeStep!(sim) |
|
139 Granular.setTotalTime!(sim, 5.) |
|
140 #Granular.setOutputFileInterval!(sim, 0.1) |
|
141 Granular.run!(sim, verbose=verbose) |
|
142 Granular.removeSimulationFiles(sim) |
|
143 @test sim.grains[1].lin_vel[1] ≈ 0. |
|
144 @test sim.grains[1].lin_vel[2] ≈ 0. |
|
145 @test sim.grains[2].lin_vel[1] ≈ 0. |
|
146 @test sim.grains[2].lin_vel[2] ≈ 0. |
|
147 @test sim.grains[1].ang_vel[3] != 0. |
|
148 @test sim.grains[2].ang_vel[3] != 0. |
|
149 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
150 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim) |
|
151 E_therm_final = Granular.totalGrainThermalEnergy(sim) |
|
152 @test E_kin_lin_init ≈ E_kin_lin_final |
|
153 @test sim.grains[1].n_contacts == 0 |
|
154 @test sim.grains[2].n_contacts == 0 |
|
155 |
|
156 @info "# Break bond through bending II" |
|
157 sim = Granular.createSimulation(id="cohesion") |
|
158 Granular.addGrainCylindrical!(sim, [0., 0.], 10.1, 1., tensile_strength=… |
|
159 shear_strength=50e3) |
|
160 Granular.addGrainCylindrical!(sim, [20., 0.], 10., 1., tensile_strength=… |
|
161 shear_strength=50e3) |
|
162 sim.grains[1].ang_vel[3] = 100 |
|
163 sim.grains[2].ang_vel[3] = 0.0 |
|
164 Granular.findContacts!(sim) # make sure contact is registered |
|
165 sim.grains[1].contact_radius=10.0 # decrease radius so there isn't compr… |
|
166 E_kin_lin_init = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
167 E_kin_rot_init = Granular.totalGrainKineticRotationalEnergy(sim) |
|
168 Granular.setTimeStep!(sim) |
|
169 Granular.setTotalTime!(sim, 5.) |
|
170 #Granular.setOutputFileInterval!(sim, 0.1) |
|
171 Granular.run!(sim, verbose=verbose) |
|
172 Granular.removeSimulationFiles(sim) |
|
173 @test sim.grains[1].lin_vel[1] ≈ 0. |
|
174 @test sim.grains[1].lin_vel[2] ≈ 0. |
|
175 @test sim.grains[2].lin_vel[1] ≈ 0. |
|
176 @test sim.grains[2].lin_vel[2] ≈ 0. |
|
177 @test sim.grains[1].ang_vel[3] != 0. |
|
178 @test sim.grains[2].ang_vel[3] != 0. |
|
179 E_kin_lin_final = Granular.totalGrainKineticTranslationalEnergy(sim) |
|
180 E_kin_rot_final = Granular.totalGrainKineticRotationalEnergy(sim) |
|
181 E_therm_final = Granular.totalGrainThermalEnergy(sim) |
|
182 @test E_kin_lin_init ≈ E_kin_lin_final |
|
183 @test sim.grains[1].n_contacts == 0 |
|
184 @test sim.grains[2].n_contacts == 0 |
