Initial public commit of this basic MD code

Simulation of interacting particles

3 files changed, 108 insertions, 0 deletions

diff --git a/examples/argon.py b/examples/argon.py
new file mode 100644
index 0000000..6eead38
--- /dev/null
+++ b/examples/argon.py
@@ -0,0 +1,39 @@
+from interacticle import MoleculeCollection, LennardJones, Coulomb, Simulation
+
+import interacticle.visualizer
+from interacticle.visual import WireBox, MolecularLinks, VelocityHistogram, TemperaturePlot
+
+from library import Argon
+
+setup = MoleculeCollection()
+
+setup.potential(LennardJones(39.948, 0.3395, 0.971))
+
+for x in range(10):
+    for y in range(10):
+        for z in range(10):
+            setup.add(Argon, (2+0.4*x,2+0.4*y,2+0.4*z))
+
+def target_temperature(step, temperature):
+    if step < 4e4: # equilibration for 20 ps
+        return None
+    elif step >= 4e4 and step <= 2e5: # Heat up to 300 Kelvin in 20 ps and keep there until t=100ps
+        return min(max(temperature, (step-4e4) / 4e4 * 300), 300)
+    else: # Cool down to 10 Kelvin within 20 ps
+        return max(min(temperature, (4e4 - (step-2e5)) / 4e4 * 300), 10)
+
+setup.domain_size = 8
+setup.tau = 0.0005
+setup.cutoff = 0.3395*2.5
+setup.target_temperature = target_temperature
+
+setup.max_lennard_jones = 200
+setup.max_coulomb = 200
+
+simulation = Simulation(setup, opengl = True)
+simulation.verbose = False
+
+histogram = VelocityHistogram(simulation, [1.1,0.5], [1,0.5])
+temperature = TemperaturePlot(simulation, [1.1,0], [1,0.5])
+
+interacticle.visualizer.simulate(setup, simulation, [ WireBox(setup.domain_size), MolecularLinks(simulation), histogram, temperature ], 100)
diff --git a/examples/dimer_only.py b/examples/dimer_only.py
new file mode 100644
index 0000000..c58eb39
--- /dev/null
+++ b/examples/dimer_only.py
@@ -0,0 +1,36 @@
+from interacticle import MoleculeCollection, LennardJones, Coulomb, Simulation
+
+import interacticle.visualizer
+from interacticle.visual import WireBox, MolecularLinks
+
+from library import Benzene
+
+setup = MoleculeCollection()
+
+setup.potential(LennardJones(12.011, 0.3550, 0.070))
+setup.potential(LennardJones( 1.008, 0.2420, 0.030))
+
+setup.potential(Coulomb(12.011, -0.115))
+setup.potential(Coulomb( 1.008,  0.115))
+
+# T-shaped dimer
+setup.add(Benzene, (1,1,1.519), rotations=[([0,1,0], 15), ([1,0,0], -9), ([1,0,0], 90), ([0,1,0], 90)])
+setup.add(Benzene, (1,1,1.0), rotations=[([0,1,0], 15), ([1,0,0], -9), ([0,0,1], 30)])
+
+# Sandwich-shaped dimer
+#setup.add(Benzene, (1,1,1), rotations=[([0,1,0], 15), ([1,0,0], -9), ([1,0,0], 90), ([0,1,0], 90)])
+#setup.add(Benzene, (1,1.377,1), rotations=[([0,1,0], 15), ([1,0,0], -9), ([1,0,0], 90), ([0,1,0], 60)])
+
+# Parallel displaced dimer
+#setup.add(Benzene, (1,1,1), rotations=[([0,1,0], 15), ([1,0,0], -9), ([1,0,0], 90), ([0,1,0], 90)])
+#setup.add(Benzene, (1.5,1.3,1), rotations=[([0,1,0], 15), ([1,0,0], -9), ([1,0,0], 90), ([0,1,0], 90)])
+
+setup.domain_size = 2
+setup.tau = 0.0005
+setup.cutoff = 0.242*2.5
+setup.statistics_step = 100
+setup.neighborhood_step = 100
+
+simulation = Simulation(setup, opengl = True)
+
+interacticle.visualizer.simulate(setup, simulation, [ WireBox(setup.domain_size), MolecularLinks(simulation) ], 100)
diff --git a/examples/water.py b/examples/water.py
new file mode 100644
index 0000000..64d62d3
--- /dev/null
+++ b/examples/water.py
@@ -0,0 +1,33 @@
+from interacticle import MoleculeCollection, LennardJones, Coulomb, Simulation
+
+import interacticle.visualizer
+from interacticle.visual import WireBox, MolecularLinks, VelocityHistogram, TemperaturePlot
+
+from library import Water
+
+setup = MoleculeCollection()
+
+setup.potential(LennardJones(15.9992, 0.31776, 0.792324))
+setup.potential(Coulomb(15.9992, -0.8450))
+setup.potential(Coulomb( 1.0079,  0.4225))
+
+for x in range(4):
+    for y in range(4):
+        for z in range(4):
+            setup.add(Water, (0.2+0.7*x,0.2+0.7*y,0.2+0.7*z))
+
+setup.domain_size = 5
+setup.tau = 0.0001
+setup.cutoff = 2.6
+
+setup.max_lennard_jones = 500
+setup.max_coulomb = 500
+
+setup.statistics_step = 100
+setup.neighborhood_step = 500
+
+simulation = Simulation(setup, opengl = True)
+
+histogram = VelocityHistogram(simulation, [1.1,0.5], [1,0.5])
+
+interacticle.visualizer.simulate(setup, simulation, [ WireBox(setup.domain_size), MolecularLinks(simulation), histogram ], 100)