thansen-zoet.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
---
thansen-zoet.py (3766B)
---
1 #!/usr/bin/env python
2
3 # Import sphere functionality
4 import sphere
5
6 ### EXPERIMENT SETUP ###
7 initialization = False
8 consolidation = True
9 shearing = True
10 rendering = False
11 plots = True
12
13 # Number of particles
14 np = 1e4
15
16 # Common simulation id
17 sim_id = "hz"
18
19 # Use specified compute device (-1 for autoselect)
20 device = 1
21
22 # Deviatoric stress [Pa]
23 #Nlist = [51e3, 101e3, 202e3, 303e3, 404e3]
24 Nlist = [404e3]
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
32 init.generateRadii(mean = 800e-5)
33
34 # Use default params
35 init.defaultParams(gamma_n = 100.0, mu_s = 0.6, mu_d = 0.6)
36
37 # Add gravity
38 init.g[2] = -9.81
39
40 # Periodic x and y boundaries
41 init.periodicBoundariesXY()
42
43 # Initialize positions in random grid (also sets world size)
44 hcells = np**(1.0/3.0)
45 init.initRandomGridPos(gridnum = [hcells, hcells, 1e9])
46
47 # Set duration of simulation
48 init.initTemporal(total = 5.0)
49
50 if (initialization == True):
51
52 # Run sphere
53 init.run(dry = True)
54 init.run(device=device)
55
56 if (plots == True):
57 # Make a graph of energies
58 init.visualize('energy')
59
60 init.writeVTKall()
61
62 if (rendering == True):
63 # Render images with raytracer
64 init.render(method = "angvel", max_val = 0.3, verbose = False)
65
66
67
68 # For each normal stress, consolidate and subsequently shear the material
69 for N in Nlist:
70
71 ### CONSOLIDATION ###
72
73 # New class
74 cons = sphere.sim(np = init.np, nw = 1, sid = sim_id +
75 "-cons-N{}".format(N))
76
77 # Read last output file of initialization step
78 lastf = sphere.status(sim_id + "-init")
79 cons.readbin("../output/" + sim_id + "-init.output{:0=5}.bin".format(lastf), verbose=False)
80
81 # Periodic x and y boundaries
82 cons.periodicBoundariesXY()
83
84 # Setup consolidation experiment
85 cons.consolidate(normal_stress = N)
86 cons.adaptiveGrid()
87
88 # Set duration of simulation
89 cons.initTemporal(total = 1.5)
90
91 """
92 cons.w_m[0] *= 0.001
93 cons.mu_s[0] = 0.0
94 cons.mu_d[0] = 0.0
95 cons.gamma_wn[0] = 1e4
96 cons.gamma_wt[0] = 1e4
97 cons.contactmodel[0] = 1
98 """
99
100 if (consolidation == True):
101
102 # Run sphere
103 cons.run(dry = True) # show values, don't run
104 cons.run(device=device) # 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 ### SHEARING ###
119
120 # New class
121 shear = sphere.sim(np = cons.np, nw = cons.nw, sid = sim_id +
122 "-shear-N{}".format(N))
123
124 # Read last output file of initialization step
125 lastf = sphere.status(sim_id + "-cons-N{}".format(N))
126 shear.readbin("../output/" + sim_id +
127 "-cons-N{}.output{:0=5}.bin".format(N, lastf),
128 verbose = False)
129
130 # Periodic x and y boundaries
131 shear.periodicBoundariesXY()
132
133 # Setup shear experiment
134 shear.shear(shear_strain_rate = 0.10)
135 shear.adaptiveGrid()
136 #shear.initFluid(mu=17.87e-4, p=0.0, hydrostatic=True, cfd_solver=1)
137
138 # Set duration of simulation
139 shear.initTemporal(total = 10.0)
140
141 if (shearing == True):
142
143 # Run sphere
144 shear.run(dry = True)
145 shear.run(device=device)
146
147 if (plots == True):
148 # Make a graph of energies
149 shear.visualize('energy')
150 shear.visualize('shear')
151
152 shear.writeVTKall()
153
154 if (rendering == True):
155 # Render images with raytracer
156 shear.render(method = "pres", max_val = 2.0*N, verbose = False)