trate-state4.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
---
trate-state4.py (3346B)
---
1 #!/usr/bin/env python
2 import sphere
3 import numpy
4 #import sys
5
6 # launch with:
7 # $ ipython sigma-sim1-starter.py <device> <fluid> <c_phi> <k_c> <sigma_0> <mu> <velfac>
8
9 # start with
10 # ipython sigma-sim1-starter.py 0 1 1.0 2.0e-16 10000.0 2.080e-7 1.0
11
12 sid_prefix = 'ratestate4'
13
14 # device = int(sys.argv[1])
15 # wet = int(sys.argv[2])
16 # c_phi = float(sys.argv[3])
17 # k_c = float(sys.argv[4])
18 # sigma0 = float(sys.argv[5])
19 # mu = float(sys.argv[6])
20 # velfac = float(sys.argv[7])
21 device = 0
22 wet = 0
23 c_phi = 1.0
24 k_c = 3.5e-13
25 #sigma0 = 80000.0
26 sigma0 = 40000.0
27 mu = 2.080e-7
28 velfac = 1.0
29
30 start_from_beginning = False
31
32 if wet == 1:
33 fluid = True
34 else:
35 fluid = False
36
37 # load consolidated granular assemblage
38 #sim = sphere.sim('halfshear-sigma0=' + str(sigma0), fluid=False)
39 sim = sphere.sim(fluid=False)
40 if start_from_beginning:
41 sim = sphere.sim('shear-sigma0=' + str(sigma0), fluid=False)
42 sim.readlast()
43 else:
44 if fluid:
45 sim.id('ratestate2-sigma0=' + str(sigma0) + '-k_c=' + str(k_c) + \
46 '-mu=' + str(mu) + '-shear')
47 else:
48 sim.id('ratestate2-sigma0=' + str(sigma0) + '-shear')
49
50 sim.readTime(4.9)
51
52 if fluid:
53 sim.id(sid_prefix + '-sigma0=' + str(sigma0) + '-k_c=' + str(k_c) + \
54 '-mu=' + str(mu) + '-shear')
55 else:
56 sim.id(sid_prefix + '-sigma0=' + str(sigma0) + '-shear')
57
58 #sim.readbin('../input/shear-sigma0=10000.0-new.bin')
59 #sim.scaleSize(0.01) # from 1 cm to 0.01 cm = 100 micro m (fine sand)
60
61
62 sim.fluid = fluid
63 if fluid:
64 sim.id(sid_prefix + '-sigma0=' + str(sigma0) + '-k_c=' + str(k_c) + \
65 '-mu=' + str(mu) + '-shear')
66 else:
67 sim.id(sid_prefix + '-sigma0=' + str(sigma0) + '-shear')
68
69 if start_from_beginning:
70 sim.checkerboardColors(nx=6,ny=3,nz=6)
71 #sim.checkerboardColors(nx=6,ny=6,nz=6)
72 sim.cleanup()
73 sim.adjustUpperWall()
74 sim.zeroKinematics()
75
76 #sim.shear(1.0/20.0 * velfac)
77 sim.shear(1.0/20.0 * velfac * 10.)
78 K_q_real = 36.4e9
79 K_w_real = 2.2e9
80
81
82 K_q_sim = 1.16e9
83 #K_q_sim = 1.16e6
84 sim.setStiffnessNormal(K_q_sim)
85 sim.setStiffnessTangential(K_q_sim)
86 K_w_sim = K_w_real/K_q_real * K_q_sim
87
88
89 if fluid:
90 #sim.num[2] *= 2
91 sim.num[:] /= 2
92 #sim.L[2] *= 2.0
93 #sim.initFluid(mu = 1.787e-6, p = 600.0e3, cfd_solver = 1)
94 sim.initFluid(mu = mu, p = 0.0, cfd_solver = 1)
95 sim.setFluidTopFixedPressure()
96 #sim.setFluidTopFixedFlow()
97 sim.setFluidBottomNoFlow()
98 #sim.setFluidBottomFixedPressure()
99 #sim.setDEMstepsPerCFDstep(10)
100 sim.setMaxIterations(2e5)
101 sim.setPermeabilityPrefactor(k_c)
102 sim.setFluidCompressibility(1.0/K_w_sim)
103
104 sim.w_sigma0[0] = sigma0
105 sim.w_m[0] = numpy.abs(sigma0*sim.L[0]*sim.L[1]/sim.g[2])
106
107 #sim.setStiffnessNormal(36.4e9 * 0.1 / 2.0)
108 #sim.setStiffnessTangential(36.4e9/3.0 * 0.1 / 2.0)
109 sim.setStiffnessNormal(K_q_sim)
110 sim.setStiffnessTangential(K_q_sim)
111 sim.mu_s[0] = 0.5
112 sim.mu_d[0] = 0.5
113 sim.setDampingNormal(0.0)
114 sim.setDampingTangential(0.0)
115
116 sim.initTemporal(total = 20.0/velfac, file_dt = 0.01/velfac, epsilon=0.07)
117
118
119 I = numpy.nonzero(sim.fixvel > 0)
120 sim.fixvel[I] = 8.0 # step-wise velocity change when fixvel in ]5.0; 10.0[
121
122 sim.run(dry=True)
123
124 sim.run(device=device)
125 sim.writeVTKall()
126 sim.visualize('shear')