tshear-test-ocr.py - sphere - GPU-based 3D discrete element method algorithm with optional fluid coupling
(HTM) git clone git://src.adamsgaard.dk/sphere
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
tshear-test-ocr.py (5194B)
---
1 #!/usr/bin/env python
2
3 # Import sphere functionality
4 import sphere
5
6 ### EXPERIMENT SETUP ###
7 initialization = True
8 consolidation = True
9 relaxation = True
10 shearing = True
11 rendering = False
12 plots = True
13
14 # Number of particles
15 np = 2e4
16
17 # Common simulation id
18 sim_id = "shear-test-ocr"
19
20 # Effective normal stresses during consolidation [Pa]
21 Nlist = [10e3, 25e3, 50e3, 100e3, 250e3, 500e3]
22
23 # Effective normal stresses during relaxation and shear [Pa]
24 Nshear = 10e3
25
26 ### INITIALIZATION ###
27
28 # New class
29 init = sphere.sim(np = np, nd = 3, nw = 0, sid = sim_id + "-init")
30
31 # Save radii with uniform size distribution
32 init.generateRadii(psd = 'uni', mean = 1e-2, variance = 2e-3, histogram = True)
33
34 # Set mechanical parameters
35 init.setYoungsModulus(7e8)
36 init.setStaticFriction(0.5)
37 init.setDynamicFriction(0.5)
38 init.setDampingNormal(5e1)
39 init.setDampingTangential(0.0)
40
41 # Add gravitational acceleration
42 init.g[0] = 0.0
43 init.g[1] = 0.0
44 init.g[2] = -9.81
45
46 # Periodic x and y boundaries
47 init.periodicBoundariesXY()
48
49 # Initialize positions in random grid (also sets world size)
50 hcells = np**(1.0/3.0)
51 init.initRandomGridPos(gridnum = [hcells, hcells, 1e9])
52
53 init.checkerboardColors(nx=init.num[0]/2, ny=init.num[1]/2, nz=init.num[2]/2)
54
55 # Set duration of simulation
56 init.initTemporal(total = 10.0, epsilon = 0.07)
57
58 if (initialization == True):
59
60 # Run sphere
61 init.run(dry = True)
62 init.run()
63
64 if (plots == True):
65 # Make a graph of energies
66 init.visualize('energy')
67
68 init.writeVTKall()
69
70 if (rendering == True):
71 # Render images with raytracer
72 init.render(method = "angvel", max_val = 0.3, verbose = False)
73
74
75 # For each normal stress, consolidate and subsequently shear the material
76 for N in Nlist:
77
78 ### CONSOLIDATION ###
79
80 # New class
81 cons = sphere.sim(np = init.np, nw = 1, sid = sim_id +
82 "-cons-N{}".format(N))
83
84 # Read last output file of initialization step
85 lastf = sphere.status(sim_id + "-init")
86 cons.readbin("../output/" + sim_id + "-init.output{:0=5}.bin".format(lastf), verbose=False)
87 cons.setDampingNormal(0.0)
88
89 # Periodic x and y boundaries
90 cons.periodicBoundariesXY()
91
92 # Setup consolidation experiment
93 cons.consolidate(normal_stress = N)
94 cons.adaptiveGrid()
95 cons.checkerboardColors(nx=cons.num[0]/2, ny=cons.num[1]/2, nz=cons.num[2]/2)
96
97 # Set duration of simulation
98 cons.initTemporal(total = 4.0)
99
100 if (consolidation == True):
101
102 # Run sphere
103 cons.run(dry = True) # show values, don't run
104 cons.run() # run
105
106 if (plots == True):
107 # Make a graph of energies
108 cons.visualize('energy')
109 cons.visualize('walls')
110
111 cons.writeVTKall()
112
113 if (rendering == True):
114 # Render images with raytracer
115 cons.render(method = "pres", max_val = 2.0*N, verbose = False)
116
117
118 ### RELAXATION at Nshear ###
119 relax = sphere.sim(np = cons.np, nw = cons.nw, sid = sim_id +
120 "-relax-from-N{}".format(N))
121 lastf = sphere.status(sim_id + "-cons-N{}".format(N))
122 relax.readbin("../output/" + sim_id +
123 "-cons-N{}.output{:0=5}.bin".format(N, lastf))
124
125 relax.periodicBoundariesXY()
126
127 # Setup relaxation experiment
128 relax.consolidate(normal_stress = Nshear)
129 relax.adaptiveGrid()
130 relax.checkerboardColors(nx=relax.num[0]/2, ny=relax.num[1]/2,
131 nz=relax.num[2]/2)
132
133 # Set duration of simulation
134 relax.initTemporal(total = 3.0)
135
136 if (relaxation == True):
137
138 # Run sphere
139 relax.run(dry = True) # show values, don't run
140 relax.run() # run
141
142 if (plots == True):
143 # Make a graph of energies
144 relax.visualize('energy')
145 relax.visualize('walls')
146
147 relax.writeVTKall()
148
149 if (rendering == True):
150 # Render images with raytracer
151 relax.render(method = "pres", max_val = 2.0*Nshear, verbose = False)
152
153 ### SHEARING ###
154
155 # New class
156 shear = sphere.sim(np = relax.np, nw = relax.nw, sid = sim_id +
157 "-shear-N{}-OCR{}".format(Nshear, N/Nshear))
158
159 # Read last output file of initialization step
160 lastf = sphere.status(sim_id + "-relax-from-N{}".format(N))
161 shear.readbin("../output/" + sim_id +
162 "-relax-from-N{}.output{:0=5}.bin".format(N, lastf),
163 verbose = False)
164
165 # Periodic x and y boundaries
166 shear.periodicBoundariesXY()
167
168 # Setup shear experiment
169 shear.shear(shear_strain_rate = 0.05)
170 shear.adaptiveGrid()
171 shear.checkerboardColors(nx=shear.num[0]/2, ny=shear.num[1]/2,
172 nz=shear.num[2]/2)
173
174 # Set duration of simulation
175 shear.initTemporal(total = 20.0)
176
177 if (shearing == True):
178
179 # Run sphere
180 shear.run(dry = True)
181 shear.run()
182
183 if (plots == True):
184 # Make a graph of energies
185 shear.visualize('energy')
186 shear.visualize('shear')
187
188 shear.writeVTKall()
189
190 if (rendering == True):
191 # Render images with raytracer
192 shear.render(method = "pres", max_val = 2.0*N, verbose = False)