tpython_api.html - 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
---
tpython_api.html (202743B)
---
1
2 <!DOCTYPE html>
3
4 <html xmlns="http://www.w3.org/1999/xhtml">
5 <head>
6 <meta charset="utf-8" />
7 <title>Python API — sphere 2.15-beta documentation</title>
8 <link rel="stylesheet" href="_static/classic.css" type="text/css" />
9 <link rel="stylesheet" href="_static/pygments.css" type="text/css" />
10
11 <script type="text/javascript" id="documentation_options" data-url_root="./" src="_static/documentation_options.js"></script>
12 <script type="text/javascript" src="_static/jquery.js"></script>
13 <script type="text/javascript" src="_static/underscore.js"></script>
14 <script type="text/javascript" src="_static/doctools.js"></script>
15 <script type="text/javascript" src="_static/language_data.js"></script>
16
17 <link rel="index" title="Index" href="genindex.html" />
18 <link rel="search" title="Search" href="search.html" />
19 <link rel="next" title="sphere internals" href="sphere_internals.html" />
20 <link rel="prev" title="Fluid simulation and particle-fluid interaction" href="cfd.html" />
21 </head><body>
22 <div class="related" role="navigation" aria-label="related navigation">
23 <h3>Navigation</h3>
24 <ul>
25 <li class="right" style="margin-right: 10px">
26 <a href="genindex.html" title="General Index"
27 accesskey="I">index</a></li>
28 <li class="right" >
29 <a href="py-modindex.html" title="Python Module Index"
30 >modules</a> |</li>
31 <li class="right" >
32 <a href="sphere_internals.html" title="sphere internals"
33 accesskey="N">next</a> |</li>
34 <li class="right" >
35 <a href="cfd.html" title="Fluid simulation and particle-fluid interaction"
36 accesskey="P">previous</a> |</li>
37 <li class="nav-item nav-item-0"><a href="index.html">sphere 2.15-beta documentation</a> »</li>
38 </ul>
39 </div>
40
41 <div class="document">
42 <div class="documentwrapper">
43 <div class="bodywrapper">
44 <div class="body" role="main">
45
46 <div class="section" id="python-api">
47 <h1>Python API<a class="headerlink" href="#python-api" title="Permalink to this headline">¶</a></h1>
48 <p>The Python module <code class="docutils literal notranslate"><span class="pre">sphere</span></code> is intended as the main interface to the <code class="docutils literal notranslate"><span class="pre">sphere</span></code>
49 application. It is recommended to use this module for simulation setup,
50 simulation execution, and analysis of the simulation output data.</p>
51 <p>In order to use the API, the file <code class="docutils literal notranslate"><span class="pre">sphere.py</span></code> must be placed in the same
52 directory as the Python files.</p>
53 <div class="section" id="sample-usage">
54 <h2>Sample usage<a class="headerlink" href="#sample-usage" title="Permalink to this headline">¶</a></h2>
55 <p>Below is a simple, annotated example of how to setup, execute, and post-process
56 a <code class="docutils literal notranslate"><span class="pre">sphere</span></code> simulation. The example is also found in the <code class="docutils literal notranslate"><span class="pre">python/</span></code> folder as
57 <code class="docutils literal notranslate"><span class="pre">collision.py</span></code>.</p>
58 <div class="highlight-python notranslate"><table class="highlighttable"><tr><td class="linenos"><div class="linenodiv"><pre> 1
59 2
60 3
61 4
62 5
63 6
64 7
65 8
66 9
67 10
68 11
69 12
70 13
71 14
72 15
73 16
74 17
75 18
76 19
77 20
78 21
79 22
80 23
81 24
82 25
83 26
84 27
85 28
86 29
87 30
88 31
89 32
90 33
91 34
92 35
93 36
94 37
95 38
96 39
97 40
98 41
99 42
100 43
101 44
102 45
103 46
104 47
105 48
106 49
107 50
108 51
109 52
110 53
111 54
112 55
113 56
114 57
115 58
116 59</pre></div></td><td class="code"><div class="highlight"><pre><span></span><span class="ch">#!/usr/bin/env python</span>
117 <span class="sd">'''</span>
118 <span class="sd">Example of two particles colliding.</span>
119 <span class="sd">Place script in sphere/python/ folder, and invoke with `python collision.py`</span>
120 <span class="sd">'''</span>
121
122 <span class="c1"># Import the sphere module for setting up, running, and analyzing the</span>
123 <span class="c1"># experiment. We also need the numpy module when setting arrays in the sphere</span>
124 <span class="c1"># object.</span>
125 <span class="kn">import</span> <span class="nn">sphere</span>
126 <span class="kn">import</span> <span class="nn">numpy</span>
127
128
129 <span class="c1">### SIMULATION SETUP</span>
130
131 <span class="c1"># Create a sphere object with two preallocated particles and a simulation ID</span>
132 <span class="n">SB</span> <span class="o">=</span> <span class="n">sphere</span><span class="o">.</span><span class="n">sim</span><span class="p">(</span><span class="n">np</span> <span class="o">=</span> <span class="mi">2</span><span class="p">,</span> <span class="n">sid</span> <span class="o">=</span> <span class="s1">'collision'</span><span class="p">)</span>
133
134 <span class="n">SB</span><span class="o">.</span><span class="n">radius</span><span class="p">[:]</span> <span class="o">=</span> <span class="mf">0.3</span> <span class="c1"># set radii to 0.3 m</span>
135
136 <span class="c1"># Define the positions of the two particles</span>
137 <span class="n">SB</span><span class="o">.</span><span class="n">x</span><span class="p">[</span><span class="mi">0</span><span class="p">,</span> <span class="p">:]</span> <span class="o">=</span> <span class="n">numpy</span><span class="o">.</span><span class="n">array</span><span class="p">([</span><span class="mf">10.0</span><span class="p">,</span> <span class="mf">5.0</span><span class="p">,</span> <span class="mf">5.0</span><span class="p">])</span> <span class="c1"># particle 1 (idx 0)</span>
138 <span class="n">SB</span><span class="o">.</span><span class="n">x</span><span class="p">[</span><span class="mi">1</span><span class="p">,</span> <span class="p">:]</span> <span class="o">=</span> <span class="n">numpy</span><span class="o">.</span><span class="n">array</span><span class="p">([</span><span class="mf">11.0</span><span class="p">,</span> <span class="mf">5.0</span><span class="p">,</span> <span class="mf">5.0</span><span class="p">])</span> <span class="c1"># particle 2 (idx 1)</span>
139
140 <span class="c1"># The default velocity is [0,0,0]. Slam particle 1 into particle 2 by defining</span>
141 <span class="c1"># a positive x velocity for particle 1.</span>
142 <span class="n">SB</span><span class="o">.</span><span class="n">vel</span><span class="p">[</span><span class="mi">0</span><span class="p">,</span> <span class="mi">0</span><span class="p">]</span> <span class="o">=</span> <span class="mf">1.0</span>
143
144 <span class="c1"># Set the world limits and the particle sorting grid. The particles need to stay</span>
145 <span class="c1"># within the world limits for the entire simulation, otherwise it will stop!</span>
146 <span class="n">SB</span><span class="o">.</span><span class="n">initGridAndWorldsize</span><span class="p">(</span><span class="n">margin</span> <span class="o">=</span> <span class="mf">5.0</span><span class="p">)</span>
147
148 <span class="c1"># Define the temporal parameters, e.g. the total time (total) and the file</span>
149 <span class="c1"># output interval (file_dt), both in seconds</span>
150 <span class="n">SB</span><span class="o">.</span><span class="n">initTemporal</span><span class="p">(</span><span class="n">total</span> <span class="o">=</span> <span class="mf">2.0</span><span class="p">,</span> <span class="n">file_dt</span> <span class="o">=</span> <span class="mf">0.1</span><span class="p">)</span>
151
152 <span class="c1"># Using a 'dry' run, the sphere main program will display important parameters.</span>
153 <span class="c1"># sphere will end after displaying these values.</span>
154 <span class="n">SB</span><span class="o">.</span><span class="n">run</span><span class="p">(</span><span class="n">dry</span> <span class="o">=</span> <span class="bp">True</span><span class="p">)</span>
155
156
157 <span class="c1">### RUNNING THE SIMULATION</span>
158
159 <span class="c1"># Start the simulation on the GPU from the sphere program</span>
160 <span class="n">SB</span><span class="o">.</span><span class="n">run</span><span class="p">()</span>
161
162
163 <span class="c1">### ANALYSIS OF SIMULATION RESULTS</span>
164
165 <span class="c1"># Plot the system energy through time, image saved as collision-energy.png</span>
166 <span class="n">SB</span><span class="o">.</span><span class="n">visualize</span><span class="p">(</span><span class="n">method</span> <span class="o">=</span> <span class="s1">'energy'</span><span class="p">)</span>
167
168 <span class="c1"># Render the particles using the built-in raytracer</span>
169 <span class="n">SB</span><span class="o">.</span><span class="n">render</span><span class="p">()</span>
170
171 <span class="c1"># Alternative visualization using ParaView. See the documentation of</span>
172 <span class="c1"># ``sim.writeVTKall()`` for more information about displaying the</span>
173 <span class="c1"># particles in ParaView.</span>
174 <span class="n">SB</span><span class="o">.</span><span class="n">writeVTKall</span><span class="p">()</span>
175 </pre></div>
176 </td></tr></table></div>
177 <p>The full documentation of the <code class="docutils literal notranslate"><span class="pre">sphere</span></code> Python API can be found below.</p>
178 </div>
179 <div class="section" id="module-sphere">
180 <span id="the-sphere-module"></span><h2>The <code class="docutils literal notranslate"><span class="pre">sphere</span></code> module<a class="headerlink" href="#module-sphere" title="Permalink to this headline">¶</a></h2>
181 <dl class="function">
182 <dt id="sphere.V_sphere">
183 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">V_sphere</code><span class="sig-paren">(</span><em class="sig-param">r</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.V_sphere" title="Permalink to this definition">¶</a></dt>
184 <dd><p>Calculates the volume of a sphere with radius r</p>
185 <dl class="field-list simple">
186 <dt class="field-odd">Returns</dt>
187 <dd class="field-odd"><p>The sphere volume [m^3]</p>
188 </dd>
189 <dt class="field-even">Return type</dt>
190 <dd class="field-even"><p>float</p>
191 </dd>
192 </dl>
193 </dd></dl>
194
195 <dl class="function">
196 <dt id="sphere.cleanup">
197 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">cleanup</code><span class="sig-paren">(</span><em class="sig-param">sim</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.cleanup" title="Permalink to this definition">¶</a></dt>
198 <dd><p>Removes the input/output files and images belonging to the object simulation
199 ID from the <code class="docutils literal notranslate"><span class="pre">input/</span></code>, <code class="docutils literal notranslate"><span class="pre">output/</span></code> and <code class="docutils literal notranslate"><span class="pre">img_out/</span></code> folders.</p>
200 <dl class="field-list simple">
201 <dt class="field-odd">Parameters</dt>
202 <dd class="field-odd"><p><strong>spherebin</strong> (<a class="reference internal" href="#sphere.sim" title="sphere.sim"><em>sim</em></a>) – A sim object</p>
203 </dd>
204 </dl>
205 </dd></dl>
206
207 <dl class="function">
208 <dt id="sphere.convert">
209 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">convert</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='png'</em>, <em class="sig-param">folder='../img_out'</em>, <em class="sig-param">remove_ppm=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.convert" title="Permalink to this definition">¶</a></dt>
210 <dd><p>Converts all PPM images in img_out to graphics_format using ImageMagick. All
211 PPM images are subsequently removed if <cite>remove_ppm</cite> is <cite>True</cite>.</p>
212 <dl class="field-list simple">
213 <dt class="field-odd">Parameters</dt>
214 <dd class="field-odd"><ul class="simple">
215 <li><p><strong>graphics_format</strong> (<em>str</em>) – Convert the images to this format</p></li>
216 <li><p><strong>folder</strong> (<em>str</em>) – The folder containing the PPM images to convert</p></li>
217 <li><p><strong>remove_ppm</strong> (<em>bool</em>) – Remove ALL ppm files in <cite>folder</cite> after conversion</p></li>
218 </ul>
219 </dd>
220 </dl>
221 </dd></dl>
222
223 <dl class="function">
224 <dt id="sphere.render">
225 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">render</code><span class="sig-paren">(</span><em class="sig-param">binary</em>, <em class="sig-param">method='pres'</em>, <em class="sig-param">max_val=1000.0</em>, <em class="sig-param">lower_cutoff=0.0</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.render" title="Permalink to this definition">¶</a></dt>
226 <dd><p>Render target binary using the <code class="docutils literal notranslate"><span class="pre">sphere</span></code> raytracer.</p>
227 <dl class="field-list simple">
228 <dt class="field-odd">Parameters</dt>
229 <dd class="field-odd"><ul class="simple">
230 <li><p><strong>method</strong> (<em>str</em>) – The color visualization method to use for the particles.
231 Possible values are: ‘normal’: color all particles with the same
232 color, ‘pres’: color by pressure, ‘vel’: color by translational
233 velocity, ‘angvel’: color by rotational velocity, ‘xdisp’: color by
234 total displacement along the x-axis, ‘angpos’: color by angular
235 position.</p></li>
236 <li><p><strong>max_val</strong> (<em>float</em>) – The maximum value of the color bar</p></li>
237 <li><p><strong>lower_cutoff</strong> (<em>float</em>) – Do not render particles with a value below this
238 value, of the field selected by <code class="docutils literal notranslate"><span class="pre">method</span></code></p></li>
239 <li><p><strong>graphics_format</strong> (<em>str</em>) – Convert the PPM images generated by the ray
240 tracer to this image format using Imagemagick</p></li>
241 <li><p><strong>verbose</strong> (<em>bool</em>) – Show verbose information during ray tracing</p></li>
242 </ul>
243 </dd>
244 </dl>
245 </dd></dl>
246
247 <dl class="function">
248 <dt id="sphere.run">
249 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">run</code><span class="sig-paren">(</span><em class="sig-param">binary</em>, <em class="sig-param">verbose=True</em>, <em class="sig-param">hideinputfile=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.run" title="Permalink to this definition">¶</a></dt>
250 <dd><p>Execute <code class="docutils literal notranslate"><span class="pre">sphere</span></code> with target binary file as input.</p>
251 <dl class="field-list simple">
252 <dt class="field-odd">Parameters</dt>
253 <dd class="field-odd"><ul class="simple">
254 <li><p><strong>binary</strong> (<em>str</em>) – Input file for <code class="docutils literal notranslate"><span class="pre">sphere</span></code></p></li>
255 <li><p><strong>verbose</strong> (<em>bool</em>) – Show <code class="docutils literal notranslate"><span class="pre">sphere</span></code> output</p></li>
256 <li><p><strong>hideinputfile</strong> (<em>bool</em>) – Hide the input file</p></li>
257 </ul>
258 </dd>
259 </dl>
260 </dd></dl>
261
262 <dl class="class">
263 <dt id="sphere.sim">
264 <em class="property">class </em><code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">sim</code><span class="sig-paren">(</span><em class="sig-param">sid='unnamed'</em>, <em class="sig-param">np=0</em>, <em class="sig-param">nd=3</em>, <em class="sig-param">nw=0</em>, <em class="sig-param">fluid=False</em>, <em class="sig-param">cfd_solver=0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim" title="Permalink to this definition">¶</a></dt>
265 <dd><p>Class containing all <code class="docutils literal notranslate"><span class="pre">sphere</span></code> data.</p>
266 <p>Contains functions for reading and writing binaries, as well as simulation
267 setup and data analysis. Most arrays are initialized to default values.</p>
268 <dl class="field-list simple">
269 <dt class="field-odd">Parameters</dt>
270 <dd class="field-odd"><ul class="simple">
271 <li><p><strong>np</strong> (<em>int</em>) – The number of particles to allocate memory for (default=1)</p></li>
272 <li><p><strong>nd</strong> (<em>int</em>) – The number of spatial dimensions (default=3). Note that 2D and
273 1D simulations currently are not possible.</p></li>
274 <li><p><strong>nw</strong> (<em>int</em>) – The number of dynamic walls (default=1)</p></li>
275 <li><p><strong>sid</strong> (<em>str</em>) – The simulation id (default=’unnamed’). The simulation files
276 will be written with this base name.</p></li>
277 <li><p><strong>fluid</strong> (<em>bool</em>) – Setup fluid simulation (default=False)</p></li>
278 <li><p><strong>cfd_solver</strong> (<em>int</em>) – Fluid solver to use if fluid == True. 0: Navier-Stokes
279 (default), 1: Darcy.</p></li>
280 </ul>
281 </dd>
282 </dl>
283 <dl class="method">
284 <dt id="sphere.sim.ReynoldsNumber">
285 <code class="sig-name descname">ReynoldsNumber</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.ReynoldsNumber" title="Permalink to this definition">¶</a></dt>
286 <dd><p>Estimate the per-cell Reynolds number by: Re=rho * ||v_f|| * dx/mu.
287 This value is returned and also stored in <cite>self.Re</cite>.</p>
288 <dl class="field-list simple">
289 <dt class="field-odd">Returns</dt>
290 <dd class="field-odd"><p>Reynolds number</p>
291 </dd>
292 <dt class="field-even">Return type</dt>
293 <dd class="field-even"><p>Numpy array with dimensions like the fluid grid</p>
294 </dd>
295 </dl>
296 </dd></dl>
297
298 <dl class="method">
299 <dt id="sphere.sim.acceleration">
300 <code class="sig-name descname">acceleration</code><span class="sig-paren">(</span><em class="sig-param">idx=-1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.acceleration" title="Permalink to this definition">¶</a></dt>
301 <dd><p>Returns the acceleration of one or more particles, selected by their
302 index. If the index is equal to -1 (default value), all accelerations
303 are returned.</p>
304 <dl class="field-list simple">
305 <dt class="field-odd">Parameters</dt>
306 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em><em>, </em><em>list</em><em> or </em><em>numpy.array</em>) – Index or index range of particles</p>
307 </dd>
308 <dt class="field-even">Returns</dt>
309 <dd class="field-even"><p>n-by-3 matrix of acceleration(s)</p>
310 </dd>
311 <dt class="field-odd">Return type</dt>
312 <dd class="field-odd"><p>numpy.array</p>
313 </dd>
314 </dl>
315 </dd></dl>
316
317 <dl class="method">
318 <dt id="sphere.sim.adaptiveGrid">
319 <code class="sig-name descname">adaptiveGrid</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.adaptiveGrid" title="Permalink to this definition">¶</a></dt>
320 <dd><p>Set the height of the fluid grid to automatically readjust to the
321 height of the granular assemblage, as dictated by the position of the
322 top wall. This will readjust <cite>self.L[2]</cite> during the simulation to
323 equal the position of the top wall <cite>self.w_x[0]</cite>.</p>
324 <p>See also <a class="reference internal" href="#sphere.sim.staticGrid" title="sphere.sim.staticGrid"><code class="xref py py-func docutils literal notranslate"><span class="pre">staticGrid()</span></code></a></p>
325 </dd></dl>
326
327 <dl class="method">
328 <dt id="sphere.sim.addParticle">
329 <code class="sig-name descname">addParticle</code><span class="sig-paren">(</span><em class="sig-param">x</em>, <em class="sig-param">radius</em>, <em class="sig-param">xyzsum=array([0.</em>, <em class="sig-param">0.</em>, <em class="sig-param">0.])</em>, <em class="sig-param">vel=array([0.</em>, <em class="sig-param">0.</em>, <em class="sig-param">0.])</em>, <em class="sig-param">fixvel=array([0.])</em>, <em class="sig-param">force=array([0.</em>, <em class="sig-param">0.</em>, <em class="sig-param">0.])</em>, <em class="sig-param">angpos=array([0.</em>, <em class="sig-param">0.</em>, <em class="sig-param">0.])</em>, <em class="sig-param">angvel=array([0.</em>, <em class="sig-param">0.</em>, <em class="sig-param">0.])</em>, <em class="sig-param">torque=array([0.</em>, <em class="sig-param">0.</em>, <em class="sig-param">0.])</em>, <em class="sig-param">es_dot=array([0.])</em>, <em class="sig-param">es=array([0.])</em>, <em class="sig-param">ev_dot=array([0.])</em>, <em class="sig-param">ev=array([0.])</em>, <em class="sig-param">p=array([0.])</em>, <em class="sig-param">color=0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.addParticle" title="Permalink to this definition">¶</a></dt>
330 <dd><p>Add a single particle to the simulation object. The only required
331 parameters are the position (x) and the radius (radius).</p>
332 <dl class="field-list simple">
333 <dt class="field-odd">Parameters</dt>
334 <dd class="field-odd"><ul class="simple">
335 <li><p><strong>x</strong> (<em>numpy.array</em>) – A vector pointing to the particle center coordinate.</p></li>
336 <li><p><strong>radius</strong> (<em>float</em>) – The particle radius</p></li>
337 <li><p><strong>vel</strong> (<em>numpy.array</em>) – The particle linear velocity (default=[0,0,0])</p></li>
338 <li><p><strong>fixvel</strong> (<em>float</em>) – 0: Do not fix particle velocity (default), 1: Fix
339 horizontal linear velocity, -1: Fix horizontal and vertical linear
340 velocity</p></li>
341 <li><p><strong>angpos</strong> (<em>numpy.array</em>) – The particle angular position (default=[0,0,0])</p></li>
342 <li><p><strong>angvel</strong> (<em>numpy.array</em>) – The particle angular velocity (default=[0,0,0])</p></li>
343 <li><p><strong>torque</strong> (<em>numpy.array</em>) – The particle torque (default=[0,0,0])</p></li>
344 <li><p><strong>es_dot</strong> (<em>float</em>) – The particle shear energy loss rate (default=0)</p></li>
345 <li><p><strong>es</strong> (<em>float</em>) – The particle shear energy loss (default=0)</p></li>
346 <li><p><strong>ev_dot</strong> (<em>float</em>) – The particle viscous energy rate loss (default=0)</p></li>
347 <li><p><strong>ev</strong> (<em>float</em>) – The particle viscous energy loss (default=0)</p></li>
348 <li><p><strong>p</strong> (<em>float</em>) – The particle pressure (default=0)</p></li>
349 </ul>
350 </dd>
351 </dl>
352 </dd></dl>
353
354 <dl class="method">
355 <dt id="sphere.sim.adjustUpperWall">
356 <code class="sig-name descname">adjustUpperWall</code><span class="sig-paren">(</span><em class="sig-param">z_adjust=1.1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.adjustUpperWall" title="Permalink to this definition">¶</a></dt>
357 <dd><p>Included for legacy purposes, calls <a class="reference internal" href="#sphere.sim.adjustWall" title="sphere.sim.adjustWall"><code class="xref py py-func docutils literal notranslate"><span class="pre">adjustWall()</span></code></a> with <code class="docutils literal notranslate"><span class="pre">idx=0</span></code>.</p>
358 <dl class="field-list simple">
359 <dt class="field-odd">Parameters</dt>
360 <dd class="field-odd"><p><strong>z_adjust</strong> (<em>float</em>) – Increase the world and grid size by this amount to
361 allow for wall movement.</p>
362 </dd>
363 </dl>
364 </dd></dl>
365
366 <dl class="method">
367 <dt id="sphere.sim.adjustWall">
368 <code class="sig-name descname">adjustWall</code><span class="sig-paren">(</span><em class="sig-param">idx</em>, <em class="sig-param">adjust=1.1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.adjustWall" title="Permalink to this definition">¶</a></dt>
369 <dd><p>Adjust grid and dynamic wall to max. particle position. The wall
370 thickness will by standard equal the maximum particle diameter. The
371 density equals the particle density, and the wall size is equal to the
372 width and depth of the simulation domain (<cite>self.L[0]</cite> and <cite>self.L[1]</cite>).</p>
373 <dl class="field-list simple">
374 <dt class="field-odd">Param</dt>
375 <dd class="field-odd"><p>idx: The wall to adjust. 0=+z, upper wall (default), 1=-x,
376 left wall, 2=+x, right wall, 3=-y, front wall, 4=+y, back
377 wall.</p>
378 </dd>
379 <dt class="field-even">Parameters</dt>
380 <dd class="field-even"><p><strong>adjust</strong> (<em>float</em>) – Increase the world and grid size by this amount to
381 allow for wall movement.</p>
382 </dd>
383 </dl>
384 </dd></dl>
385
386 <dl class="method">
387 <dt id="sphere.sim.bond">
388 <code class="sig-name descname">bond</code><span class="sig-paren">(</span><em class="sig-param">i</em>, <em class="sig-param">j</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.bond" title="Permalink to this definition">¶</a></dt>
389 <dd><p>Create a bond between particles with index i and j</p>
390 <dl class="field-list simple">
391 <dt class="field-odd">Parameters</dt>
392 <dd class="field-odd"><ul class="simple">
393 <li><p><strong>i</strong> (<em>int</em>) – Index of first particle in bond</p></li>
394 <li><p><strong>j</strong> (<em>int</em>) – Index of second particle in bond</p></li>
395 </ul>
396 </dd>
397 </dl>
398 </dd></dl>
399
400 <dl class="method">
401 <dt id="sphere.sim.bondsRose">
402 <code class="sig-name descname">bondsRose</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='pdf'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.bondsRose" title="Permalink to this definition">¶</a></dt>
403 <dd><p>Visualize the trend and plunge angles of the bond pairs in a rose plot.
404 The plot is saved in the current folder as
405 ‘bonds-<simulation id>-rose.<graphics_format>’.</p>
406 <dl class="field-list simple">
407 <dt class="field-odd">Parameters</dt>
408 <dd class="field-odd"><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p>
409 </dd>
410 </dl>
411 </dd></dl>
412
413 <dl class="method">
414 <dt id="sphere.sim.bulkPorosity">
415 <code class="sig-name descname">bulkPorosity</code><span class="sig-paren">(</span><em class="sig-param">trim=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.bulkPorosity" title="Permalink to this definition">¶</a></dt>
416 <dd><p>Calculates the bulk porosity of the particle assemblage.</p>
417 <dl class="field-list simple">
418 <dt class="field-odd">Parameters</dt>
419 <dd class="field-odd"><p><strong>trim</strong> (<em>bool</em>) – Trim the total volume to the smallest axis-parallel cube
420 containing all particles.</p>
421 </dd>
422 <dt class="field-even">Returns</dt>
423 <dd class="field-even"><p>The bulk porosity, in [0:1]</p>
424 </dd>
425 <dt class="field-odd">Return type</dt>
426 <dd class="field-odd"><p>float</p>
427 </dd>
428 </dl>
429 </dd></dl>
430
431 <dl class="method">
432 <dt id="sphere.sim.cellSize">
433 <code class="sig-name descname">cellSize</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.cellSize" title="Permalink to this definition">¶</a></dt>
434 <dd><p>Calculate the particle sorting (and fluid) cell dimensions.
435 These values are stored in <cite>self.dx</cite> and are NOT returned.</p>
436 </dd></dl>
437
438 <dl class="method">
439 <dt id="sphere.sim.checkerboardColors">
440 <code class="sig-name descname">checkerboardColors</code><span class="sig-paren">(</span><em class="sig-param">nx=6</em>, <em class="sig-param">ny=6</em>, <em class="sig-param">nz=6</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.checkerboardColors" title="Permalink to this definition">¶</a></dt>
441 <dd><p>Assign checkerboard color values to the particles in an orthogonal grid.</p>
442 <dl class="field-list simple">
443 <dt class="field-odd">Parameters</dt>
444 <dd class="field-odd"><ul class="simple">
445 <li><p><strong>nx</strong> (<em>int</em>) – Number of color values along the x axis</p></li>
446 <li><p><strong>ny</strong> (<em>int</em>) – Number of color values along the y ayis</p></li>
447 <li><p><strong>nz</strong> (<em>int</em>) – Number of color values along the z azis</p></li>
448 </ul>
449 </dd>
450 </dl>
451 </dd></dl>
452
453 <dl class="method">
454 <dt id="sphere.sim.cleanup">
455 <code class="sig-name descname">cleanup</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.cleanup" title="Permalink to this definition">¶</a></dt>
456 <dd><p>Removes the input/output files and images belonging to the object
457 simulation ID from the <code class="docutils literal notranslate"><span class="pre">input/</span></code>, <code class="docutils literal notranslate"><span class="pre">output/</span></code> and <code class="docutils literal notranslate"><span class="pre">img_out/</span></code> folders.</p>
458 </dd></dl>
459
460 <dl class="method">
461 <dt id="sphere.sim.consolidate">
462 <code class="sig-name descname">consolidate</code><span class="sig-paren">(</span><em class="sig-param">normal_stress=10000.0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.consolidate" title="Permalink to this definition">¶</a></dt>
463 <dd><p>Setup consolidation experiment. Specify the upper wall normal stress in
464 Pascal, default value is 10 kPa.</p>
465 <dl class="field-list simple">
466 <dt class="field-odd">Parameters</dt>
467 <dd class="field-odd"><p><strong>normal_stress</strong> (<em>float</em>) – The normal stress to apply from the upper wall</p>
468 </dd>
469 </dl>
470 </dd></dl>
471
472 <dl class="method">
473 <dt id="sphere.sim.contactModel">
474 <code class="sig-name descname">contactModel</code><span class="sig-paren">(</span><em class="sig-param">contactmodel</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.contactModel" title="Permalink to this definition">¶</a></dt>
475 <dd><p>Define which contact model to use for the tangential component of
476 particle-particle interactions. The elastic-viscous-frictional contact
477 model (2) is considered to be the most realistic contact model, while
478 the viscous-frictional contact model is significantly faster.</p>
479 <dl class="field-list simple">
480 <dt class="field-odd">Parameters</dt>
481 <dd class="field-odd"><p><strong>contactmodel</strong> (<em>int</em>) – The type of tangential contact model to use
482 (visco-frictional=1, elasto-visco-frictional=2)</p>
483 </dd>
484 </dl>
485 </dd></dl>
486
487 <dl class="method">
488 <dt id="sphere.sim.contactParticleArea">
489 <code class="sig-name descname">contactParticleArea</code><span class="sig-paren">(</span><em class="sig-param">i</em>, <em class="sig-param">j</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.contactParticleArea" title="Permalink to this definition">¶</a></dt>
490 <dd><p>Finds the average area of an two particles in an inter-particle contact.</p>
491 <dl class="field-list simple">
492 <dt class="field-odd">Parameters</dt>
493 <dd class="field-odd"><ul class="simple">
494 <li><p><strong>i</strong> (<em>int</em><em> or </em><em>array of ints</em>) – Index of first particle</p></li>
495 <li><p><strong>j</strong> (<em>int</em><em> or </em><em>array of ints</em>) – Index of second particle</p></li>
496 <li><p><strong>d</strong> (<em>float</em><em> or </em><em>array of floats</em>) – Overlap distance</p></li>
497 </ul>
498 </dd>
499 <dt class="field-even">Returns</dt>
500 <dd class="field-even"><p>Contact area [m*m]</p>
501 </dd>
502 <dt class="field-odd">Return type</dt>
503 <dd class="field-odd"><p>float or array of floats</p>
504 </dd>
505 </dl>
506 </dd></dl>
507
508 <dl class="method">
509 <dt id="sphere.sim.contactSurfaceArea">
510 <code class="sig-name descname">contactSurfaceArea</code><span class="sig-paren">(</span><em class="sig-param">i</em>, <em class="sig-param">j</em>, <em class="sig-param">overlap</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.contactSurfaceArea" title="Permalink to this definition">¶</a></dt>
511 <dd><p>Finds the contact surface area of an inter-particle contact.</p>
512 <dl class="field-list simple">
513 <dt class="field-odd">Parameters</dt>
514 <dd class="field-odd"><ul class="simple">
515 <li><p><strong>i</strong> (<em>int</em><em> or </em><em>array of ints</em>) – Index of first particle</p></li>
516 <li><p><strong>j</strong> (<em>int</em><em> or </em><em>array of ints</em>) – Index of second particle</p></li>
517 <li><p><strong>d</strong> (<em>float</em><em> or </em><em>array of floats</em>) – Overlap distance</p></li>
518 </ul>
519 </dd>
520 <dt class="field-even">Returns</dt>
521 <dd class="field-even"><p>Contact area [m*m]</p>
522 </dd>
523 <dt class="field-odd">Return type</dt>
524 <dd class="field-odd"><p>float or array of floats</p>
525 </dd>
526 </dl>
527 </dd></dl>
528
529 <dl class="method">
530 <dt id="sphere.sim.convergence">
531 <code class="sig-name descname">convergence</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.convergence" title="Permalink to this definition">¶</a></dt>
532 <dd><p>Read the convergence evolution in the CFD solver. The values are stored
533 in <cite>self.conv</cite> with iteration number in the first column and iteration
534 count in the second column.</p>
535 <p>See also: <a class="reference internal" href="#sphere.sim.plotConvergence" title="sphere.sim.plotConvergence"><code class="xref py py-func docutils literal notranslate"><span class="pre">plotConvergence()</span></code></a></p>
536 </dd></dl>
537
538 <dl class="method">
539 <dt id="sphere.sim.createBondPair">
540 <code class="sig-name descname">createBondPair</code><span class="sig-paren">(</span><em class="sig-param">i</em>, <em class="sig-param">j</em>, <em class="sig-param">spacing=-0.1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.createBondPair" title="Permalink to this definition">¶</a></dt>
541 <dd><p>Bond particles i and j. Particle j is moved adjacent to particle i,
542 and oriented randomly.</p>
543 <dl class="field-list simple">
544 <dt class="field-odd">Parameters</dt>
545 <dd class="field-odd"><ul class="simple">
546 <li><p><strong>i</strong> (<em>int</em>) – Index of first particle in bond</p></li>
547 <li><p><strong>j</strong> (<em>int</em>) – Index of second particle in bond</p></li>
548 <li><p><strong>spacing</strong> (<em>float</em>) – The inter-particle distance prescribed. Positive
549 values result in a inter-particle distance, negative equal an
550 overlap. The value is relative to the sum of the two radii.</p></li>
551 </ul>
552 </dd>
553 </dl>
554 </dd></dl>
555
556 <dl class="method">
557 <dt id="sphere.sim.currentNormalStress">
558 <code class="sig-name descname">currentNormalStress</code><span class="sig-paren">(</span><em class="sig-param">type='defined'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.currentNormalStress" title="Permalink to this definition">¶</a></dt>
559 <dd><p>Calculates the current magnitude of the defined or effective top wall
560 normal stress.</p>
561 <dl class="field-list simple">
562 <dt class="field-odd">Parameters</dt>
563 <dd class="field-odd"><p><strong>type</strong> (<em>str</em>) – Find the ‘defined’ (default) or ‘effective’ normal stress</p>
564 </dd>
565 <dt class="field-even">Returns</dt>
566 <dd class="field-even"><p>The current top wall normal stress in Pascal</p>
567 </dd>
568 <dt class="field-odd">Return type</dt>
569 <dd class="field-odd"><p>float</p>
570 </dd>
571 </dl>
572 </dd></dl>
573
574 <dl class="method">
575 <dt id="sphere.sim.currentTime">
576 <code class="sig-name descname">currentTime</code><span class="sig-paren">(</span><em class="sig-param">value=-1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.currentTime" title="Permalink to this definition">¶</a></dt>
577 <dd><p>Get or set the current time. If called without arguments the current
578 time is returned. If a new time is passed in the ‘value’ argument, the
579 time is written to the object.</p>
580 <dl class="field-list simple">
581 <dt class="field-odd">Parameters</dt>
582 <dd class="field-odd"><p><strong>value</strong> (<em>float</em>) – The new current time</p>
583 </dd>
584 <dt class="field-even">Returns</dt>
585 <dd class="field-even"><p>The current time</p>
586 </dd>
587 <dt class="field-odd">Return type</dt>
588 <dd class="field-odd"><p>float</p>
589 </dd>
590 </dl>
591 </dd></dl>
592
593 <dl class="method">
594 <dt id="sphere.sim.defaultParams">
595 <code class="sig-name descname">defaultParams</code><span class="sig-paren">(</span><em class="sig-param">mu_s=0.5</em>, <em class="sig-param">mu_d=0.5</em>, <em class="sig-param">mu_r=0.0</em>, <em class="sig-param">rho=2600</em>, <em class="sig-param">k_n=1160000000.0</em>, <em class="sig-param">k_t=1160000000.0</em>, <em class="sig-param">k_r=0</em>, <em class="sig-param">gamma_n=0.0</em>, <em class="sig-param">gamma_t=0.0</em>, <em class="sig-param">gamma_r=0.0</em>, <em class="sig-param">gamma_wn=0.0</em>, <em class="sig-param">gamma_wt=0.0</em>, <em class="sig-param">capillaryCohesion=0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.defaultParams" title="Permalink to this definition">¶</a></dt>
596 <dd><p>Initialize particle parameters to default values.</p>
597 <dl class="field-list simple">
598 <dt class="field-odd">Parameters</dt>
599 <dd class="field-odd"><ul class="simple">
600 <li><p><strong>mu_s</strong> (<em>float</em>) – The coefficient of static friction between particles [-]</p></li>
601 <li><p><strong>mu_d</strong> (<em>float</em>) – The coefficient of dynamic friction between particles [-]</p></li>
602 <li><p><strong>rho</strong> (<em>float</em>) – The density of the particle material [kg/(m^3)]</p></li>
603 <li><p><strong>k_n</strong> (<em>float</em>) – The normal stiffness of the particles [N/m]</p></li>
604 <li><p><strong>k_t</strong> (<em>float</em>) – The tangential stiffness of the particles [N/m]</p></li>
605 <li><p><strong>k_r</strong> (<em>float</em>) – The rolling stiffness of the particles [N/rad] <em>Parameter
606 not used</em></p></li>
607 <li><p><strong>gamma_n</strong> (<em>float</em>) – Particle-particle contact normal viscosity [Ns/m]</p></li>
608 <li><p><strong>gamma_t</strong> (<em>float</em>) – Particle-particle contact tangential viscosity [Ns/m]</p></li>
609 <li><p><strong>gamma_r</strong> (<em>float</em>) – Particle-particle contact rolling viscosity <em>Parameter
610 not used</em></p></li>
611 <li><p><strong>gamma_wn</strong> (<em>float</em>) – Wall-particle contact normal viscosity [Ns/m]</p></li>
612 <li><p><strong>gamma_wt</strong> (<em>float</em>) – Wall-particle contact tangential viscosity [Ns/m]</p></li>
613 <li><p><strong>capillaryCohesion</strong> (<em>int</em>) – Enable particle-particle capillary cohesion
614 interaction model (0=no (default), 1=yes)</p></li>
615 </ul>
616 </dd>
617 </dl>
618 </dd></dl>
619
620 <dl class="method">
621 <dt id="sphere.sim.defineWorldBoundaries">
622 <code class="sig-name descname">defineWorldBoundaries</code><span class="sig-paren">(</span><em class="sig-param">L, origo=[0.0, 0.0, 0.0], dx=-1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.defineWorldBoundaries" title="Permalink to this definition">¶</a></dt>
623 <dd><p>Set the boundaries of the world. Particles will only be able to interact
624 within this domain. With dynamic walls, allow space for expansions.
625 <em>Important</em>: The particle radii have to be set beforehand. The world
626 edges act as static walls.</p>
627 <dl class="field-list simple">
628 <dt class="field-odd">Parameters</dt>
629 <dd class="field-odd"><ul class="simple">
630 <li><p><strong>L</strong> (<em>numpy.array</em>) – The upper boundary of the domain [m]</p></li>
631 <li><p><strong>origo</strong> (<em>numpy.array</em>) – The lower boundary of the domain [m]. Negative values
632 won’t work. Default=[0.0, 0.0, 0.0].</p></li>
633 <li><p><strong>dx</strong> (<em>float</em>) – The cell width in any direction. If the default value is used
634 (-1), the cell width is calculated to fit the largest particle.</p></li>
635 </ul>
636 </dd>
637 </dl>
638 </dd></dl>
639
640 <dl class="method">
641 <dt id="sphere.sim.deleteAllParticles">
642 <code class="sig-name descname">deleteAllParticles</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.deleteAllParticles" title="Permalink to this definition">¶</a></dt>
643 <dd><p>Deletes all particles in the simulation object.</p>
644 </dd></dl>
645
646 <dl class="method">
647 <dt id="sphere.sim.deleteParticle">
648 <code class="sig-name descname">deleteParticle</code><span class="sig-paren">(</span><em class="sig-param">i</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.deleteParticle" title="Permalink to this definition">¶</a></dt>
649 <dd><p>Delete particle(s) with index <code class="docutils literal notranslate"><span class="pre">i</span></code>.</p>
650 <dl class="field-list simple">
651 <dt class="field-odd">Parameters</dt>
652 <dd class="field-odd"><p><strong>i</strong> (<em>int</em><em>, </em><em>list</em><em> or </em><em>numpy.array</em>) – One or more particle indexes to delete</p>
653 </dd>
654 </dl>
655 </dd></dl>
656
657 <dl class="method">
658 <dt id="sphere.sim.disableFluidPressureModulation">
659 <code class="sig-name descname">disableFluidPressureModulation</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.disableFluidPressureModulation" title="Permalink to this definition">¶</a></dt>
660 <dd><p>Set the parameters for the sine wave modulating the fluid pressures
661 at the top boundary to zero.</p>
662 <p>See also: <a class="reference internal" href="#sphere.sim.setFluidPressureModulation" title="sphere.sim.setFluidPressureModulation"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidPressureModulation()</span></code></a></p>
663 </dd></dl>
664
665 <dl class="method">
666 <dt id="sphere.sim.disableTopWallNormalStressModulation">
667 <code class="sig-name descname">disableTopWallNormalStressModulation</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.disableTopWallNormalStressModulation" title="Permalink to this definition">¶</a></dt>
668 <dd><p>Set the parameters for the sine wave modulating the normal stress
669 at the top dynamic wall to zero.</p>
670 <p>See also: <a class="reference internal" href="#sphere.sim.setTopWallNormalStressModulation" title="sphere.sim.setTopWallNormalStressModulation"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTopWallNormalStressModulation()</span></code></a></p>
671 </dd></dl>
672
673 <dl class="method">
674 <dt id="sphere.sim.dry">
675 <code class="sig-name descname">dry</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.dry" title="Permalink to this definition">¶</a></dt>
676 <dd><p>Set the simulation to be dry (no fluids).</p>
677 <p>See also <a class="reference internal" href="#sphere.sim.wet" title="sphere.sim.wet"><code class="xref py py-func docutils literal notranslate"><span class="pre">wet()</span></code></a></p>
678 </dd></dl>
679
680 <dl class="method">
681 <dt id="sphere.sim.energy">
682 <code class="sig-name descname">energy</code><span class="sig-paren">(</span><em class="sig-param">method</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.energy" title="Permalink to this definition">¶</a></dt>
683 <dd><p>Calculates the sum of the energy components of all particles.</p>
684 <dl class="field-list simple">
685 <dt class="field-odd">Parameters</dt>
686 <dd class="field-odd"><p><strong>method</strong> (<em>str</em>) – The type of energy to return. Possible values are ‘pot’
687 for potential energy [J], ‘kin’ for kinetic energy [J], ‘rot’ for
688 rotational energy [J], ‘shear’ for energy lost by friction,
689 ‘shearrate’ for the rate of frictional energy loss [W], ‘visc_n’ for
690 viscous losses normal to the contact [J], ‘visc_n_rate’ for the rate
691 of viscous losses normal to the contact [W], and finally ‘bondpot’
692 for the potential energy stored in bonds [J]</p>
693 </dd>
694 <dt class="field-even">Returns</dt>
695 <dd class="field-even"><p>The value of the selected energy type</p>
696 </dd>
697 <dt class="field-odd">Return type</dt>
698 <dd class="field-odd"><p>float</p>
699 </dd>
700 </dl>
701 </dd></dl>
702
703 <dl class="method">
704 <dt id="sphere.sim.findAllAverageParticlePairAreas">
705 <code class="sig-name descname">findAllAverageParticlePairAreas</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findAllAverageParticlePairAreas" title="Permalink to this definition">¶</a></dt>
706 <dd><p>Finds the average area of an inter-particle contact. This
707 function requires a prior call to <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a> as it reads
708 from the <code class="docutils literal notranslate"><span class="pre">self.pairs</span></code> and <code class="docutils literal notranslate"><span class="pre">self.overlaps</span></code> arrays.</p>
709 <dl class="field-list simple">
710 <dt class="field-odd">Returns</dt>
711 <dd class="field-odd"><p>Array of contact surface areas</p>
712 </dd>
713 <dt class="field-even">Return type</dt>
714 <dd class="field-even"><p>array of floats</p>
715 </dd>
716 </dl>
717 </dd></dl>
718
719 <dl class="method">
720 <dt id="sphere.sim.findAllContactSurfaceAreas">
721 <code class="sig-name descname">findAllContactSurfaceAreas</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findAllContactSurfaceAreas" title="Permalink to this definition">¶</a></dt>
722 <dd><p>Finds the contact surface area of an inter-particle contact. This
723 function requires a prior call to <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a> as it reads
724 from the <code class="docutils literal notranslate"><span class="pre">self.pairs</span></code> and <code class="docutils literal notranslate"><span class="pre">self.overlaps</span></code> arrays.</p>
725 <dl class="field-list simple">
726 <dt class="field-odd">Returns</dt>
727 <dd class="field-odd"><p>Array of contact surface areas</p>
728 </dd>
729 <dt class="field-even">Return type</dt>
730 <dd class="field-even"><p>array of floats</p>
731 </dd>
732 </dl>
733 </dd></dl>
734
735 <dl class="method">
736 <dt id="sphere.sim.findContactStresses">
737 <code class="sig-name descname">findContactStresses</code><span class="sig-paren">(</span><em class="sig-param">area='average'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findContactStresses" title="Permalink to this definition">¶</a></dt>
738 <dd><p>Finds all particle-particle uniaxial normal stresses (by first calling
739 <a class="reference internal" href="#sphere.sim.findNormalForces" title="sphere.sim.findNormalForces"><code class="xref py py-func docutils literal notranslate"><span class="pre">findNormalForces()</span></code></a>) and calculating the stress magnitudes by
740 dividing the normal force magnitude with the average particle area
741 (‘average’) or by the contact surface area (‘contact’).</p>
742 <p>The result is saved in <code class="docutils literal notranslate"><span class="pre">self.sigma_contacts</span></code>.</p>
743 <dl class="field-list simple">
744 <dt class="field-odd">Parameters</dt>
745 <dd class="field-odd"><p><strong>area</strong> (<em>str</em>) – Area to use: ‘average’ (default) or ‘contact’</p>
746 </dd>
747 </dl>
748 <p>See also: <a class="reference internal" href="#sphere.sim.findNormalForces" title="sphere.sim.findNormalForces"><code class="xref py py-func docutils literal notranslate"><span class="pre">findNormalForces()</span></code></a> and <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a></p>
749 </dd></dl>
750
751 <dl class="method">
752 <dt id="sphere.sim.findCoordinationNumber">
753 <code class="sig-name descname">findCoordinationNumber</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findCoordinationNumber" title="Permalink to this definition">¶</a></dt>
754 <dd><p>Finds the coordination number (the average number of contacts per
755 particle). Requires a previous call to <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a>. Values
756 are stored in <code class="docutils literal notranslate"><span class="pre">self.coordinationnumber</span></code>.</p>
757 </dd></dl>
758
759 <dl class="method">
760 <dt id="sphere.sim.findHydraulicConductivities">
761 <code class="sig-name descname">findHydraulicConductivities</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findHydraulicConductivities" title="Permalink to this definition">¶</a></dt>
762 <dd><p>Calculates the hydrological conductivities from the Kozeny-Carman
763 relationship. These values are only relevant when the Darcy solver is
764 used (<cite>self.cfd_solver=1</cite>). The permeability pre-factor <cite>self.k_c</cite>
765 and the assemblage porosities must be set beforehand. The former values
766 are set if a file from the <cite>output/</cite> folder is read using
767 <cite>self.readbin</cite>.</p>
768 </dd></dl>
769
770 <dl class="method">
771 <dt id="sphere.sim.findLoadedContacts">
772 <code class="sig-name descname">findLoadedContacts</code><span class="sig-paren">(</span><em class="sig-param">threshold</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findLoadedContacts" title="Permalink to this definition">¶</a></dt>
773 <dd><p>Finds the indices of contact pairs where the contact stress magnitude
774 exceeds or is equal to a specified threshold value. This function calls
775 <a class="reference internal" href="#sphere.sim.findContactStresses" title="sphere.sim.findContactStresses"><code class="xref py py-func docutils literal notranslate"><span class="pre">findContactStresses()</span></code></a>.</p>
776 <dl class="field-list simple">
777 <dt class="field-odd">Parameters</dt>
778 <dd class="field-odd"><p><strong>threshold</strong> (<em>float</em>) – Threshold contact stress [Pa]</p>
779 </dd>
780 <dt class="field-even">Returns</dt>
781 <dd class="field-even"><p>Array of contact indices</p>
782 </dd>
783 <dt class="field-odd">Return type</dt>
784 <dd class="field-odd"><p>array of ints</p>
785 </dd>
786 </dl>
787 </dd></dl>
788
789 <dl class="method">
790 <dt id="sphere.sim.findMeanCoordinationNumber">
791 <code class="sig-name descname">findMeanCoordinationNumber</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findMeanCoordinationNumber" title="Permalink to this definition">¶</a></dt>
792 <dd><p>Returns the coordination number (the average number of contacts per
793 particle). Requires a previous call to <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a></p>
794 <dl class="field-list simple">
795 <dt class="field-odd">Returns</dt>
796 <dd class="field-odd"><p>The mean particle coordination number</p>
797 </dd>
798 <dt class="field-even">Return type</dt>
799 <dd class="field-even"><p>float</p>
800 </dd>
801 </dl>
802 </dd></dl>
803
804 <dl class="method">
805 <dt id="sphere.sim.findNormalForces">
806 <code class="sig-name descname">findNormalForces</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findNormalForces" title="Permalink to this definition">¶</a></dt>
807 <dd><p>Finds all particle-particle overlaps (by first calling
808 <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a>) and calculating the normal magnitude by
809 multiplying the overlaps with the elastic stiffness <code class="docutils literal notranslate"><span class="pre">self.k_n</span></code>.</p>
810 <p>The result is saved in <code class="docutils literal notranslate"><span class="pre">self.f_n_magn</span></code>.</p>
811 <p>See also: <a class="reference internal" href="#sphere.sim.findOverlaps" title="sphere.sim.findOverlaps"><code class="xref py py-func docutils literal notranslate"><span class="pre">findOverlaps()</span></code></a> and <a class="reference internal" href="#sphere.sim.findContactStresses" title="sphere.sim.findContactStresses"><code class="xref py py-func docutils literal notranslate"><span class="pre">findContactStresses()</span></code></a></p>
812 </dd></dl>
813
814 <dl class="method">
815 <dt id="sphere.sim.findOverlaps">
816 <code class="sig-name descname">findOverlaps</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findOverlaps" title="Permalink to this definition">¶</a></dt>
817 <dd><p>Find all particle-particle overlaps by a n^2 contact search, which is
818 done in C++. The particle pair indexes and the distance of the overlaps
819 is saved in the object itself as the <code class="docutils literal notranslate"><span class="pre">.pairs</span></code> and <code class="docutils literal notranslate"><span class="pre">.overlaps</span></code>
820 members.</p>
821 <p>See also: <a class="reference internal" href="#sphere.sim.findNormalForces" title="sphere.sim.findNormalForces"><code class="xref py py-func docutils literal notranslate"><span class="pre">findNormalForces()</span></code></a></p>
822 </dd></dl>
823
824 <dl class="method">
825 <dt id="sphere.sim.findPermeabilities">
826 <code class="sig-name descname">findPermeabilities</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.findPermeabilities" title="Permalink to this definition">¶</a></dt>
827 <dd><p>Calculates the hydrological permeabilities from the Kozeny-Carman
828 relationship. These values are only relevant when the Darcy solver is
829 used (<cite>self.cfd_solver=1</cite>). The permeability pre-factor <cite>self.k_c</cite>
830 and the assemblage porosities must be set beforehand. The former values
831 are set if a file from the <cite>output/</cite> folder is read using
832 <cite>self.readbin</cite>.</p>
833 </dd></dl>
834
835 <dl class="method">
836 <dt id="sphere.sim.forcechains">
837 <code class="sig-name descname">forcechains</code><span class="sig-paren">(</span><em class="sig-param">lc=200.0</em>, <em class="sig-param">uc=650.0</em>, <em class="sig-param">outformat='png'</em>, <em class="sig-param">disp='2d'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.forcechains" title="Permalink to this definition">¶</a></dt>
838 <dd><p>Visualizes the force chains in the system from the magnitude of the
839 normal contact forces, and produces an image of them. Warning: Will
840 segfault if no contacts are found.</p>
841 <dl class="field-list simple">
842 <dt class="field-odd">Parameters</dt>
843 <dd class="field-odd"><ul class="simple">
844 <li><p><strong>lc</strong> (<em>float</em>) – Lower cutoff of contact forces. Contacts below are not
845 visualized</p></li>
846 <li><p><strong>uc</strong> (<em>float</em>) – Upper cutoff of contact forces. Contacts above are
847 visualized with this value</p></li>
848 <li><p><strong>outformat</strong> (<em>str</em>) – Format of output image. Possible values are
849 ‘interactive’, ‘png’, ‘epslatex’, ‘epslatex-color’</p></li>
850 <li><p><strong>disp</strong> (<em>str</em>) – Display forcechains in ‘2d’ or ‘3d’</p></li>
851 </ul>
852 </dd>
853 </dl>
854 </dd></dl>
855
856 <dl class="method">
857 <dt id="sphere.sim.forcechainsRose">
858 <code class="sig-name descname">forcechainsRose</code><span class="sig-paren">(</span><em class="sig-param">lower_limit=0.25</em>, <em class="sig-param">graphics_format='pdf'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.forcechainsRose" title="Permalink to this definition">¶</a></dt>
859 <dd><p>Visualize trend and plunge angles of the strongest force chains in a
860 rose plot. The plots are saved in the current folder with the name
861 ‘fc-<simulation id>-rose.pdf’.</p>
862 <dl class="field-list simple">
863 <dt class="field-odd">Parameters</dt>
864 <dd class="field-odd"><ul class="simple">
865 <li><p><strong>lower_limit</strong> (<em>float</em>) – Do not visualize force chains below this relative
866 contact force magnitude, in ]0;1[</p></li>
867 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
868 </ul>
869 </dd>
870 </dl>
871 </dd></dl>
872
873 <dl class="method">
874 <dt id="sphere.sim.frictionalEnergy">
875 <code class="sig-name descname">frictionalEnergy</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.frictionalEnergy" title="Permalink to this definition">¶</a></dt>
876 <dd><p>Returns the frictional dissipated energy for a particle.</p>
877 <dl class="field-list simple">
878 <dt class="field-odd">Parameters</dt>
879 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
880 </dd>
881 <dt class="field-even">Returns</dt>
882 <dd class="field-even"><p>The frictional energy lost of the particle [J]</p>
883 </dd>
884 <dt class="field-odd">Return type</dt>
885 <dd class="field-odd"><p>float</p>
886 </dd>
887 </dl>
888 </dd></dl>
889
890 <dl class="method">
891 <dt id="sphere.sim.generateBimodalRadii">
892 <code class="sig-name descname">generateBimodalRadii</code><span class="sig-paren">(</span><em class="sig-param">r_small=0.005</em>, <em class="sig-param">r_large=0.05</em>, <em class="sig-param">ratio=0.2</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.generateBimodalRadii" title="Permalink to this definition">¶</a></dt>
893 <dd><p>Draw random radii from two distinct sizes.</p>
894 <dl class="field-list simple">
895 <dt class="field-odd">Parameters</dt>
896 <dd class="field-odd"><ul class="simple">
897 <li><p><strong>r_small</strong> (<em>float</em>) – Radii of small population [m], in ]0;r_large[</p></li>
898 <li><p><strong>r_large</strong> (<em>float</em>) – Radii of large population [m], in ]r_small;inf[</p></li>
899 <li><p><strong>ratio</strong> (<em>float</em>) – Approximate volumetric ratio between the two
900 populations (large/small).</p></li>
901 </ul>
902 </dd>
903 </dl>
904 <p>See also: <a class="reference internal" href="#sphere.sim.generateRadii" title="sphere.sim.generateRadii"><code class="xref py py-func docutils literal notranslate"><span class="pre">generateRadii()</span></code></a>.</p>
905 </dd></dl>
906
907 <dl class="method">
908 <dt id="sphere.sim.generateRadii">
909 <code class="sig-name descname">generateRadii</code><span class="sig-paren">(</span><em class="sig-param">psd='logn'</em>, <em class="sig-param">mean=0.00044</em>, <em class="sig-param">variance=8.8e-09</em>, <em class="sig-param">histogram=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.generateRadii" title="Permalink to this definition">¶</a></dt>
910 <dd><p>Draw random particle radii from a selected probability distribution.
911 The larger the variance of radii is, the slower the computations will
912 run. The reason is two-fold: The smallest particle dictates the time
913 step length, where smaller particles cause shorter time steps. At the
914 same time, the largest particle determines the sorting cell size, where
915 larger particles cause larger cells. Larger cells are likely to contain
916 more particles, causing more contact checks.</p>
917 <dl class="field-list simple">
918 <dt class="field-odd">Parameters</dt>
919 <dd class="field-odd"><ul class="simple">
920 <li><p><strong>psd</strong> (<em>str</em>) – The particle side distribution. One possible value is
921 <code class="docutils literal notranslate"><span class="pre">logn</span></code>, which is a log-normal probability distribution, suitable
922 for approximating well-sorted, coarse sediments. The other possible
923 value is <code class="docutils literal notranslate"><span class="pre">uni</span></code>, which is a uniform distribution from
924 <code class="docutils literal notranslate"><span class="pre">mean</span> <span class="pre">-</span> <span class="pre">variance</span></code> to <code class="docutils literal notranslate"><span class="pre">mean</span> <span class="pre">+</span> <span class="pre">variance</span></code>.</p></li>
925 <li><p><strong>mean</strong> (<em>float</em>) – The mean radius [m] (default=440e-6 m)</p></li>
926 <li><p><strong>variance</strong> (<em>float</em>) – The variance in the probability distribution
927 [m].</p></li>
928 </ul>
929 </dd>
930 </dl>
931 <p>See also: <a class="reference internal" href="#sphere.sim.generateBimodalRadii" title="sphere.sim.generateBimodalRadii"><code class="xref py py-func docutils literal notranslate"><span class="pre">generateBimodalRadii()</span></code></a>.</p>
932 </dd></dl>
933
934 <dl class="method">
935 <dt id="sphere.sim.hydraulicConductivity">
936 <code class="sig-name descname">hydraulicConductivity</code><span class="sig-paren">(</span><em class="sig-param">phi=0.35</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.hydraulicConductivity" title="Permalink to this definition">¶</a></dt>
937 <dd><p>Determine the hydraulic conductivity (K) [m/s] from the permeability
938 prefactor and a chosen porosity. This value is stored in <cite>self.K_c</cite>.
939 This function only works for the Darcy solver (<cite>self.cfd_solver == 1</cite>)</p>
940 <dl class="field-list simple">
941 <dt class="field-odd">Parameters</dt>
942 <dd class="field-odd"><p><strong>phi</strong> (<em>float</em>) – The porosity to use in the Kozeny-Carman relationship</p>
943 </dd>
944 <dt class="field-even">Returns</dt>
945 <dd class="field-even"><p>The hydraulic conductivity [m/s]</p>
946 </dd>
947 <dt class="field-odd">Return type</dt>
948 <dd class="field-odd"><p>float</p>
949 </dd>
950 </dl>
951 </dd></dl>
952
953 <dl class="method">
954 <dt id="sphere.sim.hydraulicDiffusivity">
955 <code class="sig-name descname">hydraulicDiffusivity</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.hydraulicDiffusivity" title="Permalink to this definition">¶</a></dt>
956 <dd><p>Determine the hydraulic diffusivity (D) [m*m/s]. The result is stored in
957 <cite>self.D</cite>. This function only works for the Darcy solver
958 (<cite>self.cfd_solver[0] == 1</cite>)</p>
959 </dd></dl>
960
961 <dl class="method">
962 <dt id="sphere.sim.hydraulicPermeability">
963 <code class="sig-name descname">hydraulicPermeability</code><span class="sig-paren">(</span><em class="sig-param">phi_min=0.3</em>, <em class="sig-param">phi_max=0.99</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.hydraulicPermeability" title="Permalink to this definition">¶</a></dt>
964 <dd><p>Determine the hydraulic permeability (k) [m*m] from the Kozeny-Carman
965 relationship, using the permeability prefactor (<cite>self.k_c</cite>), and the
966 range of valid porosities set in <cite>src/darcy.cuh</cite>, by default in the
967 range 0.1 to 0.9.</p>
968 <p>This function is only valid for the Darcy solver (<cite>self.cfd_solver ==
969 1</cite>).</p>
970 </dd></dl>
971
972 <dl class="method">
973 <dt id="sphere.sim.id">
974 <code class="sig-name descname">id</code><span class="sig-paren">(</span><em class="sig-param">sid=''</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.id" title="Permalink to this definition">¶</a></dt>
975 <dd><p>Returns or sets the simulation id/name, which is used to identify
976 simulation files in the output folders.</p>
977 <dl class="field-list simple">
978 <dt class="field-odd">Parameters</dt>
979 <dd class="field-odd"><p><strong>sid</strong> (<em>str</em>) – The desired simulation id. If left blank the current
980 simulation id will be returned.</p>
981 </dd>
982 <dt class="field-even">Returns</dt>
983 <dd class="field-even"><p>The current simulation id if no new value is set.</p>
984 </dd>
985 <dt class="field-odd">Return type</dt>
986 <dd class="field-odd"><p>str</p>
987 </dd>
988 </dl>
989 </dd></dl>
990
991 <dl class="method">
992 <dt id="sphere.sim.idAppend">
993 <code class="sig-name descname">idAppend</code><span class="sig-paren">(</span><em class="sig-param">string</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.idAppend" title="Permalink to this definition">¶</a></dt>
994 <dd><p>Append a string to the simulation id/name, which is used to identify
995 simulation files in the output folders.</p>
996 <dl class="field-list simple">
997 <dt class="field-odd">Parameters</dt>
998 <dd class="field-odd"><p><strong>string</strong> (<em>str</em>) – The string to append to the simulation id (<cite>self.sid</cite>).</p>
999 </dd>
1000 </dl>
1001 </dd></dl>
1002
1003 <dl class="method">
1004 <dt id="sphere.sim.inertiaParameterPlanarShear">
1005 <code class="sig-name descname">inertiaParameterPlanarShear</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.inertiaParameterPlanarShear" title="Permalink to this definition">¶</a></dt>
1006 <dd><p>Returns the value of the inertia parameter $I$ during planar shear
1007 proposed by GDR-MiDi 2004.</p>
1008 <dl class="field-list simple">
1009 <dt class="field-odd">Returns</dt>
1010 <dd class="field-odd"><p>The value of $I$</p>
1011 </dd>
1012 <dt class="field-even">Return type</dt>
1013 <dd class="field-even"><p>float</p>
1014 </dd>
1015 </dl>
1016 <p>See also: <a class="reference internal" href="#sphere.sim.shearStrainRate" title="sphere.sim.shearStrainRate"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearStrainRate()</span></code></a> and <a class="reference internal" href="#sphere.sim.shearVel" title="sphere.sim.shearVel"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearVel()</span></code></a></p>
1017 </dd></dl>
1018
1019 <dl class="method">
1020 <dt id="sphere.sim.initFluid">
1021 <code class="sig-name descname">initFluid</code><span class="sig-paren">(</span><em class="sig-param">mu=0.00089</em>, <em class="sig-param">rho=1000.0</em>, <em class="sig-param">p=0.0</em>, <em class="sig-param">hydrostatic=False</em>, <em class="sig-param">cfd_solver=0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initFluid" title="Permalink to this definition">¶</a></dt>
1022 <dd><p>Initialize the fluid arrays and the fluid viscosity. The default value
1023 of <code class="docutils literal notranslate"><span class="pre">mu</span></code> equals the dynamic viscosity of water at 25 degrees Celcius.
1024 The value for water at 0 degrees Celcius is 17.87e-4 kg/(m*s).</p>
1025 <dl class="field-list simple">
1026 <dt class="field-odd">Parameters</dt>
1027 <dd class="field-odd"><ul class="simple">
1028 <li><p><strong>mu</strong> (<em>float</em>) – The fluid dynamic viscosity [kg/(m*s)]</p></li>
1029 <li><p><strong>rho</strong> (<em>float</em>) – The fluid density [kg/(m^3)]</p></li>
1030 <li><p><strong>p</strong> – The hydraulic pressure to initialize the cells to. If the
1031 parameter <cite>hydrostatic</cite> is set to <cite>True</cite>, this value will apply to
1032 the fluid cells at the top</p></li>
1033 <li><p><strong>hydrostatic</strong> (<em>bool</em>) – Initialize the fluid pressures to the hydrostatic
1034 pressure distribution. A pressure gradient with depth is only
1035 created if a gravitational acceleration along <img class="math" src="_images/math/8d051150f8669295ecdbe92367941012175a824d.png" alt="z"/> previously
1036 has been specified</p></li>
1037 <li><p><strong>cfd_solver</strong> (<em>int</em>) – Solver to use for the computational fluid dynamics.
1038 Accepted values: 0 (Navier Stokes, default) and 1 (Darcy).</p></li>
1039 </ul>
1040 </dd>
1041 </dl>
1042 </dd></dl>
1043
1044 <dl class="method">
1045 <dt id="sphere.sim.initGrid">
1046 <code class="sig-name descname">initGrid</code><span class="sig-paren">(</span><em class="sig-param">dx=-1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initGrid" title="Permalink to this definition">¶</a></dt>
1047 <dd><p>Initialize grid suitable for the particle positions set previously.
1048 The margin parameter adjusts the distance (in no. of max. radii)
1049 from the particle boundaries.
1050 <em>Important</em>: The particle radii have to be set beforehand if the cell
1051 width isn’t specified by <cite>dx</cite>.</p>
1052 <dl class="field-list simple">
1053 <dt class="field-odd">Parameters</dt>
1054 <dd class="field-odd"><p><strong>dx</strong> (<em>float</em>) – The cell width in any direction. If the default value is used
1055 (-1), the cell width is calculated to fit the largest particle.</p>
1056 </dd>
1057 </dl>
1058 </dd></dl>
1059
1060 <dl class="method">
1061 <dt id="sphere.sim.initGridAndWorldsize">
1062 <code class="sig-name descname">initGridAndWorldsize</code><span class="sig-paren">(</span><em class="sig-param">margin=2.0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initGridAndWorldsize" title="Permalink to this definition">¶</a></dt>
1063 <dd><p>Initialize grid suitable for the particle positions set previously.
1064 The margin parameter adjusts the distance (in no. of max. radii)
1065 from the particle boundaries. If the upper wall is dynamic, it is placed
1066 at the top boundary of the world.</p>
1067 <dl class="field-list simple">
1068 <dt class="field-odd">Parameters</dt>
1069 <dd class="field-odd"><p><strong>margin</strong> (<em>float</em>) – Distance to world boundary in no. of max. particle radii</p>
1070 </dd>
1071 </dl>
1072 </dd></dl>
1073
1074 <dl class="method">
1075 <dt id="sphere.sim.initGridPos">
1076 <code class="sig-name descname">initGridPos</code><span class="sig-paren">(</span><em class="sig-param">gridnum=array([12</em>, <em class="sig-param">12</em>, <em class="sig-param">36])</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initGridPos" title="Permalink to this definition">¶</a></dt>
1077 <dd><p>Initialize particle positions in loose, cubic configuration.
1078 <code class="docutils literal notranslate"><span class="pre">gridnum</span></code> is the number of cells in the x, y and z directions.
1079 <em>Important</em>: The particle radii and the boundary conditions (periodic or
1080 not) for the x and y boundaries have to be set beforehand.</p>
1081 <dl class="field-list simple">
1082 <dt class="field-odd">Parameters</dt>
1083 <dd class="field-odd"><p><strong>gridnum</strong> (<em>numpy.array</em>) – The number of particles in x, y and z directions</p>
1084 </dd>
1085 </dl>
1086 </dd></dl>
1087
1088 <dl class="method">
1089 <dt id="sphere.sim.initRandomGridPos">
1090 <code class="sig-name descname">initRandomGridPos</code><span class="sig-paren">(</span><em class="sig-param">gridnum=array([12</em>, <em class="sig-param">12</em>, <em class="sig-param">32])</em>, <em class="sig-param">padding=2.1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initRandomGridPos" title="Permalink to this definition">¶</a></dt>
1091 <dd><p>Initialize particle positions in loose, cubic configuration with some
1092 variance. <code class="docutils literal notranslate"><span class="pre">gridnum</span></code> is the number of cells in the x, y and z
1093 directions. <em>Important</em>: The particle radii and the boundary conditions
1094 (periodic or not) for the x and y boundaries have to be set beforehand.
1095 The world size and grid height (in the z direction) is readjusted to fit
1096 the particle positions.</p>
1097 <dl class="field-list simple">
1098 <dt class="field-odd">Parameters</dt>
1099 <dd class="field-odd"><ul class="simple">
1100 <li><p><strong>gridnum</strong> (<em>numpy.array</em>) – The number of particles in x, y and z directions</p></li>
1101 <li><p><strong>padding</strong> (<em>float</em>) – Increase distance between particles in x, y and z
1102 directions with this multiplier. Large values create more random
1103 packings.</p></li>
1104 </ul>
1105 </dd>
1106 </dl>
1107 </dd></dl>
1108
1109 <dl class="method">
1110 <dt id="sphere.sim.initRandomPos">
1111 <code class="sig-name descname">initRandomPos</code><span class="sig-paren">(</span><em class="sig-param">gridnum=array([12</em>, <em class="sig-param">12</em>, <em class="sig-param">36])</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initRandomPos" title="Permalink to this definition">¶</a></dt>
1112 <dd><p>Initialize particle positions in completely random configuration. Radii
1113 <em>must</em> be set beforehand. If the x and y boundaries are set as periodic,
1114 the particle centers will be placed all the way to the edge. On regular,
1115 non-periodic boundaries, the particles are restrained at the edges to
1116 make space for their radii within the bounding box.</p>
1117 <dl class="field-list simple">
1118 <dt class="field-odd">Parameters</dt>
1119 <dd class="field-odd"><ul class="simple">
1120 <li><p><strong>gridnum</strong> (<em>numpy.array</em>) – The number of sorting cells in each spatial direction
1121 (default=[12, 12, 36])</p></li>
1122 <li><p><strong>dx</strong> (<em>float</em>) – The cell width in any direction. If the default value is used
1123 (-1), the cell width is calculated to fit the largest particle.</p></li>
1124 </ul>
1125 </dd>
1126 </dl>
1127 </dd></dl>
1128
1129 <dl class="method">
1130 <dt id="sphere.sim.initTemporal">
1131 <code class="sig-name descname">initTemporal</code><span class="sig-paren">(</span><em class="sig-param">total</em>, <em class="sig-param">current=0.0</em>, <em class="sig-param">file_dt=0.05</em>, <em class="sig-param">step_count=0</em>, <em class="sig-param">dt=-1</em>, <em class="sig-param">epsilon=0.01</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.initTemporal" title="Permalink to this definition">¶</a></dt>
1132 <dd><p>Set temporal parameters for the simulation. <em>Important</em>: Particle radii,
1133 physical parameters, and the optional fluid grid need to be set prior to
1134 these if the computational time step (dt) isn’t set explicitly. If the
1135 parameter <cite>dt</cite> is the default value (-1), the function will estimate the
1136 best time step length. The value of the computational time step for the
1137 DEM is checked for stability in the CFD solution if fluid simulation is
1138 included.</p>
1139 <dl class="field-list simple">
1140 <dt class="field-odd">Parameters</dt>
1141 <dd class="field-odd"><ul class="simple">
1142 <li><p><strong>total</strong> (<em>float</em>) – The time at which to end the simulation [s]</p></li>
1143 <li><p><strong>current</strong> – The current time [s] (default=0.0 s)</p></li>
1144 <li><p><strong>file_dt</strong> – The interval between output files [s] (default=0.05 s)</p></li>
1145 <li><p><strong>dt</strong> – The computational time step length [s]</p></li>
1146 <li><p><strong>epsilon</strong> (<em>float</em>) – Time step multiplier (default=0.01)</p></li>
1147 </ul>
1148 </dd>
1149 <dt class="field-even">Step_count</dt>
1150 <dd class="field-even"><p>The number of the first output file (default=0)</p>
1151 </dd>
1152 </dl>
1153 </dd></dl>
1154
1155 <dl class="method">
1156 <dt id="sphere.sim.kineticEnergy">
1157 <code class="sig-name descname">kineticEnergy</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.kineticEnergy" title="Permalink to this definition">¶</a></dt>
1158 <dd><p>Returns the (linear) kinetic energy for a particle.</p>
1159 <dl class="field-list simple">
1160 <dt class="field-odd">Parameters</dt>
1161 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
1162 </dd>
1163 <dt class="field-even">Returns</dt>
1164 <dd class="field-even"><p>The kinetic energy of the particle [J]</p>
1165 </dd>
1166 <dt class="field-odd">Return type</dt>
1167 <dd class="field-odd"><p>float</p>
1168 </dd>
1169 </dl>
1170 </dd></dl>
1171
1172 <dl class="method">
1173 <dt id="sphere.sim.largestFluidTimeStep">
1174 <code class="sig-name descname">largestFluidTimeStep</code><span class="sig-paren">(</span><em class="sig-param">safety=0.5</em>, <em class="sig-param">v_max=-1.0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.largestFluidTimeStep" title="Permalink to this definition">¶</a></dt>
1175 <dd><p>Finds and returns the largest time step in the fluid phase by von
1176 Neumann and Courant-Friedrichs-Lewy analysis given the current
1177 velocities. This ensures stability in the diffusive and advective parts
1178 of the momentum equation.</p>
1179 <p>The value of the time step decreases with increasing fluid viscosity
1180 (<cite>self.mu</cite>), and increases with fluid cell size (<cite>self.L/self.num</cite>)</p>
1181 <p>and fluid velocities (<cite>self.v_f</cite>)</p>
1182 <dl class="field-list simple">
1183 <dt class="field-odd">Parameters</dt>
1184 <dd class="field-odd"><ul class="simple">
1185 <li><p><strong>safety</strong> (<em>float</em>) – Safety factor which is multiplied to the largest time
1186 step.</p></li>
1187 <li><p><strong>v_max</strong> (<em>float</em>) – The largest anticipated absolute fluid velocity [m/s]</p></li>
1188 </ul>
1189 </dd>
1190 <dt class="field-even">Returns</dt>
1191 <dd class="field-even"><p>The largest timestep stable for the current fluid state.</p>
1192 </dd>
1193 <dt class="field-odd">Return type</dt>
1194 <dd class="field-odd"><p>float</p>
1195 </dd>
1196 </dl>
1197 </dd></dl>
1198
1199 <dl class="method">
1200 <dt id="sphere.sim.largestMass">
1201 <code class="sig-name descname">largestMass</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.largestMass" title="Permalink to this definition">¶</a></dt>
1202 <dd><p>Returns the mass of the heaviest particle.</p>
1203 <dl class="field-list simple">
1204 <dt class="field-odd">Parameters</dt>
1205 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
1206 </dd>
1207 <dt class="field-even">Returns</dt>
1208 <dd class="field-even"><p>The mass of the particle [kg]</p>
1209 </dd>
1210 <dt class="field-odd">Return type</dt>
1211 <dd class="field-odd"><p>float</p>
1212 </dd>
1213 </dl>
1214 </dd></dl>
1215
1216 <dl class="method">
1217 <dt id="sphere.sim.mass">
1218 <code class="sig-name descname">mass</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.mass" title="Permalink to this definition">¶</a></dt>
1219 <dd><p>Returns the mass of a particle.</p>
1220 <dl class="field-list simple">
1221 <dt class="field-odd">Parameters</dt>
1222 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
1223 </dd>
1224 <dt class="field-even">Returns</dt>
1225 <dd class="field-even"><p>The mass of the particle [kg]</p>
1226 </dd>
1227 <dt class="field-odd">Return type</dt>
1228 <dd class="field-odd"><p>float</p>
1229 </dd>
1230 </dl>
1231 </dd></dl>
1232
1233 <dl class="method">
1234 <dt id="sphere.sim.momentOfInertia">
1235 <code class="sig-name descname">momentOfInertia</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.momentOfInertia" title="Permalink to this definition">¶</a></dt>
1236 <dd><p>Returns the moment of inertia of a particle.</p>
1237 <dl class="field-list simple">
1238 <dt class="field-odd">Parameters</dt>
1239 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
1240 </dd>
1241 <dt class="field-even">Returns</dt>
1242 <dd class="field-even"><p>The moment of inertia [kg*m^2]</p>
1243 </dd>
1244 <dt class="field-odd">Return type</dt>
1245 <dd class="field-odd"><p>float</p>
1246 </dd>
1247 </dl>
1248 </dd></dl>
1249
1250 <dl class="method">
1251 <dt id="sphere.sim.momentum">
1252 <code class="sig-name descname">momentum</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.momentum" title="Permalink to this definition">¶</a></dt>
1253 <dd><p>Returns the momentum (m*v) of a particle.</p>
1254 <dl class="field-list simple">
1255 <dt class="field-odd">Parameters</dt>
1256 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – The particle index</p>
1257 </dd>
1258 <dt class="field-even">Returns</dt>
1259 <dd class="field-even"><p>The particle momentum [N*s]</p>
1260 </dd>
1261 <dt class="field-odd">Return type</dt>
1262 <dd class="field-odd"><p>numpy.array</p>
1263 </dd>
1264 </dl>
1265 </dd></dl>
1266
1267 <dl class="method">
1268 <dt id="sphere.sim.normalBoundariesXY">
1269 <code class="sig-name descname">normalBoundariesXY</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.normalBoundariesXY" title="Permalink to this definition">¶</a></dt>
1270 <dd><p>Set the x and y boundary conditions to be static walls.</p>
1271 <p>See also <a class="reference internal" href="#sphere.sim.periodicBoundariesXY" title="sphere.sim.periodicBoundariesXY"><code class="xref py py-func docutils literal notranslate"><span class="pre">periodicBoundariesXY()</span></code></a> and
1272 <a class="reference internal" href="#sphere.sim.periodicBoundariesX" title="sphere.sim.periodicBoundariesX"><code class="xref py py-func docutils literal notranslate"><span class="pre">periodicBoundariesX()</span></code></a></p>
1273 </dd></dl>
1274
1275 <dl class="method">
1276 <dt id="sphere.sim.periodicBoundariesX">
1277 <code class="sig-name descname">periodicBoundariesX</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.periodicBoundariesX" title="Permalink to this definition">¶</a></dt>
1278 <dd><p>Set the x boundary conditions to be periodic.</p>
1279 <p>See also <a class="reference internal" href="#sphere.sim.normalBoundariesXY" title="sphere.sim.normalBoundariesXY"><code class="xref py py-func docutils literal notranslate"><span class="pre">normalBoundariesXY()</span></code></a> and
1280 <a class="reference internal" href="#sphere.sim.periodicBoundariesXY" title="sphere.sim.periodicBoundariesXY"><code class="xref py py-func docutils literal notranslate"><span class="pre">periodicBoundariesXY()</span></code></a></p>
1281 </dd></dl>
1282
1283 <dl class="method">
1284 <dt id="sphere.sim.periodicBoundariesXY">
1285 <code class="sig-name descname">periodicBoundariesXY</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.periodicBoundariesXY" title="Permalink to this definition">¶</a></dt>
1286 <dd><p>Set the x and y boundary conditions to be periodic.</p>
1287 <p>See also <a class="reference internal" href="#sphere.sim.normalBoundariesXY" title="sphere.sim.normalBoundariesXY"><code class="xref py py-func docutils literal notranslate"><span class="pre">normalBoundariesXY()</span></code></a> and
1288 <a class="reference internal" href="#sphere.sim.periodicBoundariesX" title="sphere.sim.periodicBoundariesX"><code class="xref py py-func docutils literal notranslate"><span class="pre">periodicBoundariesX()</span></code></a></p>
1289 </dd></dl>
1290
1291 <dl class="method">
1292 <dt id="sphere.sim.plotContacts">
1293 <code class="sig-name descname">plotContacts</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='png', figsize=[4, 4], title=None, lower_limit=0.0, upper_limit=1.0, alpha=1.0, return_data=False, outfolder='.', f_min=None, f_max=None, histogram=True, forcechains=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotContacts" title="Permalink to this definition">¶</a></dt>
1294 <dd><p>Plot current contact orientations on polar plot</p>
1295 <dl class="field-list simple">
1296 <dt class="field-odd">Parameters</dt>
1297 <dd class="field-odd"><ul class="simple">
1298 <li><p><strong>lower_limit</strong> (<em>float</em>) – Do not visualize force chains below this relative
1299 contact force magnitude, in ]0;1[</p></li>
1300 <li><p><strong>upper_limit</strong> (<em>float</em>) – Visualize force chains above this relative
1301 contact force magnitude but cap color bar range, in ]0;1[</p></li>
1302 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1303 </ul>
1304 </dd>
1305 </dl>
1306 </dd></dl>
1307
1308 <dl class="method">
1309 <dt id="sphere.sim.plotConvergence">
1310 <code class="sig-name descname">plotConvergence</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotConvergence" title="Permalink to this definition">¶</a></dt>
1311 <dd><p>Plot the convergence evolution in the CFD solver. The plot is saved
1312 in the output folder with the file name
1313 ‘<simulation id>-conv.<graphics_format>’.</p>
1314 <dl class="field-list simple">
1315 <dt class="field-odd">Parameters</dt>
1316 <dd class="field-odd"><ul class="simple">
1317 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1318 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1319 </ul>
1320 </dd>
1321 </dl>
1322 <p>See also: <a class="reference internal" href="#sphere.sim.convergence" title="sphere.sim.convergence"><code class="xref py py-func docutils literal notranslate"><span class="pre">convergence()</span></code></a></p>
1323 </dd></dl>
1324
1325 <dl class="method">
1326 <dt id="sphere.sim.plotFluidDiffAdvPresZ">
1327 <code class="sig-name descname">plotFluidDiffAdvPresZ</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotFluidDiffAdvPresZ" title="Permalink to this definition">¶</a></dt>
1328 <dd><p>Compare contributions to the velocity from diffusion and advection,
1329 assuming the flow is 1D along the z-axis, phi=1, and dphi=0. This
1330 solution is analog to the predicted velocity and not constrained by the
1331 conservation of mass. The plot is saved in the output folder with the
1332 name format ‘<simulation id>-diff_adv-t=<current time>s-mu=<dynamic
1333 viscosity>Pa-s.<graphics_format>’.</p>
1334 <dl class="field-list simple">
1335 <dt class="field-odd">Parameters</dt>
1336 <dd class="field-odd"><ul class="simple">
1337 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1338 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1339 </ul>
1340 </dd>
1341 </dl>
1342 </dd></dl>
1343
1344 <dl class="method">
1345 <dt id="sphere.sim.plotFluidPressuresY">
1346 <code class="sig-name descname">plotFluidPressuresY</code><span class="sig-paren">(</span><em class="sig-param">y=-1</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotFluidPressuresY" title="Permalink to this definition">¶</a></dt>
1347 <dd><p>Plot fluid pressures in a plane normal to the second axis.
1348 The plot is saved in the current folder with the format
1349 ‘p_f-<simulation id>-y<y value>.<graphics_format>’.</p>
1350 <dl class="field-list simple">
1351 <dt class="field-odd">Parameters</dt>
1352 <dd class="field-odd"><ul class="simple">
1353 <li><p><strong>y</strong> (<em>int</em>) – Plot pressures in fluid cells with these y axis values. If
1354 this value is -1, the center y position is used.</p></li>
1355 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1356 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1357 </ul>
1358 </dd>
1359 </dl>
1360 <p>See also: <a class="reference internal" href="#sphere.sim.writeFluidVTK" title="sphere.sim.writeFluidVTK"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeFluidVTK()</span></code></a> and <a class="reference internal" href="#sphere.sim.plotFluidPressuresZ" title="sphere.sim.plotFluidPressuresZ"><code class="xref py py-func docutils literal notranslate"><span class="pre">plotFluidPressuresZ()</span></code></a></p>
1361 </dd></dl>
1362
1363 <dl class="method">
1364 <dt id="sphere.sim.plotFluidPressuresZ">
1365 <code class="sig-name descname">plotFluidPressuresZ</code><span class="sig-paren">(</span><em class="sig-param">z=-1</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotFluidPressuresZ" title="Permalink to this definition">¶</a></dt>
1366 <dd><p>Plot fluid pressures in a plane normal to the third axis.
1367 The plot is saved in the current folder with the format
1368 ‘p_f-<simulation id>-z<z value>.<graphics_format>’.</p>
1369 <dl class="field-list simple">
1370 <dt class="field-odd">Parameters</dt>
1371 <dd class="field-odd"><ul class="simple">
1372 <li><p><strong>z</strong> (<em>int</em>) – Plot pressures in fluid cells with these z axis values. If
1373 this value is -1, the center z position is used.</p></li>
1374 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1375 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1376 </ul>
1377 </dd>
1378 </dl>
1379 <p>See also: <a class="reference internal" href="#sphere.sim.writeFluidVTK" title="sphere.sim.writeFluidVTK"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeFluidVTK()</span></code></a> and <a class="reference internal" href="#sphere.sim.plotFluidPressuresY" title="sphere.sim.plotFluidPressuresY"><code class="xref py py-func docutils literal notranslate"><span class="pre">plotFluidPressuresY()</span></code></a></p>
1380 </dd></dl>
1381
1382 <dl class="method">
1383 <dt id="sphere.sim.plotFluidVelocitiesY">
1384 <code class="sig-name descname">plotFluidVelocitiesY</code><span class="sig-paren">(</span><em class="sig-param">y=-1</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotFluidVelocitiesY" title="Permalink to this definition">¶</a></dt>
1385 <dd><p>Plot fluid velocities in a plane normal to the second axis.
1386 The plot is saved in the current folder with the format
1387 ‘v_f-<simulation id>-z<z value>.<graphics_format>’.</p>
1388 <dl class="field-list simple">
1389 <dt class="field-odd">Parameters</dt>
1390 <dd class="field-odd"><ul class="simple">
1391 <li><p><strong>y</strong> (<em>int</em>) – Plot velocities in fluid cells with these y axis values. If
1392 this value is -1, the center y position is used.</p></li>
1393 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1394 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1395 </ul>
1396 </dd>
1397 </dl>
1398 <p>See also: <a class="reference internal" href="#sphere.sim.writeFluidVTK" title="sphere.sim.writeFluidVTK"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeFluidVTK()</span></code></a> and <a class="reference internal" href="#sphere.sim.plotFluidVelocitiesZ" title="sphere.sim.plotFluidVelocitiesZ"><code class="xref py py-func docutils literal notranslate"><span class="pre">plotFluidVelocitiesZ()</span></code></a></p>
1399 </dd></dl>
1400
1401 <dl class="method">
1402 <dt id="sphere.sim.plotFluidVelocitiesZ">
1403 <code class="sig-name descname">plotFluidVelocitiesZ</code><span class="sig-paren">(</span><em class="sig-param">z=-1</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotFluidVelocitiesZ" title="Permalink to this definition">¶</a></dt>
1404 <dd><p>Plot fluid velocities in a plane normal to the third axis.
1405 The plot is saved in the current folder with the format
1406 ‘v_f-<simulation id>-z<z value>.<graphics_format>’.</p>
1407 <dl class="field-list simple">
1408 <dt class="field-odd">Parameters</dt>
1409 <dd class="field-odd"><ul class="simple">
1410 <li><p><strong>z</strong> (<em>int</em>) – Plot velocities in fluid cells with these z axis values. If
1411 this value is -1, the center z position is used.</p></li>
1412 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1413 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1414 </ul>
1415 </dd>
1416 </dl>
1417 <p>See also: <a class="reference internal" href="#sphere.sim.writeFluidVTK" title="sphere.sim.writeFluidVTK"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeFluidVTK()</span></code></a> and <a class="reference internal" href="#sphere.sim.plotFluidVelocitiesY" title="sphere.sim.plotFluidVelocitiesY"><code class="xref py py-func docutils literal notranslate"><span class="pre">plotFluidVelocitiesY()</span></code></a></p>
1418 </dd></dl>
1419
1420 <dl class="method">
1421 <dt id="sphere.sim.plotLoadCurve">
1422 <code class="sig-name descname">plotLoadCurve</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotLoadCurve" title="Permalink to this definition">¶</a></dt>
1423 <dd><p>Plot the load curve (log time vs. upper wall movement). The plot is
1424 saved in the current folder with the file name
1425 ‘<simulation id>-loadcurve.<graphics_format>’.
1426 The consolidation coefficient calculations are done on the base of
1427 Bowles 1992, p. 129–139, using the “Casagrande” method.
1428 It is assumed that the consolidation has stopped at the end of the
1429 simulation (i.e. flat curve).</p>
1430 <dl class="field-list simple">
1431 <dt class="field-odd">Parameters</dt>
1432 <dd class="field-odd"><ul class="simple">
1433 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1434 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1435 </ul>
1436 </dd>
1437 </dl>
1438 </dd></dl>
1439
1440 <dl class="method">
1441 <dt id="sphere.sim.plotPrescribedFluidPressures">
1442 <code class="sig-name descname">plotPrescribedFluidPressures</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotPrescribedFluidPressures" title="Permalink to this definition">¶</a></dt>
1443 <dd><p>Plot the prescribed fluid pressures through time that may be
1444 modulated through the class parameters p_mod_A, p_mod_f, and p_mod_phi.
1445 The plot is saved in the output folder with the file name
1446 ‘<simulation id>-pres.<graphics_format>’.</p>
1447 </dd></dl>
1448
1449 <dl class="method">
1450 <dt id="sphere.sim.plotSinFunction">
1451 <code class="sig-name descname">plotSinFunction</code><span class="sig-paren">(</span><em class="sig-param">baseval</em>, <em class="sig-param">A</em>, <em class="sig-param">f</em>, <em class="sig-param">phi=0.0</em>, <em class="sig-param">xlabel='$t$ [s]'</em>, <em class="sig-param">ylabel='$y$'</em>, <em class="sig-param">plotstyle='.'</em>, <em class="sig-param">outformat='png'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.plotSinFunction" title="Permalink to this definition">¶</a></dt>
1452 <dd><p>Plot the values of a sinusoidal modulated base value. Saves the output
1453 as a plot in the current folder.
1454 The time values will range from <cite>self.time_current</cite> to
1455 <cite>self.time_total</cite>.</p>
1456 <dl class="field-list simple">
1457 <dt class="field-odd">Parameters</dt>
1458 <dd class="field-odd"><ul class="simple">
1459 <li><p><strong>baseval</strong> (<em>float</em>) – The center value which the sinusoidal fluctuations are
1460 modulating</p></li>
1461 <li><p><strong>A</strong> (<em>float</em>) – The fluctuation amplitude</p></li>
1462 <li><p><strong>phi</strong> (<em>float</em>) – The phase shift [s]</p></li>
1463 <li><p><strong>xlabel</strong> (<em>str</em>) – The label for the x axis</p></li>
1464 <li><p><strong>ylabel</strong> (<em>str</em>) – The label for the y axis</p></li>
1465 <li><p><strong>plotstyle</strong> (<em>str</em>) – Matplotlib-string for specifying plotting style</p></li>
1466 <li><p><strong>outformat</strong> (<em>str</em>) – File format of the output plot</p></li>
1467 <li><p><strong>verbose</strong> (<em>bool</em>) – Print output filename after saving</p></li>
1468 </ul>
1469 </dd>
1470 </dl>
1471 </dd></dl>
1472
1473 <dl class="method">
1474 <dt id="sphere.sim.porosities">
1475 <code class="sig-name descname">porosities</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='pdf'</em>, <em class="sig-param">zslices=16</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.porosities" title="Permalink to this definition">¶</a></dt>
1476 <dd><p>Plot the averaged porosities with depth. The plot is saved in the format
1477 ‘<simulation id>-porosity.<graphics_format>’.</p>
1478 <dl class="field-list simple">
1479 <dt class="field-odd">Parameters</dt>
1480 <dd class="field-odd"><ul class="simple">
1481 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p></li>
1482 <li><p><strong>zslices</strong> (<em>int</em>) – The number of points along the vertical axis to sample
1483 the porosity in</p></li>
1484 </ul>
1485 </dd>
1486 </dl>
1487 </dd></dl>
1488
1489 <dl class="method">
1490 <dt id="sphere.sim.porosity">
1491 <code class="sig-name descname">porosity</code><span class="sig-paren">(</span><em class="sig-param">slices=10</em>, <em class="sig-param">verbose=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.porosity" title="Permalink to this definition">¶</a></dt>
1492 <dd><p>Calculates the porosity as a function of depth, by averaging values in
1493 horizontal slabs. Returns porosity values and their corresponding depth.
1494 The values are calculated using the external <code class="docutils literal notranslate"><span class="pre">porosity</span></code> program.</p>
1495 <dl class="field-list simple">
1496 <dt class="field-odd">Parameters</dt>
1497 <dd class="field-odd"><ul class="simple">
1498 <li><p><strong>slices</strong> (<em>int</em>) – The number of vertical slabs to find porosities in.</p></li>
1499 <li><p><strong>verbose</strong> (<em>bool</em>) – Show the file name of the temporary file written to
1500 disk</p></li>
1501 </ul>
1502 </dd>
1503 <dt class="field-even">Returns</dt>
1504 <dd class="field-even"><p>A 2d array of depths and their averaged porosities</p>
1505 </dd>
1506 <dt class="field-odd">Return type</dt>
1507 <dd class="field-odd"><p>numpy.array</p>
1508 </dd>
1509 </dl>
1510 </dd></dl>
1511
1512 <dl class="method">
1513 <dt id="sphere.sim.randomBondPairs">
1514 <code class="sig-name descname">randomBondPairs</code><span class="sig-paren">(</span><em class="sig-param">ratio=0.3</em>, <em class="sig-param">spacing=-0.1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.randomBondPairs" title="Permalink to this definition">¶</a></dt>
1515 <dd><p>Bond an amount of particles in two-particle clusters. The particles
1516 should be initialized beforehand. Note: The actual number of bonds is
1517 likely to be somewhat smaller than specified, due to the random
1518 selection algorithm.</p>
1519 <dl class="field-list simple">
1520 <dt class="field-odd">Parameters</dt>
1521 <dd class="field-odd"><ul class="simple">
1522 <li><p><strong>ratio</strong> (<em>float</em>) – The amount of particles to bond, values in ]0.0;1.0]</p></li>
1523 <li><p><strong>spacing</strong> (<em>float</em>) – The distance relative to the sum of radii between bonded
1524 particles, neg. values denote an overlap. Values in ]0.0,inf[.</p></li>
1525 </ul>
1526 </dd>
1527 </dl>
1528 </dd></dl>
1529
1530 <dl class="method">
1531 <dt id="sphere.sim.readTime">
1532 <code class="sig-name descname">readTime</code><span class="sig-paren">(</span><em class="sig-param">time</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.readTime" title="Permalink to this definition">¶</a></dt>
1533 <dd><p>Read the output file most closely corresponding to the time given as an
1534 argument.</p>
1535 <dl class="field-list simple">
1536 <dt class="field-odd">Parameters</dt>
1537 <dd class="field-odd"><p><strong>time</strong> (<em>float</em>) – The desired current time [s]</p>
1538 </dd>
1539 </dl>
1540 <p>See also <a class="reference internal" href="#sphere.sim.readbin" title="sphere.sim.readbin"><code class="xref py py-func docutils literal notranslate"><span class="pre">readbin()</span></code></a>, <a class="reference internal" href="#sphere.sim.readfirst" title="sphere.sim.readfirst"><code class="xref py py-func docutils literal notranslate"><span class="pre">readfirst()</span></code></a>, <a class="reference internal" href="#sphere.sim.readsecond" title="sphere.sim.readsecond"><code class="xref py py-func docutils literal notranslate"><span class="pre">readsecond()</span></code></a>, and
1541 <a class="reference internal" href="#sphere.sim.readstep" title="sphere.sim.readstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">readstep()</span></code></a>.</p>
1542 </dd></dl>
1543
1544 <dl class="method">
1545 <dt id="sphere.sim.readbin">
1546 <code class="sig-name descname">readbin</code><span class="sig-paren">(</span><em class="sig-param">targetbin</em>, <em class="sig-param">verbose=True</em>, <em class="sig-param">bonds=True</em>, <em class="sig-param">sigma0mod=True</em>, <em class="sig-param">esysparticle=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.readbin" title="Permalink to this definition">¶</a></dt>
1547 <dd><p>Reads a target <code class="docutils literal notranslate"><span class="pre">sphere</span></code> binary file.</p>
1548 <p>See also <a class="reference internal" href="#sphere.sim.writebin" title="sphere.sim.writebin"><code class="xref py py-func docutils literal notranslate"><span class="pre">writebin()</span></code></a>, <a class="reference internal" href="#sphere.sim.readfirst" title="sphere.sim.readfirst"><code class="xref py py-func docutils literal notranslate"><span class="pre">readfirst()</span></code></a>, <a class="reference internal" href="#sphere.sim.readlast" title="sphere.sim.readlast"><code class="xref py py-func docutils literal notranslate"><span class="pre">readlast()</span></code></a>,
1549 <a class="reference internal" href="#sphere.sim.readsecond" title="sphere.sim.readsecond"><code class="xref py py-func docutils literal notranslate"><span class="pre">readsecond()</span></code></a>, and <a class="reference internal" href="#sphere.sim.readstep" title="sphere.sim.readstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">readstep()</span></code></a>.</p>
1550 <dl class="field-list simple">
1551 <dt class="field-odd">Parameters</dt>
1552 <dd class="field-odd"><ul class="simple">
1553 <li><p><strong>targetbin</strong> (<em>str</em>) – The path to the binary <code class="docutils literal notranslate"><span class="pre">sphere</span></code> file</p></li>
1554 <li><p><strong>verbose</strong> (<em>bool</em>) – Show diagnostic information (default=True)</p></li>
1555 <li><p><strong>bonds</strong> (<em>bool</em>) – The input file contains bond information (default=True).
1556 This parameter should be true for all recent <code class="docutils literal notranslate"><span class="pre">sphere</span></code> versions.</p></li>
1557 <li><p><strong>sigma0mod</strong> (<em>bool</em>) – The input file contains information about modulating
1558 stresses at the top wall (default=True). This parameter should be
1559 true for all recent <code class="docutils literal notranslate"><span class="pre">sphere</span></code> versions.</p></li>
1560 <li><p><strong>esysparticle</strong> (<em>bool</em>) – Stop reading the file after reading the kinematics,
1561 which is useful for reading output files from other DEM programs.
1562 (default=False)</p></li>
1563 </ul>
1564 </dd>
1565 </dl>
1566 </dd></dl>
1567
1568 <dl class="method">
1569 <dt id="sphere.sim.readfirst">
1570 <code class="sig-name descname">readfirst</code><span class="sig-paren">(</span><em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.readfirst" title="Permalink to this definition">¶</a></dt>
1571 <dd><p>Read the first output file from the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder, corresponding
1572 to the object simulation id (<code class="docutils literal notranslate"><span class="pre">self.sid</span></code>).</p>
1573 <dl class="field-list simple">
1574 <dt class="field-odd">Parameters</dt>
1575 <dd class="field-odd"><p><strong>verbose</strong> (<em>bool</em>) – Display diagnostic information (default=True)</p>
1576 </dd>
1577 </dl>
1578 <p>See also <a class="reference internal" href="#sphere.sim.readbin" title="sphere.sim.readbin"><code class="xref py py-func docutils literal notranslate"><span class="pre">readbin()</span></code></a>, <a class="reference internal" href="#sphere.sim.readlast" title="sphere.sim.readlast"><code class="xref py py-func docutils literal notranslate"><span class="pre">readlast()</span></code></a>, <a class="reference internal" href="#sphere.sim.readsecond" title="sphere.sim.readsecond"><code class="xref py py-func docutils literal notranslate"><span class="pre">readsecond()</span></code></a>, and
1579 <a class="reference internal" href="#sphere.sim.readstep" title="sphere.sim.readstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">readstep()</span></code></a>.</p>
1580 </dd></dl>
1581
1582 <dl class="method">
1583 <dt id="sphere.sim.readlast">
1584 <code class="sig-name descname">readlast</code><span class="sig-paren">(</span><em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.readlast" title="Permalink to this definition">¶</a></dt>
1585 <dd><p>Read the last output file from the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder, corresponding
1586 to the object simulation id (<code class="docutils literal notranslate"><span class="pre">self.sid</span></code>).</p>
1587 <dl class="field-list simple">
1588 <dt class="field-odd">Parameters</dt>
1589 <dd class="field-odd"><p><strong>verbose</strong> (<em>bool</em>) – Display diagnostic information (default=True)</p>
1590 </dd>
1591 </dl>
1592 <p>See also <a class="reference internal" href="#sphere.sim.readbin" title="sphere.sim.readbin"><code class="xref py py-func docutils literal notranslate"><span class="pre">readbin()</span></code></a>, <a class="reference internal" href="#sphere.sim.readfirst" title="sphere.sim.readfirst"><code class="xref py py-func docutils literal notranslate"><span class="pre">readfirst()</span></code></a>, <a class="reference internal" href="#sphere.sim.readsecond" title="sphere.sim.readsecond"><code class="xref py py-func docutils literal notranslate"><span class="pre">readsecond()</span></code></a>, and
1593 <a class="reference internal" href="#sphere.sim.readstep" title="sphere.sim.readstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">readstep()</span></code></a>.</p>
1594 </dd></dl>
1595
1596 <dl class="method">
1597 <dt id="sphere.sim.readsecond">
1598 <code class="sig-name descname">readsecond</code><span class="sig-paren">(</span><em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.readsecond" title="Permalink to this definition">¶</a></dt>
1599 <dd><p>Read the second output file from the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder,
1600 corresponding to the object simulation id (<code class="docutils literal notranslate"><span class="pre">self.sid</span></code>).</p>
1601 <dl class="field-list simple">
1602 <dt class="field-odd">Parameters</dt>
1603 <dd class="field-odd"><p><strong>verbose</strong> (<em>bool</em>) – Display diagnostic information (default=True)</p>
1604 </dd>
1605 </dl>
1606 <p>See also <a class="reference internal" href="#sphere.sim.readbin" title="sphere.sim.readbin"><code class="xref py py-func docutils literal notranslate"><span class="pre">readbin()</span></code></a>, <a class="reference internal" href="#sphere.sim.readfirst" title="sphere.sim.readfirst"><code class="xref py py-func docutils literal notranslate"><span class="pre">readfirst()</span></code></a>, <a class="reference internal" href="#sphere.sim.readlast" title="sphere.sim.readlast"><code class="xref py py-func docutils literal notranslate"><span class="pre">readlast()</span></code></a>,
1607 and <a class="reference internal" href="#sphere.sim.readstep" title="sphere.sim.readstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">readstep()</span></code></a>.</p>
1608 </dd></dl>
1609
1610 <dl class="method">
1611 <dt id="sphere.sim.readstep">
1612 <code class="sig-name descname">readstep</code><span class="sig-paren">(</span><em class="sig-param">step</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.readstep" title="Permalink to this definition">¶</a></dt>
1613 <dd><p>Read a output file from the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder, corresponding
1614 to the object simulation id (<code class="docutils literal notranslate"><span class="pre">self.sid</span></code>).</p>
1615 <dl class="field-list simple">
1616 <dt class="field-odd">Parameters</dt>
1617 <dd class="field-odd"><ul class="simple">
1618 <li><p><strong>step</strong> (<em>int</em>) – The output file number to read, starting from 0.</p></li>
1619 <li><p><strong>verbose</strong> (<em>bool</em>) – Display diagnostic information (default=True)</p></li>
1620 </ul>
1621 </dd>
1622 </dl>
1623 <p>See also <a class="reference internal" href="#sphere.sim.readbin" title="sphere.sim.readbin"><code class="xref py py-func docutils literal notranslate"><span class="pre">readbin()</span></code></a>, <a class="reference internal" href="#sphere.sim.readfirst" title="sphere.sim.readfirst"><code class="xref py py-func docutils literal notranslate"><span class="pre">readfirst()</span></code></a>, <a class="reference internal" href="#sphere.sim.readlast" title="sphere.sim.readlast"><code class="xref py py-func docutils literal notranslate"><span class="pre">readlast()</span></code></a>,
1624 and <a class="reference internal" href="#sphere.sim.readsecond" title="sphere.sim.readsecond"><code class="xref py py-func docutils literal notranslate"><span class="pre">readsecond()</span></code></a>.</p>
1625 </dd></dl>
1626
1627 <dl class="method">
1628 <dt id="sphere.sim.render">
1629 <code class="sig-name descname">render</code><span class="sig-paren">(</span><em class="sig-param">method='pres'</em>, <em class="sig-param">max_val=1000.0</em>, <em class="sig-param">lower_cutoff=0.0</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.render" title="Permalink to this definition">¶</a></dt>
1630 <dd><p>Using the built-in ray tracer, render all output files that belong to
1631 the simulation, determined by the simulation id (<code class="docutils literal notranslate"><span class="pre">sid</span></code>).</p>
1632 <dl class="field-list simple">
1633 <dt class="field-odd">Parameters</dt>
1634 <dd class="field-odd"><ul class="simple">
1635 <li><p><strong>method</strong> (<em>str</em>) – The color visualization method to use for the particles.
1636 Possible values are: ‘normal’: color all particles with the same
1637 color, ‘pres’: color by pressure, ‘vel’: color by translational
1638 velocity, ‘angvel’: color by rotational velocity, ‘xdisp’: color by
1639 total displacement along the x-axis, ‘angpos’: color by angular
1640 position.</p></li>
1641 <li><p><strong>max_val</strong> (<em>float</em>) – The maximum value of the color bar</p></li>
1642 <li><p><strong>lower_cutoff</strong> (<em>float</em>) – Do not render particles with a value below this
1643 value, of the field selected by <code class="docutils literal notranslate"><span class="pre">method</span></code></p></li>
1644 <li><p><strong>graphics_format</strong> (<em>str</em>) – Convert the PPM images generated by the ray
1645 tracer to this image format using Imagemagick</p></li>
1646 <li><p><strong>verbose</strong> (<em>bool</em>) – Show verbose information during ray tracing</p></li>
1647 </ul>
1648 </dd>
1649 </dl>
1650 </dd></dl>
1651
1652 <dl class="method">
1653 <dt id="sphere.sim.rotationalEnergy">
1654 <code class="sig-name descname">rotationalEnergy</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.rotationalEnergy" title="Permalink to this definition">¶</a></dt>
1655 <dd><p>Returns the rotational energy for a particle.</p>
1656 <dl class="field-list simple">
1657 <dt class="field-odd">Parameters</dt>
1658 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
1659 </dd>
1660 <dt class="field-even">Returns</dt>
1661 <dd class="field-even"><p>The rotational kinetic energy of the particle [J]</p>
1662 </dd>
1663 <dt class="field-odd">Return type</dt>
1664 <dd class="field-odd"><p>float</p>
1665 </dd>
1666 </dl>
1667 </dd></dl>
1668
1669 <dl class="method">
1670 <dt id="sphere.sim.run">
1671 <code class="sig-name descname">run</code><span class="sig-paren">(</span><em class="sig-param">verbose=True</em>, <em class="sig-param">hideinputfile=False</em>, <em class="sig-param">dry=False</em>, <em class="sig-param">valgrind=False</em>, <em class="sig-param">cudamemcheck=False</em>, <em class="sig-param">device=-1</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.run" title="Permalink to this definition">¶</a></dt>
1672 <dd><p>Start <code class="docutils literal notranslate"><span class="pre">sphere</span></code> calculations on the <code class="docutils literal notranslate"><span class="pre">sim</span></code> object</p>
1673 <dl class="field-list simple">
1674 <dt class="field-odd">Parameters</dt>
1675 <dd class="field-odd"><ul class="simple">
1676 <li><p><strong>verbose</strong> (<em>bool</em>) – Show <code class="docutils literal notranslate"><span class="pre">sphere</span></code> output</p></li>
1677 <li><p><strong>hideinputfile</strong> (<em>bool</em>) – Hide the file name of the <code class="docutils literal notranslate"><span class="pre">sphere</span></code> input file</p></li>
1678 <li><p><strong>dry</strong> (<em>bool</em>) – Perform a dry run. Important parameter values are shown by
1679 the <code class="docutils literal notranslate"><span class="pre">sphere</span></code> program, and it exits afterwards.</p></li>
1680 <li><p><strong>valgrind</strong> (<em>bool</em>) – Run the program with <code class="docutils literal notranslate"><span class="pre">valgrind</span></code> in order to check
1681 memory leaks in the host code. This causes a significant increase in
1682 computational time.</p></li>
1683 <li><p><strong>cudamemcheck</strong> (<em>bool</em>) – Run the program with <code class="docutils literal notranslate"><span class="pre">cudamemcheck</span></code> in order to
1684 check for device memory leaks and errors. This causes a significant
1685 increase in computational time.</p></li>
1686 <li><p><strong>device</strong> (<em>int</em>) – Specify the GPU device to execute the program on.
1687 If not specified, sphere will use the device with the most CUDA cores.
1688 To see a list of devices, run <code class="docutils literal notranslate"><span class="pre">nvidia-smi</span></code> in the system shell.</p></li>
1689 </ul>
1690 </dd>
1691 </dl>
1692 </dd></dl>
1693
1694 <dl class="method">
1695 <dt id="sphere.sim.scaleSize">
1696 <code class="sig-name descname">scaleSize</code><span class="sig-paren">(</span><em class="sig-param">factor</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.scaleSize" title="Permalink to this definition">¶</a></dt>
1697 <dd><p>Scale the positions, linear velocities, forces, torques and radii of all
1698 particles and mobile walls.</p>
1699 <dl class="field-list simple">
1700 <dt class="field-odd">Parameters</dt>
1701 <dd class="field-odd"><p><strong>factor</strong> (<em>float</em>) – Spatial scaling factor ]0;inf[</p>
1702 </dd>
1703 </dl>
1704 </dd></dl>
1705
1706 <dl class="method">
1707 <dt id="sphere.sim.setBeta">
1708 <code class="sig-name descname">setBeta</code><span class="sig-paren">(</span><em class="sig-param">beta</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setBeta" title="Permalink to this definition">¶</a></dt>
1709 <dd><p>Beta is a fluid solver parameter, used in velocity prediction and
1710 pressure iteration 1.0: Use old pressures for fluid velocity prediction
1711 (see Langtangen et al. 2002) 0.0: Do not use old pressures for fluid
1712 velocity prediction (Chorin’s original projection method, see Chorin
1713 (1968) and “Projection method (fluid dynamics)” page on Wikipedia. The
1714 best results precision and performance-wise are obtained by using a beta
1715 of 0 and a low tolerance criteria value.</p>
1716 <p>The default and recommended value is 0.0.</p>
1717 <p>Other solver parameter setting functions: <a class="reference internal" href="#sphere.sim.setGamma" title="sphere.sim.setGamma"><code class="xref py py-func docutils literal notranslate"><span class="pre">setGamma()</span></code></a>,
1718 <a class="reference internal" href="#sphere.sim.setTheta" title="sphere.sim.setTheta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTheta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setTolerance" title="sphere.sim.setTolerance"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTolerance()</span></code></a>,
1719 <a class="reference internal" href="#sphere.sim.setDEMstepsPerCFDstep" title="sphere.sim.setDEMstepsPerCFDstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">setDEMstepsPerCFDstep()</span></code></a> and
1720 <a class="reference internal" href="#sphere.sim.setMaxIterations" title="sphere.sim.setMaxIterations"><code class="xref py py-func docutils literal notranslate"><span class="pre">setMaxIterations()</span></code></a></p>
1721 </dd></dl>
1722
1723 <dl class="method">
1724 <dt id="sphere.sim.setDEMstepsPerCFDstep">
1725 <code class="sig-name descname">setDEMstepsPerCFDstep</code><span class="sig-paren">(</span><em class="sig-param">ndem</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setDEMstepsPerCFDstep" title="Permalink to this definition">¶</a></dt>
1726 <dd><p>A fluid solver parameter, the value of the maxiter parameter denotes the
1727 number of DEM time steps to be performed per CFD time step.</p>
1728 <p>The default value is 1.</p>
1729 <dl class="field-list simple">
1730 <dt class="field-odd">Parameters</dt>
1731 <dd class="field-odd"><p><strong>ndem</strong> (<em>int</em>) – The DEM/CFD time step ratio</p>
1732 </dd>
1733 </dl>
1734 <p>Other solver parameter setting functions: <a class="reference internal" href="#sphere.sim.setGamma" title="sphere.sim.setGamma"><code class="xref py py-func docutils literal notranslate"><span class="pre">setGamma()</span></code></a>,
1735 <a class="reference internal" href="#sphere.sim.setTheta" title="sphere.sim.setTheta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTheta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setBeta" title="sphere.sim.setBeta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setBeta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setTolerance" title="sphere.sim.setTolerance"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTolerance()</span></code></a> and
1736 <a class="reference internal" href="#sphere.sim.setMaxIterations" title="sphere.sim.setMaxIterations"><code class="xref py py-func docutils literal notranslate"><span class="pre">setMaxIterations()</span></code></a>.</p>
1737 </dd></dl>
1738
1739 <dl class="method">
1740 <dt id="sphere.sim.setDampingNormal">
1741 <code class="sig-name descname">setDampingNormal</code><span class="sig-paren">(</span><em class="sig-param">gamma</em>, <em class="sig-param">over_damping=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setDampingNormal" title="Permalink to this definition">¶</a></dt>
1742 <dd><p>Set the dampening coefficient (gamma) in the normal direction of the
1743 particle-particle contact model. The function will print the fraction
1744 between the chosen damping and the critical damping value.</p>
1745 <dl class="field-list simple">
1746 <dt class="field-odd">Parameters</dt>
1747 <dd class="field-odd"><ul class="simple">
1748 <li><p><strong>gamma</strong> (<em>float</em>) – The viscous damping constant [N/(m/s)]</p></li>
1749 <li><p><strong>over_damping</strong> (<em>boolean</em>) – Accept overdampening</p></li>
1750 </ul>
1751 </dd>
1752 </dl>
1753 <p>See also: <code class="xref py py-func docutils literal notranslate"><span class="pre">setDampingTangential(gamma)()</span></code></p>
1754 </dd></dl>
1755
1756 <dl class="method">
1757 <dt id="sphere.sim.setDampingTangential">
1758 <code class="sig-name descname">setDampingTangential</code><span class="sig-paren">(</span><em class="sig-param">gamma</em>, <em class="sig-param">over_damping=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setDampingTangential" title="Permalink to this definition">¶</a></dt>
1759 <dd><p>Set the dampening coefficient (gamma) in the tangential direction of the
1760 particle-particle contact model. The function will print the fraction
1761 between the chosen damping and the critical damping value.</p>
1762 <dl class="field-list simple">
1763 <dt class="field-odd">Parameters</dt>
1764 <dd class="field-odd"><ul class="simple">
1765 <li><p><strong>gamma</strong> (<em>float</em>) – The viscous damping constant [N/(m/s)]</p></li>
1766 <li><p><strong>over_damping</strong> (<em>boolean</em>) – Accept overdampening</p></li>
1767 </ul>
1768 </dd>
1769 </dl>
1770 <p>See also: <code class="xref py py-func docutils literal notranslate"><span class="pre">setDampingNormal(gamma)()</span></code></p>
1771 </dd></dl>
1772
1773 <dl class="method">
1774 <dt id="sphere.sim.setDynamicFriction">
1775 <code class="sig-name descname">setDynamicFriction</code><span class="sig-paren">(</span><em class="sig-param">mu_d</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setDynamicFriction" title="Permalink to this definition">¶</a></dt>
1776 <dd><p>Set the dynamic friction coefficient for particle-particle interactions
1777 (<cite>self.mu_d</cite>). This value describes the resistance to a shearing motion
1778 while it is happening (contact tangential velocity larger than 0).
1779 Strain softening can be introduced by having a smaller dynamic
1780 frictional coefficient than the static fricion coefficient. Usually this
1781 value is identical to the static friction coefficient.</p>
1782 <dl class="field-list simple">
1783 <dt class="field-odd">Parameters</dt>
1784 <dd class="field-odd"><p><strong>mu_d</strong> (<em>float</em>) – Value of the dynamic friction coefficient, in [0;inf[.
1785 Usually between 0 and 1.</p>
1786 </dd>
1787 </dl>
1788 <p>See also: <code class="xref py py-func docutils literal notranslate"><span class="pre">setStaticFriction(mu_s)()</span></code></p>
1789 </dd></dl>
1790
1791 <dl class="method">
1792 <dt id="sphere.sim.setFluidBottomFixedFlux">
1793 <code class="sig-name descname">setFluidBottomFixedFlux</code><span class="sig-paren">(</span><em class="sig-param">specific_flux</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidBottomFixedFlux" title="Permalink to this definition">¶</a></dt>
1794 <dd><p>Define a constant fluid flux normal to the boundary.</p>
1795 <p>The default behavior for the boundary is fixed value (Dirichlet), see
1796 <a class="reference internal" href="#sphere.sim.setFluidBottomFixedPressure" title="sphere.sim.setFluidBottomFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidBottomFixedPressure()</span></code></a>.</p>
1797 <dl class="field-list simple">
1798 <dt class="field-odd">Parameters</dt>
1799 <dd class="field-odd"><p><strong>specific_flux</strong> – Specific flux values across boundary (positive
1800 values upwards), [m/s]</p>
1801 </dd>
1802 </dl>
1803 </dd></dl>
1804
1805 <dl class="method">
1806 <dt id="sphere.sim.setFluidBottomFixedPressure">
1807 <code class="sig-name descname">setFluidBottomFixedPressure</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidBottomFixedPressure" title="Permalink to this definition">¶</a></dt>
1808 <dd><p>Set the lower boundary of the fluid domain to follow the fixed pressure
1809 value (Dirichlet) boundary condition.</p>
1810 <p>This is the default behavior for the boundary. See also
1811 <a class="reference internal" href="#sphere.sim.setFluidBottomNoFlow" title="sphere.sim.setFluidBottomNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidBottomNoFlow()</span></code></a></p>
1812 </dd></dl>
1813
1814 <dl class="method">
1815 <dt id="sphere.sim.setFluidBottomNoFlow">
1816 <code class="sig-name descname">setFluidBottomNoFlow</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidBottomNoFlow" title="Permalink to this definition">¶</a></dt>
1817 <dd><p>Set the lower boundary of the fluid domain to follow the no-flow
1818 (Neumann) boundary condition with free slip parallel to the boundary.</p>
1819 <p>The default behavior for the boundary is fixed value (Dirichlet), see
1820 <a class="reference internal" href="#sphere.sim.setFluidBottomFixedPressure" title="sphere.sim.setFluidBottomFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidBottomFixedPressure()</span></code></a>.</p>
1821 </dd></dl>
1822
1823 <dl class="method">
1824 <dt id="sphere.sim.setFluidBottomNoFlowNoSlip">
1825 <code class="sig-name descname">setFluidBottomNoFlowNoSlip</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidBottomNoFlowNoSlip" title="Permalink to this definition">¶</a></dt>
1826 <dd><p>Set the lower boundary of the fluid domain to follow the no-flow
1827 (Neumann) boundary condition with no slip parallel to the boundary.</p>
1828 <p>The default behavior for the boundary is fixed value (Dirichlet), see
1829 <a class="reference internal" href="#sphere.sim.setFluidBottomFixedPressure" title="sphere.sim.setFluidBottomFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidBottomFixedPressure()</span></code></a>.</p>
1830 </dd></dl>
1831
1832 <dl class="method">
1833 <dt id="sphere.sim.setFluidCompressibility">
1834 <code class="sig-name descname">setFluidCompressibility</code><span class="sig-paren">(</span><em class="sig-param">beta_f</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidCompressibility" title="Permalink to this definition">¶</a></dt>
1835 <dd><p>Set the fluid adiabatic compressibility [1/Pa]. This value is equal to
1836 <cite>1/K</cite> where <cite>K</cite> is the bulk modulus [Pa]. The value for water is 5.1e-10
1837 for water at 0 degrees Celcius. This parameter is used for the Darcy
1838 solver exclusively.</p>
1839 <dl class="field-list simple">
1840 <dt class="field-odd">Parameters</dt>
1841 <dd class="field-odd"><p><strong>beta_f</strong> (<em>float</em>) – The fluid compressibility [1/Pa]</p>
1842 </dd>
1843 </dl>
1844 <p>See also: <a class="reference internal" href="#sphere.sim.setFluidDensity" title="sphere.sim.setFluidDensity"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidDensity()</span></code></a> and <a class="reference internal" href="#sphere.sim.setFluidViscosity" title="sphere.sim.setFluidViscosity"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidViscosity()</span></code></a></p>
1845 </dd></dl>
1846
1847 <dl class="method">
1848 <dt id="sphere.sim.setFluidDensity">
1849 <code class="sig-name descname">setFluidDensity</code><span class="sig-paren">(</span><em class="sig-param">rho_f</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidDensity" title="Permalink to this definition">¶</a></dt>
1850 <dd><p>Set the fluid density [kg/(m*m*m)]. The value for water is 1000. This
1851 parameter is used for the Navier-Stokes fluid solver exclusively.</p>
1852 <dl class="field-list simple">
1853 <dt class="field-odd">Parameters</dt>
1854 <dd class="field-odd"><p><strong>rho_f</strong> (<em>float</em>) – The fluid density [kg/(m*m*m)]</p>
1855 </dd>
1856 </dl>
1857 <dl class="simple">
1858 <dt>See also: <a class="reference internal" href="#sphere.sim.setFluidViscosity" title="sphere.sim.setFluidViscosity"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidViscosity()</span></code></a> and</dt><dd><p><a class="reference internal" href="#sphere.sim.setFluidCompressibility" title="sphere.sim.setFluidCompressibility"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidCompressibility()</span></code></a></p>
1859 </dd>
1860 </dl>
1861 </dd></dl>
1862
1863 <dl class="method">
1864 <dt id="sphere.sim.setFluidPressureModulation">
1865 <code class="sig-name descname">setFluidPressureModulation</code><span class="sig-paren">(</span><em class="sig-param">A</em>, <em class="sig-param">f</em>, <em class="sig-param">phi=0.0</em>, <em class="sig-param">plot=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidPressureModulation" title="Permalink to this definition">¶</a></dt>
1866 <dd><p>Set the parameters for the sine wave modulating the fluid pressures
1867 at the top boundary. Note that a cos-wave is obtained with phi=pi/2.</p>
1868 <dl class="field-list simple">
1869 <dt class="field-odd">Parameters</dt>
1870 <dd class="field-odd"><ul class="simple">
1871 <li><p><strong>A</strong> (<em>float</em>) – Fluctuation amplitude [Pa]</p></li>
1872 <li><p><strong>f</strong> (<em>float</em>) – Fluctuation frequency [Hz]</p></li>
1873 <li><p><strong>phi</strong> (<em>float</em>) – Fluctuation phase shift (default=0.0) [rad]</p></li>
1874 <li><p><strong>plot</strong> (<em>bool</em>) – Show a plot of the resulting modulation</p></li>
1875 </ul>
1876 </dd>
1877 </dl>
1878 <p>See also: <a class="reference internal" href="#sphere.sim.setTopWallNormalStressModulation" title="sphere.sim.setTopWallNormalStressModulation"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTopWallNormalStressModulation()</span></code></a> and
1879 <a class="reference internal" href="#sphere.sim.disableFluidPressureModulation" title="sphere.sim.disableFluidPressureModulation"><code class="xref py py-func docutils literal notranslate"><span class="pre">disableFluidPressureModulation()</span></code></a></p>
1880 </dd></dl>
1881
1882 <dl class="method">
1883 <dt id="sphere.sim.setFluidTopFixedFlux">
1884 <code class="sig-name descname">setFluidTopFixedFlux</code><span class="sig-paren">(</span><em class="sig-param">specific_flux</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidTopFixedFlux" title="Permalink to this definition">¶</a></dt>
1885 <dd><p>Define a constant fluid flux normal to the boundary.</p>
1886 <p>The default behavior for the boundary is fixed value (Dirichlet), see
1887 <a class="reference internal" href="#sphere.sim.setFluidBottomFixedPressure" title="sphere.sim.setFluidBottomFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidBottomFixedPressure()</span></code></a>.</p>
1888 <dl class="field-list simple">
1889 <dt class="field-odd">Parameters</dt>
1890 <dd class="field-odd"><p><strong>specific_flux</strong> – Specific flux values across boundary (positive
1891 values upwards), [m/s]</p>
1892 </dd>
1893 </dl>
1894 </dd></dl>
1895
1896 <dl class="method">
1897 <dt id="sphere.sim.setFluidTopFixedPressure">
1898 <code class="sig-name descname">setFluidTopFixedPressure</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidTopFixedPressure" title="Permalink to this definition">¶</a></dt>
1899 <dd><p>Set the upper boundary of the fluid domain to follow the fixed pressure
1900 value (Dirichlet) boundary condition.</p>
1901 <p>This is the default behavior for the boundary. See also
1902 <a class="reference internal" href="#sphere.sim.setFluidTopNoFlow" title="sphere.sim.setFluidTopNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidTopNoFlow()</span></code></a></p>
1903 </dd></dl>
1904
1905 <dl class="method">
1906 <dt id="sphere.sim.setFluidTopNoFlow">
1907 <code class="sig-name descname">setFluidTopNoFlow</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidTopNoFlow" title="Permalink to this definition">¶</a></dt>
1908 <dd><p>Set the upper boundary of the fluid domain to follow the no-flow
1909 (Neumann) boundary condition with free slip parallel to the boundary.</p>
1910 <p>The default behavior for the boundary is fixed value (Dirichlet), see
1911 <a class="reference internal" href="#sphere.sim.setFluidTopFixedPressure" title="sphere.sim.setFluidTopFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidTopFixedPressure()</span></code></a>.</p>
1912 </dd></dl>
1913
1914 <dl class="method">
1915 <dt id="sphere.sim.setFluidTopNoFlowNoSlip">
1916 <code class="sig-name descname">setFluidTopNoFlowNoSlip</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidTopNoFlowNoSlip" title="Permalink to this definition">¶</a></dt>
1917 <dd><p>Set the upper boundary of the fluid domain to follow the no-flow
1918 (Neumann) boundary condition with no slip parallel to the boundary.</p>
1919 <p>The default behavior for the boundary is fixed value (Dirichlet), see
1920 <a class="reference internal" href="#sphere.sim.setFluidTopFixedPressure" title="sphere.sim.setFluidTopFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidTopFixedPressure()</span></code></a>.</p>
1921 </dd></dl>
1922
1923 <dl class="method">
1924 <dt id="sphere.sim.setFluidViscosity">
1925 <code class="sig-name descname">setFluidViscosity</code><span class="sig-paren">(</span><em class="sig-param">mu</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidViscosity" title="Permalink to this definition">¶</a></dt>
1926 <dd><p>Set the fluid dynamic viscosity [Pa*s]. The value for water is
1927 1.797e-3 at 0 degrees Celcius. This parameter is used for both the Darcy
1928 and Navier-Stokes fluid solver.</p>
1929 <dl class="field-list simple">
1930 <dt class="field-odd">Parameters</dt>
1931 <dd class="field-odd"><p><strong>mu</strong> (<em>float</em>) – The fluid dynamic viscosity [Pa*s]</p>
1932 </dd>
1933 </dl>
1934 <dl class="simple">
1935 <dt>See also: <a class="reference internal" href="#sphere.sim.setFluidDensity" title="sphere.sim.setFluidDensity"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidDensity()</span></code></a> and</dt><dd><p><a class="reference internal" href="#sphere.sim.setFluidCompressibility" title="sphere.sim.setFluidCompressibility"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidCompressibility()</span></code></a></p>
1936 </dd>
1937 </dl>
1938 </dd></dl>
1939
1940 <dl class="method">
1941 <dt id="sphere.sim.setFluidXFixedPressure">
1942 <code class="sig-name descname">setFluidXFixedPressure</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidXFixedPressure" title="Permalink to this definition">¶</a></dt>
1943 <dd><p>Set the X boundaries of the fluid domain to follow the fixed pressure
1944 value (Dirichlet) boundary condition.</p>
1945 <p>This is not the default behavior for the boundary. See also
1946 <a class="reference internal" href="#sphere.sim.setFluidXFixedPressure" title="sphere.sim.setFluidXFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXFixedPressure()</span></code></a>,
1947 <a class="reference internal" href="#sphere.sim.setFluidXNoFlow" title="sphere.sim.setFluidXNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXNoFlow()</span></code></a>, and
1948 <a class="reference internal" href="#sphere.sim.setFluidXPeriodic" title="sphere.sim.setFluidXPeriodic"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXPeriodic()</span></code></a> (default)</p>
1949 </dd></dl>
1950
1951 <dl class="method">
1952 <dt id="sphere.sim.setFluidXNoFlow">
1953 <code class="sig-name descname">setFluidXNoFlow</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidXNoFlow" title="Permalink to this definition">¶</a></dt>
1954 <dd><p>Set the X boundaries of the fluid domain to follow the no-flow
1955 (Neumann) boundary condition.</p>
1956 <p>This is not the default behavior for the boundary. See also
1957 <a class="reference internal" href="#sphere.sim.setFluidXFixedPressure" title="sphere.sim.setFluidXFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXFixedPressure()</span></code></a>,
1958 <a class="reference internal" href="#sphere.sim.setFluidXNoFlow" title="sphere.sim.setFluidXNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXNoFlow()</span></code></a>, and
1959 <a class="reference internal" href="#sphere.sim.setFluidXPeriodic" title="sphere.sim.setFluidXPeriodic"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXPeriodic()</span></code></a> (default)</p>
1960 </dd></dl>
1961
1962 <dl class="method">
1963 <dt id="sphere.sim.setFluidXPeriodic">
1964 <code class="sig-name descname">setFluidXPeriodic</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidXPeriodic" title="Permalink to this definition">¶</a></dt>
1965 <dd><p>Set the X boundaries of the fluid domain to follow the periodic
1966 (cyclic) boundary condition.</p>
1967 <p>This is the default behavior for the boundary. See also
1968 <a class="reference internal" href="#sphere.sim.setFluidXFixedPressure" title="sphere.sim.setFluidXFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXFixedPressure()</span></code></a> and
1969 <a class="reference internal" href="#sphere.sim.setFluidXNoFlow" title="sphere.sim.setFluidXNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidXNoFlow()</span></code></a></p>
1970 </dd></dl>
1971
1972 <dl class="method">
1973 <dt id="sphere.sim.setFluidYFixedPressure">
1974 <code class="sig-name descname">setFluidYFixedPressure</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidYFixedPressure" title="Permalink to this definition">¶</a></dt>
1975 <dd><p>Set the Y boundaries of the fluid domain to follow the fixed pressure
1976 value (Dirichlet) boundary condition.</p>
1977 <p>This is not the default behavior for the boundary. See also
1978 <a class="reference internal" href="#sphere.sim.setFluidYNoFlow" title="sphere.sim.setFluidYNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidYNoFlow()</span></code></a> and
1979 <a class="reference internal" href="#sphere.sim.setFluidYPeriodic" title="sphere.sim.setFluidYPeriodic"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidYPeriodic()</span></code></a> (default)</p>
1980 </dd></dl>
1981
1982 <dl class="method">
1983 <dt id="sphere.sim.setFluidYNoFlow">
1984 <code class="sig-name descname">setFluidYNoFlow</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidYNoFlow" title="Permalink to this definition">¶</a></dt>
1985 <dd><p>Set the Y boundaries of the fluid domain to follow the no-flow
1986 (Neumann) boundary condition.</p>
1987 <p>This is not the default behavior for the boundary. See also
1988 <a class="reference internal" href="#sphere.sim.setFluidYFixedPressure" title="sphere.sim.setFluidYFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidYFixedPressure()</span></code></a> and
1989 <a class="reference internal" href="#sphere.sim.setFluidYPeriodic" title="sphere.sim.setFluidYPeriodic"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidYPeriodic()</span></code></a> (default)</p>
1990 </dd></dl>
1991
1992 <dl class="method">
1993 <dt id="sphere.sim.setFluidYPeriodic">
1994 <code class="sig-name descname">setFluidYPeriodic</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setFluidYPeriodic" title="Permalink to this definition">¶</a></dt>
1995 <dd><p>Set the Y boundaries of the fluid domain to follow the periodic
1996 (cyclic) boundary condition.</p>
1997 <p>This is the default behavior for the boundary. See also
1998 <a class="reference internal" href="#sphere.sim.setFluidYFixedPressure" title="sphere.sim.setFluidYFixedPressure"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidYFixedPressure()</span></code></a> and
1999 <a class="reference internal" href="#sphere.sim.setFluidYNoFlow" title="sphere.sim.setFluidYNoFlow"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidYNoFlow()</span></code></a></p>
2000 </dd></dl>
2001
2002 <dl class="method">
2003 <dt id="sphere.sim.setGamma">
2004 <code class="sig-name descname">setGamma</code><span class="sig-paren">(</span><em class="sig-param">gamma</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setGamma" title="Permalink to this definition">¶</a></dt>
2005 <dd><p>Gamma is a fluid solver parameter, used for smoothing the pressure
2006 values. The epsilon (pressure) values are smoothed by including the
2007 average epsilon value of the six closest (face) neighbor cells. This
2008 parameter should be in the range [0.0;1.0[. The higher the value, the
2009 more averaging is introduced. A value of 0.0 disables all averaging.</p>
2010 <p>The default and recommended value is 0.0.</p>
2011 <dl class="field-list simple">
2012 <dt class="field-odd">Parameters</dt>
2013 <dd class="field-odd"><p><strong>theta</strong> (<em>float</em>) – The smoothing parameter value</p>
2014 </dd>
2015 </dl>
2016 <p>Other solver parameter setting functions: <a class="reference internal" href="#sphere.sim.setTheta" title="sphere.sim.setTheta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTheta()</span></code></a>,
2017 <a class="reference internal" href="#sphere.sim.setBeta" title="sphere.sim.setBeta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setBeta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setTolerance" title="sphere.sim.setTolerance"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTolerance()</span></code></a>,
2018 <a class="reference internal" href="#sphere.sim.setDEMstepsPerCFDstep" title="sphere.sim.setDEMstepsPerCFDstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">setDEMstepsPerCFDstep()</span></code></a> and <a class="reference internal" href="#sphere.sim.setMaxIterations" title="sphere.sim.setMaxIterations"><code class="xref py py-func docutils literal notranslate"><span class="pre">setMaxIterations()</span></code></a></p>
2019 </dd></dl>
2020
2021 <dl class="method">
2022 <dt id="sphere.sim.setMaxIterations">
2023 <code class="sig-name descname">setMaxIterations</code><span class="sig-paren">(</span><em class="sig-param">maxiter</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setMaxIterations" title="Permalink to this definition">¶</a></dt>
2024 <dd><p>A fluid solver parameter, the value of the maxiter parameter denotes the
2025 maximal allowed number of fluid solver iterations before ending the
2026 fluid solver loop prematurely. The residual values are at that point not
2027 fulfilling the tolerance criteria. The parameter is included to avoid
2028 infinite hangs.</p>
2029 <p>The default and recommended value is 1e4.</p>
2030 <dl class="field-list simple">
2031 <dt class="field-odd">Parameters</dt>
2032 <dd class="field-odd"><p><strong>maxiter</strong> (<em>int</em>) – The maximum number of Jacobi iterations in the fluid
2033 solver</p>
2034 </dd>
2035 </dl>
2036 <p>Other solver parameter setting functions: <a class="reference internal" href="#sphere.sim.setGamma" title="sphere.sim.setGamma"><code class="xref py py-func docutils literal notranslate"><span class="pre">setGamma()</span></code></a>,
2037 <a class="reference internal" href="#sphere.sim.setTheta" title="sphere.sim.setTheta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTheta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setBeta" title="sphere.sim.setBeta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setBeta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setDEMstepsPerCFDstep" title="sphere.sim.setDEMstepsPerCFDstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">setDEMstepsPerCFDstep()</span></code></a>
2038 and <a class="reference internal" href="#sphere.sim.setTolerance" title="sphere.sim.setTolerance"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTolerance()</span></code></a></p>
2039 </dd></dl>
2040
2041 <dl class="method">
2042 <dt id="sphere.sim.setPermeabilityGrainSize">
2043 <code class="sig-name descname">setPermeabilityGrainSize</code><span class="sig-paren">(</span><em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setPermeabilityGrainSize" title="Permalink to this definition">¶</a></dt>
2044 <dd><p>Set the permeability prefactor based on the mean grain size (Damsgaard
2045 et al., 2015, eq. 10).</p>
2046 <dl class="field-list simple">
2047 <dt class="field-odd">Parameters</dt>
2048 <dd class="field-odd"><p><strong>verbose</strong> (<em>bool</em>) – Print information about the realistic permeabilities
2049 hydraulic conductivities to expect with the chosen permeability
2050 prefactor.</p>
2051 </dd>
2052 </dl>
2053 </dd></dl>
2054
2055 <dl class="method">
2056 <dt id="sphere.sim.setPermeabilityPrefactor">
2057 <code class="sig-name descname">setPermeabilityPrefactor</code><span class="sig-paren">(</span><em class="sig-param">k_c</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setPermeabilityPrefactor" title="Permalink to this definition">¶</a></dt>
2058 <dd><p>Set the permeability prefactor from Goren et al 2011, eq. 24. The
2059 function will print the limits of permeabilities to be simulated. This
2060 parameter is only used in the Darcy solver.</p>
2061 <dl class="field-list simple">
2062 <dt class="field-odd">Parameters</dt>
2063 <dd class="field-odd"><ul class="simple">
2064 <li><p><strong>k_c</strong> (<em>float</em>) – Permeability prefactor value [m*m]</p></li>
2065 <li><p><strong>verbose</strong> (<em>bool</em>) – Print information about the realistic permeabilities and
2066 hydraulic conductivities to expect with the chosen permeability
2067 prefactor.</p></li>
2068 </ul>
2069 </dd>
2070 </dl>
2071 </dd></dl>
2072
2073 <dl class="method">
2074 <dt id="sphere.sim.setStaticFriction">
2075 <code class="sig-name descname">setStaticFriction</code><span class="sig-paren">(</span><em class="sig-param">mu_s</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setStaticFriction" title="Permalink to this definition">¶</a></dt>
2076 <dd><p>Set the static friction coefficient for particle-particle interactions
2077 (<cite>self.mu_s</cite>). This value describes the resistance to a shearing motion
2078 while it is not happenind (contact tangential velocity zero).</p>
2079 <dl class="field-list simple">
2080 <dt class="field-odd">Parameters</dt>
2081 <dd class="field-odd"><p><strong>mu_s</strong> (<em>float</em>) – Value of the static friction coefficient, in [0;inf[.
2082 Usually between 0 and 1.</p>
2083 </dd>
2084 </dl>
2085 <p>See also: <code class="xref py py-func docutils literal notranslate"><span class="pre">setDynamicFriction(mu_d)()</span></code></p>
2086 </dd></dl>
2087
2088 <dl class="method">
2089 <dt id="sphere.sim.setStiffnessNormal">
2090 <code class="sig-name descname">setStiffnessNormal</code><span class="sig-paren">(</span><em class="sig-param">k_n</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setStiffnessNormal" title="Permalink to this definition">¶</a></dt>
2091 <dd><p>Set the elastic stiffness (<cite>k_n</cite>) in the normal direction of the
2092 contact.</p>
2093 <dl class="field-list simple">
2094 <dt class="field-odd">Parameters</dt>
2095 <dd class="field-odd"><p><strong>k_n</strong> (<em>float</em>) – The elastic stiffness coefficient [N/m]</p>
2096 </dd>
2097 </dl>
2098 </dd></dl>
2099
2100 <dl class="method">
2101 <dt id="sphere.sim.setStiffnessTangential">
2102 <code class="sig-name descname">setStiffnessTangential</code><span class="sig-paren">(</span><em class="sig-param">k_t</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setStiffnessTangential" title="Permalink to this definition">¶</a></dt>
2103 <dd><p>Set the elastic stiffness (<cite>k_t</cite>) in the tangential direction of the
2104 contact.</p>
2105 <dl class="field-list simple">
2106 <dt class="field-odd">Parameters</dt>
2107 <dd class="field-odd"><p><strong>k_t</strong> (<em>float</em>) – The elastic stiffness coefficient [N/m]</p>
2108 </dd>
2109 </dl>
2110 </dd></dl>
2111
2112 <dl class="method">
2113 <dt id="sphere.sim.setTheta">
2114 <code class="sig-name descname">setTheta</code><span class="sig-paren">(</span><em class="sig-param">theta</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setTheta" title="Permalink to this definition">¶</a></dt>
2115 <dd><p>Theta is a fluid solver under-relaxation parameter, used in solution of
2116 Poisson equation. The value should be within the range ]0.0;1.0]. At a
2117 value of 1.0, the new estimate of epsilon values is used exclusively. At
2118 lower values, a linear interpolation between new and old values is used.
2119 The solution typically converges faster with a value of 1.0, but
2120 instabilities may be avoided with lower values.</p>
2121 <p>The default and recommended value is 1.0.</p>
2122 <dl class="field-list simple">
2123 <dt class="field-odd">Parameters</dt>
2124 <dd class="field-odd"><p><strong>theta</strong> (<em>float</em>) – The under-relaxation parameter value</p>
2125 </dd>
2126 </dl>
2127 <p>Other solver parameter setting functions: <a class="reference internal" href="#sphere.sim.setGamma" title="sphere.sim.setGamma"><code class="xref py py-func docutils literal notranslate"><span class="pre">setGamma()</span></code></a>,
2128 <a class="reference internal" href="#sphere.sim.setBeta" title="sphere.sim.setBeta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setBeta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setTolerance" title="sphere.sim.setTolerance"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTolerance()</span></code></a>,
2129 <a class="reference internal" href="#sphere.sim.setDEMstepsPerCFDstep" title="sphere.sim.setDEMstepsPerCFDstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">setDEMstepsPerCFDstep()</span></code></a> and <a class="reference internal" href="#sphere.sim.setMaxIterations" title="sphere.sim.setMaxIterations"><code class="xref py py-func docutils literal notranslate"><span class="pre">setMaxIterations()</span></code></a></p>
2130 </dd></dl>
2131
2132 <dl class="method">
2133 <dt id="sphere.sim.setTolerance">
2134 <code class="sig-name descname">setTolerance</code><span class="sig-paren">(</span><em class="sig-param">tolerance</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setTolerance" title="Permalink to this definition">¶</a></dt>
2135 <dd><p>A fluid solver parameter, the value of the tolerance parameter denotes
2136 the required value of the maximum normalized residual for the fluid
2137 solver.</p>
2138 <p>The default and recommended value is 1.0e-3.</p>
2139 <dl class="field-list simple">
2140 <dt class="field-odd">Parameters</dt>
2141 <dd class="field-odd"><p><strong>tolerance</strong> (<em>float</em>) – The tolerance criteria for the maximal normalized
2142 residual</p>
2143 </dd>
2144 </dl>
2145 <p>Other solver parameter setting functions: <a class="reference internal" href="#sphere.sim.setGamma" title="sphere.sim.setGamma"><code class="xref py py-func docutils literal notranslate"><span class="pre">setGamma()</span></code></a>,
2146 <a class="reference internal" href="#sphere.sim.setTheta" title="sphere.sim.setTheta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setTheta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setBeta" title="sphere.sim.setBeta"><code class="xref py py-func docutils literal notranslate"><span class="pre">setBeta()</span></code></a>, <a class="reference internal" href="#sphere.sim.setDEMstepsPerCFDstep" title="sphere.sim.setDEMstepsPerCFDstep"><code class="xref py py-func docutils literal notranslate"><span class="pre">setDEMstepsPerCFDstep()</span></code></a> and
2147 <a class="reference internal" href="#sphere.sim.setMaxIterations" title="sphere.sim.setMaxIterations"><code class="xref py py-func docutils literal notranslate"><span class="pre">setMaxIterations()</span></code></a></p>
2148 </dd></dl>
2149
2150 <dl class="method">
2151 <dt id="sphere.sim.setTopWallNormalStressModulation">
2152 <code class="sig-name descname">setTopWallNormalStressModulation</code><span class="sig-paren">(</span><em class="sig-param">A</em>, <em class="sig-param">f</em>, <em class="sig-param">plot=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setTopWallNormalStressModulation" title="Permalink to this definition">¶</a></dt>
2153 <dd><p>Set the parameters for the sine wave modulating the normal stress
2154 at the top wall. Note that a cos-wave is obtained with phi=pi/2.</p>
2155 <dl class="field-list simple">
2156 <dt class="field-odd">Parameters</dt>
2157 <dd class="field-odd"><ul class="simple">
2158 <li><p><strong>A</strong> (<em>float</em>) – Fluctuation amplitude [Pa]</p></li>
2159 <li><p><strong>f</strong> (<em>float</em>) – Fluctuation frequency [Hz]</p></li>
2160 <li><p><strong>plot</strong> (<em>bool</em>) – Show a plot of the resulting modulation</p></li>
2161 </ul>
2162 </dd>
2163 </dl>
2164 <p>See also: <a class="reference internal" href="#sphere.sim.setFluidPressureModulation" title="sphere.sim.setFluidPressureModulation"><code class="xref py py-func docutils literal notranslate"><span class="pre">setFluidPressureModulation()</span></code></a> and
2165 <a class="reference internal" href="#sphere.sim.disableTopWallNormalStressModulation" title="sphere.sim.disableTopWallNormalStressModulation"><code class="xref py py-func docutils literal notranslate"><span class="pre">disableTopWallNormalStressModulation()</span></code></a></p>
2166 </dd></dl>
2167
2168 <dl class="method">
2169 <dt id="sphere.sim.setYoungsModulus">
2170 <code class="sig-name descname">setYoungsModulus</code><span class="sig-paren">(</span><em class="sig-param">E</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.setYoungsModulus" title="Permalink to this definition">¶</a></dt>
2171 <dd><p>Set the elastic Young’s modulus (<cite>E</cite>) for the contact model. This
2172 parameter is used over normal stiffness (<cite>k_n</cite>) and tangential
2173 stiffness (<cite>k_t</cite>) when its value is greater than zero. Using this
2174 parameter produces size-invariant behavior.</p>
2175 <p>Example values are ~70e9 Pa for quartz,
2176 <a class="reference external" href="http://www.engineeringtoolbox.com/young-modulus-d_417.html">http://www.engineeringtoolbox.com/young-modulus-d_417.html</a></p>
2177 <dl class="field-list simple">
2178 <dt class="field-odd">Parameters</dt>
2179 <dd class="field-odd"><p><strong>E</strong> (<em>float</em>) – The elastic modulus [Pa]</p>
2180 </dd>
2181 </dl>
2182 </dd></dl>
2183
2184 <dl class="method">
2185 <dt id="sphere.sim.shear">
2186 <code class="sig-name descname">shear</code><span class="sig-paren">(</span><em class="sig-param">shear_strain_rate=1.0</em>, <em class="sig-param">shear_stress=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shear" title="Permalink to this definition">¶</a></dt>
2187 <dd><p>Setup shear experiment either by a constant shear rate or a constant
2188 shear stress. The shear strain rate is the shear velocity divided by
2189 the initial height per second. The shear movement is along the positive
2190 x axis. The function zeroes the tangential wall viscosity (gamma_wt) and
2191 the wall friction coefficients (mu_ws, mu_wn).</p>
2192 <dl class="field-list simple">
2193 <dt class="field-odd">Parameters</dt>
2194 <dd class="field-odd"><ul class="simple">
2195 <li><p><strong>shear_strain_rate</strong> (<em>float</em>) – The shear strain rate [-] to use if
2196 shear_stress isn’t False.</p></li>
2197 <li><p><strong>shear_stress</strong> (<em>float</em><em> or </em><em>bool</em>) – The shear stress value to use [Pa].</p></li>
2198 </ul>
2199 </dd>
2200 </dl>
2201 </dd></dl>
2202
2203 <dl class="method">
2204 <dt id="sphere.sim.shearDisplacement">
2205 <code class="sig-name descname">shearDisplacement</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shearDisplacement" title="Permalink to this definition">¶</a></dt>
2206 <dd><p>Calculates and returns the current shear displacement. The displacement
2207 is found by determining the total x-axis displacement of the upper,
2208 fixed particles.</p>
2209 <dl class="field-list simple">
2210 <dt class="field-odd">Returns</dt>
2211 <dd class="field-odd"><p>The total shear displacement [m]</p>
2212 </dd>
2213 <dt class="field-even">Return type</dt>
2214 <dd class="field-even"><p>float</p>
2215 </dd>
2216 </dl>
2217 <p>See also: <a class="reference internal" href="#sphere.sim.shearStrain" title="sphere.sim.shearStrain"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearStrain()</span></code></a> and <a class="reference internal" href="#sphere.sim.shearVelocity" title="sphere.sim.shearVelocity"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearVelocity()</span></code></a></p>
2218 </dd></dl>
2219
2220 <dl class="method">
2221 <dt id="sphere.sim.shearStrain">
2222 <code class="sig-name descname">shearStrain</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shearStrain" title="Permalink to this definition">¶</a></dt>
2223 <dd><p>Calculates and returns the current shear strain (gamma) value of the
2224 experiment. The shear strain is found by determining the total x-axis
2225 displacement of the upper, fixed particles.</p>
2226 <dl class="field-list simple">
2227 <dt class="field-odd">Returns</dt>
2228 <dd class="field-odd"><p>The total shear strain [-]</p>
2229 </dd>
2230 <dt class="field-even">Return type</dt>
2231 <dd class="field-even"><p>float</p>
2232 </dd>
2233 </dl>
2234 <p>See also: <a class="reference internal" href="#sphere.sim.shearStrainRate" title="sphere.sim.shearStrainRate"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearStrainRate()</span></code></a> and <a class="reference internal" href="#sphere.sim.shearVel" title="sphere.sim.shearVel"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearVel()</span></code></a></p>
2235 </dd></dl>
2236
2237 <dl class="method">
2238 <dt id="sphere.sim.shearStrainRate">
2239 <code class="sig-name descname">shearStrainRate</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shearStrainRate" title="Permalink to this definition">¶</a></dt>
2240 <dd><p>Calculates the shear strain rate (dot(gamma)) value of the experiment.</p>
2241 <dl class="field-list simple">
2242 <dt class="field-odd">Returns</dt>
2243 <dd class="field-odd"><p>The value of dot(gamma)</p>
2244 </dd>
2245 <dt class="field-even">Return type</dt>
2246 <dd class="field-even"><p>float</p>
2247 </dd>
2248 </dl>
2249 <p>See also: <a class="reference internal" href="#sphere.sim.shearStrain" title="sphere.sim.shearStrain"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearStrain()</span></code></a> and <a class="reference internal" href="#sphere.sim.shearVel" title="sphere.sim.shearVel"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearVel()</span></code></a></p>
2250 </dd></dl>
2251
2252 <dl class="method">
2253 <dt id="sphere.sim.shearStress">
2254 <code class="sig-name descname">shearStress</code><span class="sig-paren">(</span><em class="sig-param">type='effective'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shearStress" title="Permalink to this definition">¶</a></dt>
2255 <dd><p>Calculates the sum of shear stress values measured on any moving
2256 particles with a finite and fixed velocity.</p>
2257 <dl class="field-list simple">
2258 <dt class="field-odd">Parameters</dt>
2259 <dd class="field-odd"><p><strong>type</strong> (<em>str</em>) – Find the ‘defined’ or ‘effective’ (default) shear stress</p>
2260 </dd>
2261 <dt class="field-even">Returns</dt>
2262 <dd class="field-even"><p>The shear stress in Pa</p>
2263 </dd>
2264 <dt class="field-odd">Return type</dt>
2265 <dd class="field-odd"><p>numpy.array</p>
2266 </dd>
2267 </dl>
2268 </dd></dl>
2269
2270 <dl class="method">
2271 <dt id="sphere.sim.shearVel">
2272 <code class="sig-name descname">shearVel</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shearVel" title="Permalink to this definition">¶</a></dt>
2273 <dd><p>Alias of <a class="reference internal" href="#sphere.sim.shearVelocity" title="sphere.sim.shearVelocity"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearVelocity()</span></code></a></p>
2274 </dd></dl>
2275
2276 <dl class="method">
2277 <dt id="sphere.sim.shearVelocity">
2278 <code class="sig-name descname">shearVelocity</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.shearVelocity" title="Permalink to this definition">¶</a></dt>
2279 <dd><p>Calculates and returns the current shear velocity. The displacement
2280 is found by determining the total x-axis velocity of the upper,
2281 fixed particles.</p>
2282 <dl class="field-list simple">
2283 <dt class="field-odd">Returns</dt>
2284 <dd class="field-odd"><p>The shear velocity [m/s]</p>
2285 </dd>
2286 <dt class="field-even">Return type</dt>
2287 <dd class="field-even"><p>float</p>
2288 </dd>
2289 </dl>
2290 <p>See also: <a class="reference internal" href="#sphere.sim.shearStrainRate" title="sphere.sim.shearStrainRate"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearStrainRate()</span></code></a> and <a class="reference internal" href="#sphere.sim.shearDisplacement" title="sphere.sim.shearDisplacement"><code class="xref py py-func docutils literal notranslate"><span class="pre">shearDisplacement()</span></code></a></p>
2291 </dd></dl>
2292
2293 <dl class="method">
2294 <dt id="sphere.sim.sheardisp">
2295 <code class="sig-name descname">sheardisp</code><span class="sig-paren">(</span><em class="sig-param">graphics_format='pdf'</em>, <em class="sig-param">zslices=32</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.sheardisp" title="Permalink to this definition">¶</a></dt>
2296 <dd><p>Plot the particle x-axis displacement against the original vertical
2297 particle position. The plot is saved in the current directory with the
2298 file name ‘<simulation id>-sheardisp.<graphics_format>’.</p>
2299 <dl class="field-list simple">
2300 <dt class="field-odd">Parameters</dt>
2301 <dd class="field-odd"><p><strong>graphics_format</strong> (<em>str</em>) – Save the plot in this format</p>
2302 </dd>
2303 </dl>
2304 </dd></dl>
2305
2306 <dl class="method">
2307 <dt id="sphere.sim.show">
2308 <code class="sig-name descname">show</code><span class="sig-paren">(</span><em class="sig-param">coloring=array([]</em>, <em class="sig-param">dtype=float64)</em>, <em class="sig-param">resolution=6</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.show" title="Permalink to this definition">¶</a></dt>
2309 <dd><p>Show a rendering of all particles in a window.</p>
2310 <dl class="field-list simple">
2311 <dt class="field-odd">Parameters</dt>
2312 <dd class="field-odd"><ul class="simple">
2313 <li><p><strong>coloring</strong> (<em>numpy.array</em>) – Color the particles from red to white to blue according
2314 to the values in this array.</p></li>
2315 <li><p><strong>resolution</strong> (<em>int</em>) – The resolution of the rendered spheres. Larger values
2316 increase the performance requirements.</p></li>
2317 </ul>
2318 </dd>
2319 </dl>
2320 </dd></dl>
2321
2322 <dl class="method">
2323 <dt id="sphere.sim.smallestMass">
2324 <code class="sig-name descname">smallestMass</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.smallestMass" title="Permalink to this definition">¶</a></dt>
2325 <dd><p>Returns the mass of the leightest particle.</p>
2326 <dl class="field-list simple">
2327 <dt class="field-odd">Parameters</dt>
2328 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
2329 </dd>
2330 <dt class="field-even">Returns</dt>
2331 <dd class="field-even"><p>The mass of the particle [kg]</p>
2332 </dd>
2333 <dt class="field-odd">Return type</dt>
2334 <dd class="field-odd"><p>float</p>
2335 </dd>
2336 </dl>
2337 </dd></dl>
2338
2339 <dl class="method">
2340 <dt id="sphere.sim.staticGrid">
2341 <code class="sig-name descname">staticGrid</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.staticGrid" title="Permalink to this definition">¶</a></dt>
2342 <dd><p>Set the height of the fluid grid to be constant as set in <cite>self.L[2]</cite>.</p>
2343 <p>See also <a class="reference internal" href="#sphere.sim.adaptiveGrid" title="sphere.sim.adaptiveGrid"><code class="xref py py-func docutils literal notranslate"><span class="pre">adaptiveGrid()</span></code></a></p>
2344 </dd></dl>
2345
2346 <dl class="method">
2347 <dt id="sphere.sim.status">
2348 <code class="sig-name descname">status</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.status" title="Permalink to this definition">¶</a></dt>
2349 <dd><p>Returns the current simulation status by using the simulation id
2350 (<code class="docutils literal notranslate"><span class="pre">sid</span></code>) as an identifier.</p>
2351 <dl class="field-list simple">
2352 <dt class="field-odd">Returns</dt>
2353 <dd class="field-odd"><p>The number of the last output file written</p>
2354 </dd>
2355 <dt class="field-even">Return type</dt>
2356 <dd class="field-even"><p>int</p>
2357 </dd>
2358 </dl>
2359 </dd></dl>
2360
2361 <dl class="method">
2362 <dt id="sphere.sim.surfaceArea">
2363 <code class="sig-name descname">surfaceArea</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.surfaceArea" title="Permalink to this definition">¶</a></dt>
2364 <dd><p>Returns the surface area of a particle.</p>
2365 <dl class="field-list simple">
2366 <dt class="field-odd">Parameters</dt>
2367 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
2368 </dd>
2369 <dt class="field-even">Returns</dt>
2370 <dd class="field-even"><p>The surface area of the particle [m^2]</p>
2371 </dd>
2372 <dt class="field-odd">Return type</dt>
2373 <dd class="field-odd"><p>float</p>
2374 </dd>
2375 </dl>
2376 </dd></dl>
2377
2378 <dl class="method">
2379 <dt id="sphere.sim.thinsection_x1x3">
2380 <code class="sig-name descname">thinsection_x1x3</code><span class="sig-paren">(</span><em class="sig-param">x2='center'</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">cbmax=None</em>, <em class="sig-param">arrowscale=0.01</em>, <em class="sig-param">velarrowscale=1.0</em>, <em class="sig-param">slipscale=1.0</em>, <em class="sig-param">verbose=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.thinsection_x1x3" title="Permalink to this definition">¶</a></dt>
2381 <dd><p>Produce a 2D image of particles on a x1,x3 plane, intersecting the
2382 second axis at x2. Output is saved as ‘<sid>-ts-x1x3.txt’ in the
2383 current folder.</p>
2384 <p>An upper limit to the pressure color bar range can be set by the
2385 cbmax parameter.</p>
2386 <dl class="simple">
2387 <dt>The data can be plotted in gnuplot with:</dt><dd><p>gnuplot> set size ratio -1
2388 gnuplot> set palette defined (0 “blue”, 0.5 “gray”, 1 “red”)
2389 gnuplot> plot ‘<sid>-ts-x1x3.txt’ with circles palette fs transparent solid 0.4 noborder</p>
2390 </dd>
2391 </dl>
2392 <p>This function also saves a plot of the inter-particle slip angles.</p>
2393 <dl class="field-list simple">
2394 <dt class="field-odd">Parameters</dt>
2395 <dd class="field-odd"><ul class="simple">
2396 <li><p><strong>x2</strong> (<em>foat</em>) – The position along the second axis of the intersecting plane</p></li>
2397 <li><p><strong>graphics_format</strong> (<em>str</em>) – Save the slip angle plot in this format</p></li>
2398 <li><p><strong>cbmax</strong> (<em>float</em>) – The maximal value of the pressure color bar range</p></li>
2399 <li><p><strong>arrowscale</strong> (<em>float</em>) – Scale the rotational arrows by this value</p></li>
2400 <li><p><strong>velarrowscale</strong> (<em>float</em>) – Scale the translational arrows by this value</p></li>
2401 <li><p><strong>slipscale</strong> (<em>float</em>) – Scale the slip arrows by this value</p></li>
2402 <li><p><strong>verbose</strong> (<em>bool</em>) – Show function output during calculations</p></li>
2403 </ul>
2404 </dd>
2405 </dl>
2406 </dd></dl>
2407
2408 <dl class="method">
2409 <dt id="sphere.sim.torqueScript">
2410 <code class="sig-name descname">torqueScript</code><span class="sig-paren">(</span><em class="sig-param">email='adc@geo.au.dk'</em>, <em class="sig-param">email_alerts='ae'</em>, <em class="sig-param">walltime='24:00:00'</em>, <em class="sig-param">queue='qfermi'</em>, <em class="sig-param">cudapath='/com/cuda/4.0.17/cuda'</em>, <em class="sig-param">spheredir='/home/adc/code/sphere'</em>, <em class="sig-param">use_workdir=False</em>, <em class="sig-param">workdir='/scratch'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.torqueScript" title="Permalink to this definition">¶</a></dt>
2411 <dd><p>Creates a job script for the Torque queue manager for the simulation
2412 object.</p>
2413 <dl class="field-list simple">
2414 <dt class="field-odd">Parameters</dt>
2415 <dd class="field-odd"><ul class="simple">
2416 <li><p><strong>email</strong> (<em>str</em>) – The e-mail address that Torque messages should be sent to</p></li>
2417 <li><p><strong>email_alerts</strong> (<em>str</em>) – The type of Torque messages to send to the e-mail
2418 address. The character ‘b’ causes a mail to be sent when the
2419 execution begins. The character ‘e’ causes a mail to be sent when
2420 the execution ends normally. The character ‘a’ causes a mail to be
2421 sent if the execution ends abnormally. The characters can be written
2422 in any order.</p></li>
2423 <li><p><strong>walltime</strong> (<em>str</em>) – The maximal allowed time for the job, in the format
2424 ‘HH:MM:SS’.</p></li>
2425 <li><p><strong>queue</strong> (<em>str</em>) – The Torque queue to schedule the job for</p></li>
2426 <li><p><strong>cudapath</strong> (<em>str</em>) – The path of the CUDA library on the cluster compute
2427 nodes</p></li>
2428 <li><p><strong>spheredir</strong> (<em>str</em>) – The path to the root directory of sphere on the
2429 cluster</p></li>
2430 <li><p><strong>use_workdir</strong> (<em>bool</em>) – Use a different working directory than the sphere
2431 folder</p></li>
2432 <li><p><strong>workdir</strong> (<em>str</em>) – The working directory during the calculations, if
2433 <cite>use_workdir=True</cite></p></li>
2434 </ul>
2435 </dd>
2436 </dl>
2437 </dd></dl>
2438
2439 <dl class="method">
2440 <dt id="sphere.sim.totalFrictionalEnergy">
2441 <code class="sig-name descname">totalFrictionalEnergy</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.totalFrictionalEnergy" title="Permalink to this definition">¶</a></dt>
2442 <dd><p>Returns the total frictional dissipated energy for all particles.</p>
2443 <dl class="field-list simple">
2444 <dt class="field-odd">Returns</dt>
2445 <dd class="field-odd"><p>The total frictional energy lost of all particles [J]</p>
2446 </dd>
2447 <dt class="field-even">Return type</dt>
2448 <dd class="field-even"><p>float</p>
2449 </dd>
2450 </dl>
2451 </dd></dl>
2452
2453 <dl class="method">
2454 <dt id="sphere.sim.totalKineticEnergy">
2455 <code class="sig-name descname">totalKineticEnergy</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.totalKineticEnergy" title="Permalink to this definition">¶</a></dt>
2456 <dd><p>Returns the total linear kinetic energy for all particles.</p>
2457 <dl class="field-list simple">
2458 <dt class="field-odd">Returns</dt>
2459 <dd class="field-odd"><p>The kinetic energy of all particles [J]</p>
2460 </dd>
2461 </dl>
2462 </dd></dl>
2463
2464 <dl class="method">
2465 <dt id="sphere.sim.totalMass">
2466 <code class="sig-name descname">totalMass</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.totalMass" title="Permalink to this definition">¶</a></dt>
2467 <dd><p>Returns the total mass of all particles.</p>
2468 <dl class="field-list simple">
2469 <dt class="field-odd">Returns</dt>
2470 <dd class="field-odd"><p>The total mass in [kg]</p>
2471 </dd>
2472 </dl>
2473 </dd></dl>
2474
2475 <dl class="method">
2476 <dt id="sphere.sim.totalMomentum">
2477 <code class="sig-name descname">totalMomentum</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.totalMomentum" title="Permalink to this definition">¶</a></dt>
2478 <dd><p>Returns the sum of particle momentums.</p>
2479 <dl class="field-list simple">
2480 <dt class="field-odd">Returns</dt>
2481 <dd class="field-odd"><p>The sum of particle momentums (m*v) [N*s]</p>
2482 </dd>
2483 <dt class="field-even">Return type</dt>
2484 <dd class="field-even"><p>numpy.array</p>
2485 </dd>
2486 </dl>
2487 </dd></dl>
2488
2489 <dl class="method">
2490 <dt id="sphere.sim.totalRotationalEnergy">
2491 <code class="sig-name descname">totalRotationalEnergy</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.totalRotationalEnergy" title="Permalink to this definition">¶</a></dt>
2492 <dd><p>Returns the total rotational kinetic energy for all particles.</p>
2493 <dl class="field-list simple">
2494 <dt class="field-odd">Returns</dt>
2495 <dd class="field-odd"><p>The rotational energy of all particles [J]</p>
2496 </dd>
2497 </dl>
2498 </dd></dl>
2499
2500 <dl class="method">
2501 <dt id="sphere.sim.totalViscousEnergy">
2502 <code class="sig-name descname">totalViscousEnergy</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.totalViscousEnergy" title="Permalink to this definition">¶</a></dt>
2503 <dd><p>Returns the total viscous dissipated energy for all particles.</p>
2504 <dl class="field-list simple">
2505 <dt class="field-odd">Returns</dt>
2506 <dd class="field-odd"><p>The normal viscous energy lost by all particles [J]</p>
2507 </dd>
2508 <dt class="field-even">Return type</dt>
2509 <dd class="field-even"><p>float</p>
2510 </dd>
2511 </dl>
2512 </dd></dl>
2513
2514 <dl class="method">
2515 <dt id="sphere.sim.triaxial">
2516 <code class="sig-name descname">triaxial</code><span class="sig-paren">(</span><em class="sig-param">wvel=-0.001</em>, <em class="sig-param">normal_stress=10000.0</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.triaxial" title="Permalink to this definition">¶</a></dt>
2517 <dd><p>Setup triaxial experiment. The upper wall is moved at a fixed velocity
2518 in m/s, default values is -0.001 m/s (i.e. downwards). The side walls
2519 are exerting a defined normal stress.</p>
2520 <dl class="field-list simple">
2521 <dt class="field-odd">Parameters</dt>
2522 <dd class="field-odd"><ul class="simple">
2523 <li><p><strong>wvel</strong> (<em>float</em>) – Upper wall velocity. Negative values mean that the wall
2524 moves downwards.</p></li>
2525 <li><p><strong>normal_stress</strong> (<em>float</em>) – The normal stress to apply from the upper wall.</p></li>
2526 </ul>
2527 </dd>
2528 </dl>
2529 </dd></dl>
2530
2531 <dl class="method">
2532 <dt id="sphere.sim.uniaxialStrainRate">
2533 <code class="sig-name descname">uniaxialStrainRate</code><span class="sig-paren">(</span><em class="sig-param">wvel=-0.001</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.uniaxialStrainRate" title="Permalink to this definition">¶</a></dt>
2534 <dd><p>Setup consolidation experiment. Specify the upper wall velocity in m/s,
2535 default value is -0.001 m/s (i.e. downwards).</p>
2536 <dl class="field-list simple">
2537 <dt class="field-odd">Parameters</dt>
2538 <dd class="field-odd"><p><strong>wvel</strong> (<em>float</em>) – Upper wall velocity. Negative values mean that the wall
2539 moves downwards.</p>
2540 </dd>
2541 </dl>
2542 </dd></dl>
2543
2544 <dl class="method">
2545 <dt id="sphere.sim.video">
2546 <code class="sig-name descname">video</code><span class="sig-paren">(</span><em class="sig-param">out_folder='./'</em>, <em class="sig-param">video_format='mp4'</em>, <em class="sig-param">graphics_folder='../img_out/'</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">fps=25</em>, <em class="sig-param">qscale=1</em>, <em class="sig-param">bitrate=1800</em>, <em class="sig-param">verbose=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.video" title="Permalink to this definition">¶</a></dt>
2547 <dd><p>Uses ffmpeg to combine images to animation. All images should be
2548 rendered beforehand using <a class="reference internal" href="#sphere.render" title="sphere.render"><code class="xref py py-func docutils literal notranslate"><span class="pre">render()</span></code></a>.</p>
2549 <dl class="field-list simple">
2550 <dt class="field-odd">Parameters</dt>
2551 <dd class="field-odd"><ul class="simple">
2552 <li><p><strong>out_folder</strong> (<em>str</em>) – The output folder for the video file</p></li>
2553 <li><p><strong>video_format</strong> (<em>str</em>) – The format of the output video</p></li>
2554 <li><p><strong>graphics_folder</strong> (<em>str</em>) – The folder containing the rendered images</p></li>
2555 <li><p><strong>graphics_format</strong> (<em>str</em>) – The format of the rendered images</p></li>
2556 <li><p><strong>fps</strong> (<em>int</em>) – The number of frames per second to use in the video</p></li>
2557 <li><p><strong>qscale</strong> (<em>float</em>) – The output video quality, in ]0;1]</p></li>
2558 <li><p><strong>bitrate</strong> (<em>int</em>) – The bitrate to use in the output video</p></li>
2559 <li><p><strong>verbose</strong> (<em>bool</em>) – Show ffmpeg output</p></li>
2560 </ul>
2561 </dd>
2562 </dl>
2563 </dd></dl>
2564
2565 <dl class="method">
2566 <dt id="sphere.sim.viscousEnergy">
2567 <code class="sig-name descname">viscousEnergy</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.viscousEnergy" title="Permalink to this definition">¶</a></dt>
2568 <dd><p>Returns the viscous dissipated energy for a particle.</p>
2569 <dl class="field-list simple">
2570 <dt class="field-odd">Parameters</dt>
2571 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
2572 </dd>
2573 <dt class="field-even">Returns</dt>
2574 <dd class="field-even"><p>The energy lost by the particle by viscous dissipation [J]</p>
2575 </dd>
2576 <dt class="field-odd">Return type</dt>
2577 <dd class="field-odd"><p>float</p>
2578 </dd>
2579 </dl>
2580 </dd></dl>
2581
2582 <dl class="method">
2583 <dt id="sphere.sim.visualize">
2584 <code class="sig-name descname">visualize</code><span class="sig-paren">(</span><em class="sig-param">method='energy'</em>, <em class="sig-param">savefig=True</em>, <em class="sig-param">outformat='png'</em>, <em class="sig-param">figsize=False</em>, <em class="sig-param">pickle=False</em>, <em class="sig-param">xlim=False</em>, <em class="sig-param">firststep=0</em>, <em class="sig-param">f_min=None</em>, <em class="sig-param">f_max=None</em>, <em class="sig-param">cmap=None</em>, <em class="sig-param">smoothing=0</em>, <em class="sig-param">smoothing_window='hanning'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.visualize" title="Permalink to this definition">¶</a></dt>
2585 <dd><p>Visualize output from the simulation, where the temporal progress is
2586 of interest. The output will be saved in the current folder with a name
2587 combining the simulation id of the simulation, and the visualization
2588 method.</p>
2589 <dl class="field-list simple">
2590 <dt class="field-odd">Parameters</dt>
2591 <dd class="field-odd"><ul class="simple">
2592 <li><p><strong>method</strong> (<em>str</em>) – The type of plot to render. Possible values are ‘energy’,
2593 ‘walls’, ‘triaxial’, ‘inertia’, ‘mean-fluid-pressure’,
2594 ‘fluid-pressure’, ‘shear’, ‘shear-displacement’, ‘porosity’,
2595 ‘rate-dependence’, ‘contacts’</p></li>
2596 <li><p><strong>savefig</strong> (<em>bool</em>) – Save the image instead of showing it on screen</p></li>
2597 <li><p><strong>outformat</strong> – The output format of the plot data. This can be an
2598 image format, or in text (‘txt’).</p></li>
2599 <li><p><strong>figsize</strong> (<em>array</em>) – Specify output figure size in inches</p></li>
2600 <li><p><strong>pickle</strong> (<em>bool</em>) – Save all figure content as a Python pickle file. It can
2601 be opened later using <cite>fig=pickle.load(open(‘file.pickle’,’rb’))</cite>.</p></li>
2602 <li><p><strong>xlim</strong> (<em>array</em>) – Set custom limits to the x axis. If not specified, the x
2603 range will correspond to the entire data interval.</p></li>
2604 <li><p><strong>firststep</strong> (<em>int</em>) – The first output file step to read (default: 0)</p></li>
2605 <li><p><strong>cmap</strong> (<em>matplotlib.colors.LinearSegmentedColormap</em>) – Choose custom color map, e.g.
2606 <cite>cmap=matplotlib.cm.get_cmap(‘afmhot’)</cite></p></li>
2607 <li><p><strong>smoothing</strong> (<em>int</em>) – Apply smoothing across a number of output files to the
2608 <cite>method=’shear’</cite> plot. A value of less than 3 means that no
2609 smoothing occurs.</p></li>
2610 <li><p><strong>smoothing_window</strong> (<em>str</em>) – Type of smoothing to use when <cite>smoothing >= 3</cite>.
2611 Valid values are ‘flat’, ‘hanning’ (default), ‘hamming’, ‘bartlett’,
2612 and ‘blackman’.</p></li>
2613 </ul>
2614 </dd>
2615 </dl>
2616 </dd></dl>
2617
2618 <dl class="method">
2619 <dt id="sphere.sim.voidRatio">
2620 <code class="sig-name descname">voidRatio</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.voidRatio" title="Permalink to this definition">¶</a></dt>
2621 <dd><p>Calculates the current void ratio</p>
2622 <dl class="field-list simple">
2623 <dt class="field-odd">Returns</dt>
2624 <dd class="field-odd"><p>The void ratio (pore volume relative to solid volume)</p>
2625 </dd>
2626 <dt class="field-even">Return type</dt>
2627 <dd class="field-even"><p>float</p>
2628 </dd>
2629 </dl>
2630 </dd></dl>
2631
2632 <dl class="method">
2633 <dt id="sphere.sim.volume">
2634 <code class="sig-name descname">volume</code><span class="sig-paren">(</span><em class="sig-param">idx</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.volume" title="Permalink to this definition">¶</a></dt>
2635 <dd><p>Returns the volume of a particle.</p>
2636 <dl class="field-list simple">
2637 <dt class="field-odd">Parameters</dt>
2638 <dd class="field-odd"><p><strong>idx</strong> (<em>int</em>) – Particle index</p>
2639 </dd>
2640 <dt class="field-even">Returns</dt>
2641 <dd class="field-even"><p>The volume of the particle [m^3]</p>
2642 </dd>
2643 <dt class="field-odd">Return type</dt>
2644 <dd class="field-odd"><p>float</p>
2645 </dd>
2646 </dl>
2647 </dd></dl>
2648
2649 <dl class="method">
2650 <dt id="sphere.sim.wall0iz">
2651 <code class="sig-name descname">wall0iz</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.wall0iz" title="Permalink to this definition">¶</a></dt>
2652 <dd><p>Returns the cell index of wall 0 along z.</p>
2653 <dl class="field-list simple">
2654 <dt class="field-odd">Returns</dt>
2655 <dd class="field-odd"><p>z cell index</p>
2656 </dd>
2657 <dt class="field-even">Return type</dt>
2658 <dd class="field-even"><p>int</p>
2659 </dd>
2660 </dl>
2661 </dd></dl>
2662
2663 <dl class="method">
2664 <dt id="sphere.sim.wet">
2665 <code class="sig-name descname">wet</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.wet" title="Permalink to this definition">¶</a></dt>
2666 <dd><p>Set the simulation to be wet (total fluid saturation).</p>
2667 <p>See also <a class="reference internal" href="#sphere.sim.dry" title="sphere.sim.dry"><code class="xref py py-func docutils literal notranslate"><span class="pre">dry()</span></code></a></p>
2668 </dd></dl>
2669
2670 <dl class="method">
2671 <dt id="sphere.sim.writeFluidVTK">
2672 <code class="sig-name descname">writeFluidVTK</code><span class="sig-paren">(</span><em class="sig-param">folder='../output/'</em>, <em class="sig-param">cell_centered=True</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.writeFluidVTK" title="Permalink to this definition">¶</a></dt>
2673 <dd><p>Writes a VTK file for the fluid grid to the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder by
2674 default. The file name will be in the format <code class="docutils literal notranslate"><span class="pre">fluid-<self.sid>.vti</span></code>.
2675 The vti files can be used for visualizing the fluid in ParaView.</p>
2676 <p>The scalars (pressure, porosity, porosity change) and the velocity
2677 vectors are either placed in a grid where the grid corners correspond to
2678 the computational grid center (cell_centered=False). This results in a
2679 grid that doesn’t appears to span the simulation domain, and values are
2680 smoothly interpolated on the cell faces. Alternatively, the
2681 visualization grid is equal to the computational grid, and cells face
2682 colors are not interpolated (cell_centered=True, default behavior).</p>
2683 <p>The fluid grid is visualized by opening the vti files, and pressing
2684 “Apply” to import all fluid field properties. To visualize the scalar
2685 fields, such as the pressure, the porosity, the porosity change or the
2686 velocity magnitude, choose “Surface” or “Surface With Edges” as the
2687 “Representation”. Choose the desired property as the “Coloring” field.
2688 It may be desirable to show the color bar by pressing the “Show” button,
2689 and “Rescale” to fit the color range limits to the current file. The
2690 coordinate system can be displayed by checking the “Show Axis” field.
2691 All adjustments by default require the “Apply” button to be pressed
2692 before regenerating the view.</p>
2693 <p>The fluid vector fields (e.g. the fluid velocity) can be visualizing by
2694 e.g. arrows. To do this, select the fluid data in the “Pipeline
2695 Browser”. Press “Glyph” from the “Common” toolbar, or go to the
2696 “Filters” mennu, and press “Glyph” from the “Common” list. Make sure
2697 that “Arrow” is selected as the “Glyph type”, and “Velocity” as the
2698 “Vectors” value. Adjust the “Maximum Number of Points” to be at least as
2699 big as the number of fluid cells in the grid. Press “Apply” to visualize
2700 the arrows.</p>
2701 <p>To visualize the cell-centered data with smooth interpolation, and in
2702 order to visualize fluid vector fields, the cell-centered mesh is
2703 selected in the “Pipeline Browser”, and is filtered using “Filters” ->
2704 “Alphabetical” -> “Cell Data to Point Data”.</p>
2705 <p>If several data files are generated for the same simulation (e.g. using
2706 the <a class="reference internal" href="#sphere.sim.writeVTKall" title="sphere.sim.writeVTKall"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeVTKall()</span></code></a> function), it is able to step the
2707 visualization through time by using the ParaView controls.</p>
2708 <dl class="field-list simple">
2709 <dt class="field-odd">Parameters</dt>
2710 <dd class="field-odd"><ul class="simple">
2711 <li><p><strong>folder</strong> (<em>str</em>) – The folder where to place the output binary file (default
2712 (default=’../output/’)</p></li>
2713 <li><p><strong>cell_centered</strong> (<em>bool</em>) – put scalars and vectors at cell centers (True) or
2714 cell corners (False), (default=True)</p></li>
2715 <li><p><strong>verbose</strong> (<em>bool</em>) – Show diagnostic information (default=True)</p></li>
2716 </ul>
2717 </dd>
2718 </dl>
2719 </dd></dl>
2720
2721 <dl class="method">
2722 <dt id="sphere.sim.writeVTK">
2723 <code class="sig-name descname">writeVTK</code><span class="sig-paren">(</span><em class="sig-param">folder='../output/'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.writeVTK" title="Permalink to this definition">¶</a></dt>
2724 <dd><p>Writes a VTK file with particle information to the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder
2725 by default. The file name will be in the format <code class="docutils literal notranslate"><span class="pre"><self.sid>.vtu</span></code>.
2726 The vtu files can be used to visualize the particles in ParaView.</p>
2727 <p>After opening the vtu files, the particle fields will show up in the
2728 “Properties” list. Press “Apply” to import all fields into the ParaView
2729 session. The particles are visualized by selecting the imported data in
2730 the “Pipeline Browser”. Afterwards, click the “Glyph” button in the
2731 “Common” toolbar, or go to the “Filters” menu, and press “Glyph” from
2732 the “Common” list. Choose “Sphere” as the “Glyph Type”, choose “scalar”
2733 as the “Scale Mode”. Check the “Edit” checkbox, and set the “Set Scale
2734 Factor” to 1.0. The field “Maximum Number of Points” may be increased if
2735 the number of particles exceed the default value. Finally press “Apply”,
2736 and the particles will appear in the main window.</p>
2737 <p>The sphere resolution may be adjusted (“Theta resolution”, “Phi
2738 resolution”) to increase the quality and the computational requirements
2739 of the rendering. All adjustments by default require the “Apply” button
2740 to be pressed before regenerating the view.</p>
2741 <p>If several vtu files are generated for the same simulation (e.g. using
2742 the <a class="reference internal" href="#sphere.sim.writeVTKall" title="sphere.sim.writeVTKall"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeVTKall()</span></code></a> function), it is able to step the
2743 visualization through time by using the ParaView controls.</p>
2744 <dl class="field-list simple">
2745 <dt class="field-odd">Parameters</dt>
2746 <dd class="field-odd"><ul class="simple">
2747 <li><p><strong>folder</strong> (<em>str</em>) – The folder where to place the output binary file (default
2748 (default=’../output/’)</p></li>
2749 <li><p><strong>verbose</strong> (<em>bool</em>) – Show diagnostic information (default=True)</p></li>
2750 </ul>
2751 </dd>
2752 </dl>
2753 </dd></dl>
2754
2755 <dl class="method">
2756 <dt id="sphere.sim.writeVTKall">
2757 <code class="sig-name descname">writeVTKall</code><span class="sig-paren">(</span><em class="sig-param">cell_centered=True</em>, <em class="sig-param">verbose=True</em>, <em class="sig-param">forces=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.writeVTKall" title="Permalink to this definition">¶</a></dt>
2758 <dd><p>Writes a VTK file for each simulation output file with particle
2759 information and the fluid grid to the <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder by default.
2760 The file name will be in the format <code class="docutils literal notranslate"><span class="pre"><self.sid>.vtu</span></code> and
2761 <code class="docutils literal notranslate"><span class="pre">fluid-<self.sid>.vti</span></code>. The vtu files can be used to visualize the
2762 particles, and the vti files for visualizing the fluid in ParaView.</p>
2763 <p>After opening the vtu files, the particle fields will show up in the
2764 “Properties” list. Press “Apply” to import all fields into the ParaView
2765 session. The particles are visualized by selecting the imported data in
2766 the “Pipeline Browser”. Afterwards, click the “Glyph” button in the
2767 “Common” toolbar, or go to the “Filters” menu, and press “Glyph” from
2768 the “Common” list. Choose “Sphere” as the “Glyph Type”, set “Radius” to
2769 1.0, choose “scalar” as the “Scale Mode”. Check the “Edit” checkbox, and
2770 set the “Set Scale Factor” to 1.0. The field “Maximum Number of Points”
2771 may be increased if the number of particles exceed the default value.
2772 Finally press “Apply”, and the particles will appear in the main window.</p>
2773 <p>The sphere resolution may be adjusted (“Theta resolution”, “Phi
2774 resolution”) to increase the quality and the computational requirements
2775 of the rendering.</p>
2776 <p>The fluid grid is visualized by opening the vti files, and pressing
2777 “Apply” to import all fluid field properties. To visualize the scalar
2778 fields, such as the pressure, the porosity, the porosity change or the
2779 velocity magnitude, choose “Surface” or “Surface With Edges” as the
2780 “Representation”. Choose the desired property as the “Coloring” field.
2781 It may be desirable to show the color bar by pressing the “Show” button,
2782 and “Rescale” to fit the color range limits to the current file. The
2783 coordinate system can be displayed by checking the “Show Axis” field.
2784 All adjustments by default require the “Apply” button to be pressed
2785 before regenerating the view.</p>
2786 <p>The fluid vector fields (e.g. the fluid velocity) can be visualizing by
2787 e.g. arrows. To do this, select the fluid data in the “Pipeline
2788 Browser”. Press “Glyph” from the “Common” toolbar, or go to the
2789 “Filters” mennu, and press “Glyph” from the “Common” list. Make sure
2790 that “Arrow” is selected as the “Glyph type”, and “Velocity” as the
2791 “Vectors” value. Adjust the “Maximum Number of Points” to be at least as
2792 big as the number of fluid cells in the grid. Press “Apply” to visualize
2793 the arrows.</p>
2794 <p>If several data files are generated for the same simulation (e.g. using
2795 the <a class="reference internal" href="#sphere.sim.writeVTKall" title="sphere.sim.writeVTKall"><code class="xref py py-func docutils literal notranslate"><span class="pre">writeVTKall()</span></code></a> function), it is able to step the
2796 visualization through time by using the ParaView controls.</p>
2797 <dl class="field-list simple">
2798 <dt class="field-odd">Parameters</dt>
2799 <dd class="field-odd"><ul class="simple">
2800 <li><p><strong>verbose</strong> (<em>bool</em>) – Show diagnostic information (default=True)</p></li>
2801 <li><p><strong>cell_centered</strong> (<em>bool</em>) – Write fluid values to cell centered positions
2802 (default=true)</p></li>
2803 <li><p><strong>forces</strong> (<em>bool</em>) – Write contact force files (slow) (default=False)</p></li>
2804 </ul>
2805 </dd>
2806 </dl>
2807 </dd></dl>
2808
2809 <dl class="method">
2810 <dt id="sphere.sim.writeVTKforces">
2811 <code class="sig-name descname">writeVTKforces</code><span class="sig-paren">(</span><em class="sig-param">folder='../output/'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.writeVTKforces" title="Permalink to this definition">¶</a></dt>
2812 <dd><p>Writes a VTK file with particle-interaction information to the
2813 <code class="docutils literal notranslate"><span class="pre">../output/</span></code> folder by default. The file name will be in the format
2814 <code class="docutils literal notranslate"><span class="pre"><self.sid>.vtp</span></code>. The vtp files can be used to visualize the
2815 particle interactions in ParaView. First use the “Cell Data to Point
2816 Data” filter, and afterwards show the contact network with the “Tube”
2817 filter.</p>
2818 <dl class="field-list simple">
2819 <dt class="field-odd">Parameters</dt>
2820 <dd class="field-odd"><ul class="simple">
2821 <li><p><strong>folder</strong> (<em>str</em>) – The folder where to place the output file (default
2822 (default=’../output/’)</p></li>
2823 <li><p><strong>verbose</strong> (<em>bool</em>) – Show diagnostic information (default=True)</p></li>
2824 </ul>
2825 </dd>
2826 </dl>
2827 </dd></dl>
2828
2829 <dl class="method">
2830 <dt id="sphere.sim.writebin">
2831 <code class="sig-name descname">writebin</code><span class="sig-paren">(</span><em class="sig-param">folder='../input/'</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.writebin" title="Permalink to this definition">¶</a></dt>
2832 <dd><p>Writes a <code class="docutils literal notranslate"><span class="pre">sphere</span></code> binary file to the <code class="docutils literal notranslate"><span class="pre">../input/</span></code> folder by default.
2833 The file name will be in the format <code class="docutils literal notranslate"><span class="pre"><self.sid>.bin</span></code>.</p>
2834 <p>See also <a class="reference internal" href="#sphere.sim.readbin" title="sphere.sim.readbin"><code class="xref py py-func docutils literal notranslate"><span class="pre">readbin()</span></code></a>.</p>
2835 <dl class="field-list simple">
2836 <dt class="field-odd">Parameters</dt>
2837 <dd class="field-odd"><ul class="simple">
2838 <li><p><strong>folder</strong> (<em>str</em>) – The folder where to place the output binary file</p></li>
2839 <li><p><strong>verbose</strong> (<em>bool</em>) – Show diagnostic information (default=True)</p></li>
2840 </ul>
2841 </dd>
2842 </dl>
2843 </dd></dl>
2844
2845 <dl class="method">
2846 <dt id="sphere.sim.zeroKinematics">
2847 <code class="sig-name descname">zeroKinematics</code><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="headerlink" href="#sphere.sim.zeroKinematics" title="Permalink to this definition">¶</a></dt>
2848 <dd><p>Zero all kinematic parameters of the particles. This function is useful
2849 when output from one simulation is reused in another simulation.</p>
2850 </dd></dl>
2851
2852 </dd></dl>
2853
2854 <dl class="function">
2855 <dt id="sphere.status">
2856 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">status</code><span class="sig-paren">(</span><em class="sig-param">project</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.status" title="Permalink to this definition">¶</a></dt>
2857 <dd><p>Check the status.dat file for the target project, and return the last output
2858 file number.</p>
2859 <dl class="field-list simple">
2860 <dt class="field-odd">Parameters</dt>
2861 <dd class="field-odd"><p><strong>project</strong> (<em>str</em>) – The simulation id of the target project</p>
2862 </dd>
2863 <dt class="field-even">Returns</dt>
2864 <dd class="field-even"><p>The last output file written in the simulation calculations</p>
2865 </dd>
2866 <dt class="field-odd">Return type</dt>
2867 <dd class="field-odd"><p>int</p>
2868 </dd>
2869 </dl>
2870 </dd></dl>
2871
2872 <dl class="function">
2873 <dt id="sphere.thinsectionVideo">
2874 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">thinsectionVideo</code><span class="sig-paren">(</span><em class="sig-param">project</em>, <em class="sig-param">out_folder='./'</em>, <em class="sig-param">video_format='mp4'</em>, <em class="sig-param">fps=25</em>, <em class="sig-param">qscale=1</em>, <em class="sig-param">bitrate=1800</em>, <em class="sig-param">verbose=False</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.thinsectionVideo" title="Permalink to this definition">¶</a></dt>
2875 <dd><p>Uses ffmpeg to combine thin section images to an animation. This function
2876 will implicity render the thin section images beforehand.</p>
2877 <dl class="field-list simple">
2878 <dt class="field-odd">Parameters</dt>
2879 <dd class="field-odd"><ul class="simple">
2880 <li><p><strong>project</strong> (<em>str</em>) – The simulation id of the project to render</p></li>
2881 <li><p><strong>out_folder</strong> (<em>str</em>) – The output folder for the video file</p></li>
2882 <li><p><strong>video_format</strong> (<em>str</em>) – The format of the output video</p></li>
2883 <li><p><strong>fps</strong> (<em>int</em>) – The number of frames per second to use in the video</p></li>
2884 <li><p><strong>qscale</strong> (<em>float</em>) – The output video quality, in ]0;1]</p></li>
2885 <li><p><strong>bitrate</strong> (<em>int</em>) – The bitrate to use in the output video</p></li>
2886 <li><p><strong>verbose</strong> (<em>bool</em>) – Show ffmpeg output</p></li>
2887 </ul>
2888 </dd>
2889 </dl>
2890 </dd></dl>
2891
2892 <dl class="function">
2893 <dt id="sphere.torqueScriptParallel3">
2894 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">torqueScriptParallel3</code><span class="sig-paren">(</span><em class="sig-param">obj1</em>, <em class="sig-param">obj2</em>, <em class="sig-param">obj3</em>, <em class="sig-param">email='adc@geo.au.dk'</em>, <em class="sig-param">email_alerts='ae'</em>, <em class="sig-param">walltime='24:00:00'</em>, <em class="sig-param">queue='qfermi'</em>, <em class="sig-param">cudapath='/com/cuda/4.0.17/cuda'</em>, <em class="sig-param">spheredir='/home/adc/code/sphere'</em>, <em class="sig-param">use_workdir=False</em>, <em class="sig-param">workdir='/scratch'</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.torqueScriptParallel3" title="Permalink to this definition">¶</a></dt>
2895 <dd><p>Create job script for the Torque queue manager for three binaries,
2896 executed in parallel, ideally on three GPUs.</p>
2897 <dl class="field-list simple">
2898 <dt class="field-odd">Parameters</dt>
2899 <dd class="field-odd"><ul class="simple">
2900 <li><p><strong>email</strong> (<em>str</em>) – The e-mail address that Torque messages should be sent to</p></li>
2901 <li><p><strong>email_alerts</strong> (<em>str</em>) – The type of Torque messages to send to the e-mail
2902 address. The character ‘b’ causes a mail to be sent when the
2903 execution begins. The character ‘e’ causes a mail to be sent when
2904 the execution ends normally. The character ‘a’ causes a mail to be
2905 sent if the execution ends abnormally. The characters can be written
2906 in any order.</p></li>
2907 <li><p><strong>walltime</strong> (<em>str</em>) – The maximal allowed time for the job, in the format
2908 ‘HH:MM:SS’.</p></li>
2909 <li><p><strong>queue</strong> (<em>str</em>) – The Torque queue to schedule the job for</p></li>
2910 <li><p><strong>cudapath</strong> (<em>str</em>) – The path of the CUDA library on the cluster compute nodes</p></li>
2911 <li><p><strong>spheredir</strong> (<em>str</em>) – The path to the root directory of sphere on the cluster</p></li>
2912 <li><p><strong>use_workdir</strong> (<em>bool</em>) – Use a different working directory than the sphere folder</p></li>
2913 <li><p><strong>workdir</strong> (<em>str</em>) – The working directory during the calculations, if
2914 <cite>use_workdir=True</cite></p></li>
2915 </ul>
2916 </dd>
2917 <dt class="field-even">Returns</dt>
2918 <dd class="field-even"><p>The filename of the script</p>
2919 </dd>
2920 <dt class="field-odd">Return type</dt>
2921 <dd class="field-odd"><p>str</p>
2922 </dd>
2923 </dl>
2924 <p>See also <code class="xref py py-func docutils literal notranslate"><span class="pre">torqueScript()</span></code></p>
2925 </dd></dl>
2926
2927 <dl class="function">
2928 <dt id="sphere.video">
2929 <code class="sig-prename descclassname">sphere.</code><code class="sig-name descname">video</code><span class="sig-paren">(</span><em class="sig-param">project</em>, <em class="sig-param">out_folder='./'</em>, <em class="sig-param">video_format='mp4'</em>, <em class="sig-param">graphics_folder='../img_out/'</em>, <em class="sig-param">graphics_format='png'</em>, <em class="sig-param">fps=25</em>, <em class="sig-param">verbose=True</em><span class="sig-paren">)</span><a class="headerlink" href="#sphere.video" title="Permalink to this definition">¶</a></dt>
2930 <dd><p>Uses ffmpeg to combine images to animation. All images should be
2931 rendered beforehand using <a class="reference internal" href="#sphere.render" title="sphere.render"><code class="xref py py-func docutils literal notranslate"><span class="pre">render()</span></code></a>.</p>
2932 <dl class="field-list simple">
2933 <dt class="field-odd">Parameters</dt>
2934 <dd class="field-odd"><ul class="simple">
2935 <li><p><strong>project</strong> (<em>str</em>) – The simulation id of the project to render</p></li>
2936 <li><p><strong>out_folder</strong> (<em>str</em>) – The output folder for the video file</p></li>
2937 <li><p><strong>video_format</strong> (<em>str</em>) – The format of the output video</p></li>
2938 <li><p><strong>graphics_folder</strong> (<em>str</em>) – The folder containing the rendered images</p></li>
2939 <li><p><strong>graphics_format</strong> (<em>str</em>) – The format of the rendered images</p></li>
2940 <li><p><strong>fps</strong> (<em>int</em>) – The number of frames per second to use in the video</p></li>
2941 <li><p><strong>qscale</strong> (<em>float</em>) – The output video quality, in ]0;1]</p></li>
2942 <li><p><strong>bitrate</strong> (<em>int</em>) – The bitrate to use in the output video</p></li>
2943 <li><p><strong>verbose</strong> (<em>bool</em>) – Show ffmpeg output</p></li>
2944 </ul>
2945 </dd>
2946 </dl>
2947 </dd></dl>
2948
2949 </div>
2950 </div>
2951
2952
2953 </div>
2954 </div>
2955 </div>
2956 <div class="sphinxsidebar" role="navigation" aria-label="main navigation">
2957 <div class="sphinxsidebarwrapper">
2958 <h3><a href="index.html">Table of Contents</a></h3>
2959 <ul>
2960 <li><a class="reference internal" href="#">Python API</a><ul>
2961 <li><a class="reference internal" href="#sample-usage">Sample usage</a></li>
2962 <li><a class="reference internal" href="#module-sphere">The <code class="docutils literal notranslate"><span class="pre">sphere</span></code> module</a></li>
2963 </ul>
2964 </li>
2965 </ul>
2966
2967 <h4>Previous topic</h4>
2968 <p class="topless"><a href="cfd.html"
2969 title="previous chapter">Fluid simulation and particle-fluid interaction</a></p>
2970 <h4>Next topic</h4>
2971 <p class="topless"><a href="sphere_internals.html"
2972 title="next chapter">sphere internals</a></p>
2973 <div role="note" aria-label="source link">
2974 <h3>This Page</h3>
2975 <ul class="this-page-menu">
2976 <li><a href="_sources/python_api.rst.txt"
2977 rel="nofollow">Show Source</a></li>
2978 </ul>
2979 </div>
2980 <div id="searchbox" style="display: none" role="search">
2981 <h3 id="searchlabel">Quick search</h3>
2982 <div class="searchformwrapper">
2983 <form class="search" action="search.html" method="get">
2984 <input type="text" name="q" aria-labelledby="searchlabel" />
2985 <input type="submit" value="Go" />
2986 </form>
2987 </div>
2988 </div>
2989 <script type="text/javascript">$('#searchbox').show(0);</script>
2990 </div>
2991 </div>
2992 <div class="clearer"></div>
2993 </div>
2994 <div class="related" role="navigation" aria-label="related navigation">
2995 <h3>Navigation</h3>
2996 <ul>
2997 <li class="right" style="margin-right: 10px">
2998 <a href="genindex.html" title="General Index"
2999 >index</a></li>
3000 <li class="right" >
3001 <a href="py-modindex.html" title="Python Module Index"
3002 >modules</a> |</li>
3003 <li class="right" >
3004 <a href="sphere_internals.html" title="sphere internals"
3005 >next</a> |</li>
3006 <li class="right" >
3007 <a href="cfd.html" title="Fluid simulation and particle-fluid interaction"
3008 >previous</a> |</li>
3009 <li class="nav-item nav-item-0"><a href="index.html">sphere 2.15-beta documentation</a> »</li>
3010 </ul>
3011 </div>
3012 <div class="footer" role="contentinfo">
3013 © Copyright 2014, Anders Damsgaard.
3014 Created using <a href="http://sphinx-doc.org/">Sphinx</a> 2.2.0.
3015 </div>
3016 </body>
3017