forked from gfrd/gfrd.github.com
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathsimulation.html
More file actions
416 lines (306 loc) · 36.5 KB
/
simulation.html
File metadata and controls
416 lines (306 loc) · 36.5 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
<!DOCTYPE html>
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8"/>
<title>eGFRD simulation</title>
<link rel="stylesheet" href="/css/reset.css" type="text/css"/>
<link rel="stylesheet" href="/css/style.css" type="text/css"/>
<link rel="stylesheet" href="/css/syntax.css" type="text/css"/>
<link rel="stylesheet" href="/css/menu.css" type="text/css"/>
<link rel="stylesheet" href="/css/footer.css" type="text/css"/>
<!--[if lt IE 9]>
<script src="/ie7/IE9.js"></script>
<![endif]-->
</head>
<body>
<div id="globalheader">
<!--googleoff: all-->
<ul id="globalnav">
<li id="gn-home"><a href="/">Home</a></li>
<li id="gn-algorithm"><a href="/algorithm">Algorithm</a></li>
<li id="gn-background"><a href="/background">Background</a></li>
<li id="gn-downloads"><a href="/downloads">Downloads</a></li>
<li id="gn-documentation"><a href="/documentation">Documentation</a></li>
<li id="gn-support"><a href="/support">Support</a></li>
<li class="gn-empty"><a></a></li>
</ul>
<!--googleon: all-->
<div id="lastbutton"> </div>
</div>
<h1>Running a simulation</h1>
<div id="container">
<ul>
<li>Go to the directory containing your simulation script.</li>
</ul>
<div class="highlight"><pre><code class="bash"><span class="nv">$ </span><span class="nb">cd</span> <egfrd_directory>/samples/example
</code></pre>
</div>
<ul>
<li>Run the simulation script writing output to stdout for each step.</li>
</ul>
<div class="highlight"><pre><code class="bash"><span class="nv">$ PYTHONPATH</span><span class="o">=</span>../../ python example.py
</code></pre>
</div>
<ul>
<li>Run the simulation writing output to stdout only if an error occurs.</li>
</ul>
<div class="highlight"><pre><code class="bash"><span class="nv">$ LOGLEVEL</span><span class="o">=</span>ERROR <span class="nv">PYTHONPATH</span><span class="o">=</span>../../ python -O example.py
</code></pre>
</div>
<p>Below you find a simulation
<a href="http://github.com/gfrd/gfrd/tree/develop/samples/example/example.py">script</a>,
with comments explaining each step. More example scripts can be found in the
<a href="http://github.com/gfrd/egfrd/tree/develop/samples">samples</a> directory of the
code.</p>
<div class="highlight"><pre><code class="python"><span class="c">#!/usr/bin/env python</span>
<span class="sd">'''</span>
<span class="sd">To run this script:</span>
<span class="sd">PYTHONPATH=../../ python example.py</span>
<span class="sd">'''</span>
<span class="kn">from</span> <span class="nn">egfrd</span> <span class="kn">import</span> <span class="o">*</span>
<span class="kn">from</span> <span class="nn">vtklogger</span> <span class="kn">import</span> <span class="n">VTKLogger</span>
<span class="sd">'''</span>
<span class="sd">The reaction network:</span>
<span class="sd">A <-> B</span>
<span class="sd">A + B <-> C</span>
<span class="sd">C -> 0</span>
<span class="sd">'''</span>
<span class="k">def</span> <span class="nf">run</span><span class="p">(</span> <span class="p">):</span>
<span class="sd">''' Dimensions.</span>
<span class="sd"> '''</span>
<span class="n">L</span> <span class="o">=</span> <span class="mi">1e-6</span> <span class="c"># Size of simulation box.</span>
<span class="n">D</span> <span class="o">=</span> <span class="mi">1e-12</span> <span class="c"># Diffusion constant.</span>
<span class="n">radius</span> <span class="o">=</span> <span class="mi">2.5e-8</span> <span class="c"># Radius of particles.</span>
<span class="sd">''' Simulator.</span>
<span class="sd"> A cube of dimension (L x L x L) is created, with periodic boundary </span>
<span class="sd"> conditions.</span>
<span class="sd"> '''</span>
<span class="n">s</span> <span class="o">=</span> <span class="n">EGFRDSimulator</span><span class="p">(</span> <span class="n">worldSize</span><span class="o">=</span><span class="n">L</span> <span class="p">)</span>
<span class="sd">''' Adding surfaces. </span>
<span class="sd"> Define the surfaces to add to the simulator.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addPlanarSurface( origin, vectorX, vectorY, Lx, Ly, Lz, name )</span>
<span class="sd"> s.addCylindricalSurface( origin, radius, orientation, size, name )</span>
<span class="sd"> See the docstrings. </span>
<span class="sd"> Note that surface are not allowed to touch or overlap. Surfaces should </span>
<span class="sd"> extend over the whole simulation box, so finite surfaces are not </span>
<span class="sd"> supported.</span>
<span class="sd"> Both methods return the surface they create. Assign it </span>
<span class="sd"> to some variable, so it can be used when adding species.</span>
<span class="sd"> '''</span>
<span class="n">WORLD</span> <span class="o">=</span> <span class="bp">True</span>
<span class="n">MEMBRANE1</span> <span class="o">=</span> <span class="bp">True</span>
<span class="n">MEMBRANE2</span> <span class="o">=</span> <span class="bp">True</span>
<span class="n">DNA</span> <span class="o">=</span> <span class="bp">True</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span><span class="p">:</span>
<span class="n">m1</span> <span class="o">=</span> <span class="n">s</span><span class="o">.</span><span class="n">addPlanarSurface</span><span class="p">(</span><span class="n">origin</span><span class="o">=</span><span class="p">[</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">,</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">,</span> <span class="mi">2</span> <span class="o">*</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">10</span> <span class="p">],</span>
<span class="n">vectorX</span><span class="o">=</span><span class="p">[</span> <span class="mi">1</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="n">vectorY</span><span class="o">=</span><span class="p">[</span> <span class="mi">0</span><span class="p">,</span> <span class="mi">1</span><span class="p">,</span> <span class="mi">0</span> <span class="p">],</span>
<span class="n">Lx</span><span class="o">=</span><span class="p">(</span><span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">),</span>
<span class="n">Ly</span><span class="o">=</span><span class="p">(</span><span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">),</span>
<span class="n">Lz</span><span class="o">=</span><span class="n">radius</span><span class="p">,</span>
<span class="n">name</span><span class="o">=</span><span class="s">'m1'</span><span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE2</span><span class="p">:</span>
<span class="n">m2</span> <span class="o">=</span> <span class="n">s</span><span class="o">.</span><span class="n">addPlanarSurface</span><span class="p">(</span> <span class="n">origin</span><span class="o">=</span><span class="p">[</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">,</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">,</span> <span class="mi">8</span> <span class="o">*</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">10</span> <span class="p">],</span>
<span class="n">vectorX</span><span class="o">=</span><span class="p">[</span> <span class="mi">1</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="n">vectorY</span><span class="o">=</span><span class="p">[</span> <span class="mi">0</span><span class="p">,</span> <span class="mi">1</span><span class="p">,</span> <span class="mi">0</span> <span class="p">],</span>
<span class="n">Lx</span><span class="o">=</span><span class="p">(</span><span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">),</span>
<span class="n">Ly</span><span class="o">=</span><span class="p">(</span><span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">),</span>
<span class="n">Lz</span><span class="o">=</span><span class="n">radius</span><span class="p">,</span>
<span class="n">name</span><span class="o">=</span><span class="s">'m2'</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">DNA</span><span class="p">:</span>
<span class="n">d</span> <span class="o">=</span> <span class="n">s</span><span class="o">.</span><span class="n">addCylindricalSurface</span><span class="p">(</span> <span class="n">origin</span><span class="o">=</span><span class="p">[</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">,</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">,</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">2</span> <span class="p">],</span>
<span class="n">radius</span><span class="o">=</span><span class="n">radius</span><span class="p">,</span>
<span class="n">orientation</span><span class="o">=</span><span class="p">[</span> <span class="mi">0</span><span class="p">,</span> <span class="mi">1</span><span class="p">,</span> <span class="mi">0</span> <span class="p">],</span>
<span class="n">size</span><span class="o">=</span><span class="p">(</span><span class="n">L</span> <span class="o">/</span> <span class="mi">2</span><span class="p">),</span>
<span class="n">name</span><span class="o">=</span><span class="s">'d'</span> <span class="p">)</span>
<span class="sd">''' Defining species.</span>
<span class="sd"> Define the type of particles that can exist.</span>
<span class="sd"> Usage:</span>
<span class="sd"> species = Species( 'name', D, radius )</span>
<span class="sd"> If no D or radius is specified, it has to be specified when adding it to a </span>
<span class="sd"> specific surface.</span>
<span class="sd"> '''</span>
<span class="n">A</span> <span class="o">=</span> <span class="n">Species</span><span class="p">(</span> <span class="s">'A'</span><span class="p">,</span> <span class="n">D</span><span class="p">,</span> <span class="n">radius</span> <span class="p">)</span>
<span class="n">B</span> <span class="o">=</span> <span class="n">Species</span><span class="p">(</span> <span class="s">'B'</span><span class="p">,</span> <span class="n">D</span><span class="p">,</span> <span class="n">radius</span> <span class="p">)</span>
<span class="n">C</span> <span class="o">=</span> <span class="n">Species</span><span class="p">(</span> <span class="s">'C'</span><span class="p">,</span> <span class="n">D</span><span class="p">,</span> <span class="n">radius</span> <span class="p">)</span>
<span class="sd">''' Adding species to 3D space.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addSpecies( species )</span>
<span class="sd"> When adding a species to the system without explicitly assigning the </span>
<span class="sd"> surface it can live on, it is assumed to be in the 'world'. The 'world' is </span>
<span class="sd"> a 3D space, also referred to as the 'bulk' or 'cytoplasm'.</span>
<span class="sd"> '''</span>
<span class="k">if</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">A</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">B</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">C</span> <span class="p">)</span>
<span class="sd">''' Adding species to surfaces.</span>
<span class="sd"> Specify which species can exist on which surface.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addSpecies( species, surface, D, radius )</span>
<span class="sd"> See the docstring.</span>
<span class="sd"> Per surface a different diffusion constant D and radius can be specified. </span>
<span class="sd"> By default the ones for the 'world' are used. </span>
<span class="sd"> Note: if particles of species 'A' should be able to exist in the world as </span>
<span class="sd"> well as on one of the previously added surfaces, then it should be added </span>
<span class="sd"> twice. Ones with and ones without an explicit surface argument.</span>
<span class="sd"> '''</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span><span class="p">:</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">A</span><span class="p">,</span> <span class="n">m1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">B</span><span class="p">,</span> <span class="n">m1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">C</span><span class="p">,</span> <span class="n">m1</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE2</span><span class="p">:</span>
<span class="c"># No species can live on membrane 2.</span>
<span class="k">pass</span>
<span class="k">if</span> <span class="n">DNA</span><span class="p">:</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">A</span><span class="p">,</span> <span class="n">d</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">B</span><span class="p">,</span> <span class="n">d</span><span class="p">,</span> <span class="n">D</span> <span class="o">*</span> <span class="mi">2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addSpecies</span><span class="p">(</span> <span class="n">C</span><span class="p">,</span> <span class="n">d</span><span class="p">,</span> <span class="n">D</span> <span class="o">*</span> <span class="mi">2</span><span class="p">,</span> <span class="n">radius</span> <span class="o">*</span> <span class="mi">2</span> <span class="p">)</span>
<span class="sd">''' Adding reactions in 3D.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addReaction( [reactants], [products], rate )</span>
<span class="sd"> For now: a bimolecular reaction can only have 1 product species.</span>
<span class="sd"> Units for the rate of a bimolecular reaction:</span>
<span class="sd"> [kf] = meters^3 / (molecules * second)</span>
<span class="sd"> (Order of magnitude kf: 1e-18)</span>
<span class="sd"> '''</span>
<span class="n">sigma</span> <span class="o">=</span> <span class="mi">2</span> <span class="o">*</span> <span class="n">radius</span>
<span class="n">D_tot</span> <span class="o">=</span> <span class="n">D</span>
<span class="n">tau</span> <span class="o">=</span> <span class="n">sigma</span> <span class="o">*</span> <span class="n">sigma</span> <span class="o">/</span> <span class="n">D_tot</span>
<span class="c"># Bimolecular. </span>
<span class="n">kf_2</span> <span class="o">=</span> <span class="mi">100</span> <span class="o">*</span> <span class="n">sigma</span> <span class="o">*</span> <span class="n">D_tot</span>
<span class="n">kb_2</span> <span class="o">=</span> <span class="mi">0.1</span> <span class="o">/</span> <span class="n">tau</span>
<span class="c"># Unimolecular.</span>
<span class="n">kf_1</span> <span class="o">=</span> <span class="mi">0.1</span> <span class="o">/</span> <span class="n">tau</span>
<span class="n">kb_1</span> <span class="o">=</span> <span class="mi">0.1</span> <span class="o">/</span> <span class="n">tau</span>
<span class="k">if</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># A <-> B</span>
<span class="c"># A + B <-> C</span>
<span class="c"># C -> 0</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">A</span> <span class="p">],</span> <span class="p">[</span> <span class="n">B</span> <span class="p">],</span> <span class="n">kf_1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">B</span> <span class="p">],</span> <span class="p">[</span> <span class="n">A</span> <span class="p">],</span> <span class="n">kb_1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">A</span><span class="p">,</span> <span class="n">B</span> <span class="p">],</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="n">kf_2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="p">[</span> <span class="n">A</span><span class="p">,</span> <span class="n">B</span> <span class="p">],</span> <span class="n">kb_2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="p">[</span> <span class="p">],</span> <span class="n">kf_1</span> <span class="p">)</span>
<span class="sd">''' Adding reactions on surfaces.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addReaction( [reactants], [products], rate )</span>
<span class="sd"> , where one reactant or product is a tuple: (species, surface).</span>
<span class="sd"> Combinations of species which do not appear together as reactants in any </span>
<span class="sd"> reaction are made repulsive by default on every surface.</span>
<span class="sd"> '''</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span><span class="p">:</span>
<span class="c"># A <-> B</span>
<span class="c"># A + B <-> C</span>
<span class="c"># C -> 0</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">m1</span><span class="p">),</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="n">kf_1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">m1</span><span class="p">),</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="n">kb_1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">m1</span><span class="p">),</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="n">kf_2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">m1</span><span class="p">),</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="n">kb_2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">m1</span><span class="p">),</span> <span class="p">],</span> <span class="p">[</span> <span class="p">],</span> <span class="n">kf_1</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE2</span><span class="p">:</span>
<span class="c"># No particles can live on membrane2.</span>
<span class="k">pass</span>
<span class="k">if</span> <span class="n">DNA</span><span class="p">:</span>
<span class="c"># A <-> B</span>
<span class="c"># A + B <-> C</span>
<span class="c"># C -> 0</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">d</span><span class="p">),</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="n">kf_1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">d</span><span class="p">),</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="n">kb_1</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">d</span><span class="p">),</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="n">kf_2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">A</span><span class="p">,</span> <span class="n">d</span><span class="p">),</span> <span class="p">(</span><span class="n">B</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="n">kb_2</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">d</span><span class="p">),</span> <span class="p">],</span> <span class="p">[</span> <span class="p">],</span> <span class="n">kf_1</span> <span class="p">)</span>
<span class="sd">''' Surface binding reactions.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addReaction( [reactant], [product], rate ) )</span>
<span class="sd"> , where product is a tuple: (species, surface).</span>
<span class="sd"> The reactant species for every surface binding reaction is always a </span>
<span class="sd"> species that can only exist in the 'world', so no surface has to be </span>
<span class="sd"> specified there.</span>
<span class="sd"> If a surface should be absorbive, specify (0, surface) as the product.</span>
<span class="sd"> When no surface binding reaction is defined for a combination of a species </span>
<span class="sd"> and a surface, they are made repulsive by default.</span>
<span class="sd"> '''</span>
<span class="c"># Molecule + surface.</span>
<span class="n">kon</span> <span class="o">=</span> <span class="n">kf_2</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span> <span class="ow">and</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># Species C can bind to the membrane. The membrane is reflective, by </span>
<span class="c"># default, to species A and B.</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="n">kon</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE2</span> <span class="ow">and</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># Membrane 2 absorbs all particles.</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">A</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="mi">0</span><span class="p">,</span> <span class="n">m2</span><span class="p">)</span> <span class="p">],</span> <span class="n">kon</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">B</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="mi">0</span><span class="p">,</span> <span class="n">m2</span><span class="p">)</span> <span class="p">],</span> <span class="n">kon</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="mi">0</span><span class="p">,</span> <span class="n">m2</span><span class="p">)</span> <span class="p">],</span> <span class="n">kon</span> <span class="p">)</span>
<span class="k">pass</span>
<span class="k">if</span> <span class="n">DNA</span> <span class="ow">and</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># Species C can bind to the dna. The dna is reflective, by default, to </span>
<span class="c"># species A and B.</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="n">kon</span> <span class="p">)</span>
<span class="sd">''' Surface unbinding reactions.</span>
<span class="sd"> Usage:</span>
<span class="sd"> s.addReaction( [reactant], [product], k ) )</span>
<span class="sd"> , where reactant is a tuple: (species, surface).</span>
<span class="sd"> Unbinding is a Poissonian reaction, from a surface to the 'world'.</span>
<span class="sd"> The product species for every surface binding reaction is always a </span>
<span class="sd"> species that can only exist in the 'world', so no surface has to be </span>
<span class="sd"> specified there.</span>
<span class="sd"> '''</span>
<span class="c"># Unimolecular.</span>
<span class="n">koff</span> <span class="o">=</span> <span class="n">kb_1</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span> <span class="ow">and</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># Species C can unbind to the 'world'.</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">m1</span><span class="p">)</span> <span class="p">],</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="n">koff</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE2</span> <span class="ow">and</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># No particles can live on membrane2.</span>
<span class="k">pass</span>
<span class="k">if</span> <span class="n">DNA</span> <span class="ow">and</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="c"># Species C can unbind to the 'world'.</span>
<span class="n">s</span><span class="o">.</span><span class="n">addReaction</span><span class="p">(</span> <span class="p">[</span> <span class="p">(</span><span class="n">C</span><span class="p">,</span> <span class="n">d</span><span class="p">)</span> <span class="p">],</span> <span class="p">[</span> <span class="n">C</span> <span class="p">],</span> <span class="n">koff</span> <span class="p">)</span>
<span class="sd">'''</span>
<span class="sd"> Particles.</span>
<span class="sd"> '''</span>
<span class="k">if</span> <span class="n">WORLD</span><span class="p">:</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span> <span class="ow">and</span> <span class="n">MEMBRANE2</span><span class="p">:</span>
<span class="c"># Add world particles inside the two planes.</span>
<span class="c"># Note that a CuboidalRegion is defined by 2 corners.</span>
<span class="n">box1</span> <span class="o">=</span> <span class="n">CuboidalRegion</span><span class="p">(</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="mi">2</span> <span class="o">*</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">10</span> <span class="p">],</span> <span class="p">[</span> <span class="n">L</span><span class="p">,</span> <span class="n">L</span><span class="p">,</span> <span class="mi">8</span> <span class="o">*</span> <span class="n">L</span> <span class="o">/</span> <span class="mi">10</span> <span class="p">]</span> <span class="p">)</span>
<span class="n">s</span><span class="o">.</span><span class="n">throwInParticles</span><span class="p">(</span> <span class="n">C</span><span class="p">,</span> <span class="mi">2</span><span class="p">,</span> <span class="n">box1</span> <span class="p">)</span>
<span class="k">else</span><span class="p">:</span>
<span class="c"># Particles are added to world (3D) by default.</span>
<span class="n">s</span><span class="o">.</span><span class="n">throwInParticles</span><span class="p">(</span> <span class="n">C</span><span class="p">,</span> <span class="mi">2</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE1</span><span class="p">:</span> <span class="n">s</span><span class="o">.</span><span class="n">throwInParticles</span><span class="p">(</span> <span class="n">C</span><span class="p">,</span> <span class="mi">2</span><span class="p">,</span> <span class="n">m1</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">MEMBRANE2</span><span class="p">:</span> <span class="k">pass</span>
<span class="k">if</span> <span class="n">DNA</span><span class="p">:</span> <span class="n">s</span><span class="o">.</span><span class="n">throwInParticles</span><span class="p">(</span> <span class="n">C</span><span class="p">,</span> <span class="mi">2</span><span class="p">,</span> <span class="n">d</span> <span class="p">)</span>
<span class="sd">'''</span>
<span class="sd"> Simulation.</span>
<span class="sd"> '''</span>
<span class="n">s</span><span class="o">.</span><span class="n">initialize</span><span class="p">()</span>
<span class="k">print</span> <span class="n">s</span><span class="o">.</span><span class="n">dumpReactions</span><span class="p">()</span>
<span class="c"># Write vtk files for last 100 steps only. Use VTK or Paraview to</span>
<span class="c"># visualize.</span>
<span class="n">vtklogger</span> <span class="o">=</span> <span class="n">VTKLogger</span><span class="p">(</span> <span class="n">s</span><span class="p">,</span> <span class="s">'run'</span><span class="p">,</span> <span class="mi">100</span> <span class="p">)</span>
<span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="nb">range</span><span class="p">(</span> <span class="mi">200</span> <span class="p">):</span>
<span class="k">try</span><span class="p">:</span>
<span class="n">vtklogger</span><span class="o">.</span><span class="n">log</span><span class="p">()</span>
<span class="n">s</span><span class="o">.</span><span class="n">step</span><span class="p">()</span>
<span class="k">except</span> <span class="ne">RuntimeError</span><span class="p">,</span> <span class="n">message</span><span class="p">:</span>
<span class="k">print</span> <span class="n">message</span>
<span class="k">break</span>
<span class="n">vtklogger</span><span class="o">.</span><span class="n">stop</span><span class="p">()</span>
<span class="n">s</span><span class="o">.</span><span class="n">stop</span><span class="p">(</span> <span class="n">s</span><span class="o">.</span><span class="n">t</span> <span class="p">)</span>
<span class="k">if</span> <span class="n">__name__</span> <span class="o">==</span> <span class="s">'__main__'</span><span class="p">:</span>
<span class="n">run</span><span class="p">(</span> <span class="p">)</span>
</code></pre>
</div>
</div>
<div id="globalfooter">
<p>Copyright © 2012 FOM institute AMOLF</p>
<ul class="piped">
<li><a href="/about">About</a></li>
<li><a href="/credits">Credits</a></li>
</ul>
</div>
</body>
</html>