In [1]:

import readdy
print(readdy.__version__)

v2.0.3-27

Step 1: Set up a reaction diffusion system

In [2]:

system = readdy.ReactionDiffusionSystem([25.,25.,25.], temperature=300.*readdy.units.kelvin)
system.add_species("A", 1.0)

Define reactions

In [3]:

system.reactions.add("myfusion: A +(2) A -> A", rate=10.)
system.reactions.add("myfission: A -> A +(2) A", rate=3.)

Define potentials

In [4]:

system.potentials.add_harmonic_repulsion("A", "A", force_constant=10., 
                                         interaction_distance=2.)

Step 2: Create a corresponding simulation

In [5]:

simulation = system.simulation(kernel="CPU")

simulation.output_file = "out.h5"
simulation.reaction_handler = "UncontrolledApproximation"

simulation.add_particle("A", [0.,0.,0.])

simulation.record_trajectory(stride=1)
simulation.observe.number_of_particles(stride=100, callback=lambda n: print("#A:", n))

Step 3: Run the simulation

In [6]:

simulation.run(n_steps=3000, timestep=1e-2)

  0%|          | 0/300 [00:00<?, ?it/s]

Configured kernel context with:
--------------------------------
 - kBT = 2.49434
 - periodic b.c. = (true, true, true)
 - box size = (25.0, 25.0, 25.0)
 - particle types:
     *  particle type "A" with D=1.0
 - potentials of order 2:
     * for types "A" and "A"
         * Harmonic repulsion with Force constant k=10.0
 - unimolecular reactions:
     * Fission A -> A + A with a rate of 3.0, a product distance of 2.0, and weights 0.5 and 0.5
 - bimolecular reactions:
     * Fusion A + A -> A with a rate of 10.0, an educt distance of 2.0, and weights 0.5 and 0.5

Configured simulation loop with:
--------------------------------
 - timeStep = 0.01
 - evaluateObservables = true
 - progressOutputStride = 100
 - context written to file = true
 - Performing actions:
   * Initialize neighbor list? true
   * Update neighbor list? true
   * Clear neighbor list? true
   * Integrate diffusion? true
   * Calculate forces? true
   * Handle reactions? true
   * Handle topology reactions? true

#A: [1]
#A: [7]

  5%|▌         | 16/300 [00:00<00:01, 153.23it/s]

#A: [16]

 13%|█▎        | 39/300 [00:00<00:02, 116.80it/s]

#A: [36]
#A: [94]

 21%|██        | 63/300 [00:00<00:02, 80.95it/s]

#A: [175]
#A: [242]

 26%|██▌       | 78/300 [00:00<00:03, 68.90it/s]

#A: [328]
#A: [403]

 33%|███▎      | 99/300 [00:01<00:03, 65.20it/s]

#A: [450]
#A: [512]

 40%|████      | 120/300 [00:01<00:02, 60.11it/s]

#A: [529]
#A: [525]

 47%|████▋     | 141/300 [00:01<00:02, 58.78it/s]

#A: [502]
#A: [507]

 54%|█████▎    | 161/300 [00:02<00:02, 61.82it/s]

#A: [475]
#A: [504]

 61%|██████    | 182/300 [00:02<00:01, 63.12it/s]

#A: [476]
#A: [517]

 65%|██████▌   | 196/300 [00:02<00:01, 60.30it/s]

#A: [485]
#A: [469]

 72%|███████▏  | 217/300 [00:03<00:01, 61.38it/s]

#A: [516]
#A: [503]

 79%|███████▉  | 237/300 [00:03<00:01, 56.02it/s]

#A: [487]
#A: [497]

 86%|████████▌ | 257/300 [00:03<00:00, 62.03it/s]

#A: [483]
#A: [531]

 93%|█████████▎| 278/300 [00:04<00:00, 62.94it/s]

#A: [482]
#A: [489]

100%|██████████| 300/300 [00:04<00:00, 65.86it/s]

#A: [496]
#A: [533]

Step 4: Look at results

In [7]:

trajectory = readdy.Trajectory('out.h5')
trajectory.convert_to_xyz(particle_radii={'A': 1.})

Visualize in VMD

In [8]:

!vmd -e out.xyz.tcl

/usr/lib/vmd/vmd_LINUXAMD64: /usr/lib/libGL.so.1: no version information available (required by /usr/lib/vmd/vmd_LINUXAMD64)
Info) VMD for LINUXAMD64, version 1.9.4a38 (October 20, 2019)
Info) http://www.ks.uiuc.edu/Research/vmd/                         
Info) Email questions and bug reports to vmd@ks.uiuc.edu           
Info) Please include this reference in published work using VMD:   
Info)    Humphrey, W., Dalke, A. and Schulten, K., `VMD - Visual   
Info)    Molecular Dynamics', J. Molec. Graphics 1996, 14.1, 33-38.
Info) -------------------------------------------------------------
Info) Multithreading available, 12 CPUs detected.
Info)   CPU features: SSE2 AVX AVX2 FMA 
Info) Free system memory: 25GB (78%)
Info) Creating CUDA device pool and initializing hardware...
Info) Detected 1 available CUDA accelerator::
Info) [0] GeForce GTX 1080     20 SM_6.1 1.9 GHz, 7.9GB RAM SP32 KT AE2 ZC
Warning) Detected X11 'Composite' extension: if incorrect display occurs
Warning) try disabling this X server option.  Most OpenGL drivers
Warning) disable stereoscopic display when 'Composite' is enabled.
Info) OpenGL renderer: GeForce GTX 1080/PCIe/SSE2
Info)   Features: STENCIL MSAA(4) MDE CVA MTX NPOT PP PS GLSL(OVFGS) 
Info)   Full GLSL rendering mode is available.
Info)   Textures: 2-D (32768x32768), 3-D (16384x16384x16384), Multitexture (4)
Info) Detected 1 available TachyonL/OptiX ray tracing accelerator
Info)   Compiling 1 OptiX shaders on 1 target GPU...
Info) Dynamically loaded 3 plugins in directory:
Info) /usr/lib/vmd/plugins/LINUXAMD64/molfile
ERROR) Error opening file out.xyz.tcl
ERROR) 1.0couldn't open "out.xyz.tcl": no such file or directory
vmd > Info) VMD for LINUXAMD64, version 1.9.4a38 (October 20, 2019)
Info) Exiting normally.

Simple example