Skip to content
Snippets Groups Projects
Commit ab6011eb authored by cmaffeo2's avatar cmaffeo2
Browse files

Update

parent c45980ef
No related branches found
No related tags found
No related merge requests found
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
# Basic usage of the `arbdmodel` package # Basic usage of the `arbdmodel` package
We'll learn how to build simple models using the arbdmodel Python package that will run on the Atomic Resolution Brownian Dynamics (ARBD/arbd) simulation engine. ARBD requires a GPU and linux operating system. We'll learn how to build simple models using the arbdmodel Python package that will run on the Atomic Resolution Brownian Dynamics (ARBD/arbd) simulation engine. ARBD requires a GPU and linux operating system.
This tutorial assumes you have a working knowledge of Python, including its object-oriented features, and the numpy package. This tutorial assumes you have a working knowledge of Python, including its object-oriented features, and the numpy package.
Also, a word of caution: at the time of writing ARBD will not prevent you from simulating with poorly selected runtime parameters (e.g. too large a timestep; or too large a decomposition period with too small a pairlist distance; or simply potentials that prescribe unreasonably large forces). Hence it is crucial to have a good understanding of the MD/BD algorithm if you are using it in production settings. Also, a word of caution: at the time of writing ARBD will not prevent you from simulating with poorly selected runtime parameters (e.g. too large a timestep; or too large a decomposition period with too small a pairlist distance; or simply potentials that prescribe unreasonably large forces). Hence it is crucial to have a good understanding of the MD/BD algorithm if you are using it in production settings.
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 1: Create particle types ### Step 1: Create particle types
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
import ipdb # in case we need to debug, after a crash, run `ipdb.pm()` import ipdb # in case we need to debug, after a crash, run `ipdb.pm()`
import numpy as np import numpy as np
from arbdmodel import ParticleType, PointParticle, ArbdModel from arbdmodel import ParticleType, PointParticle, ArbdModel
system_size = 620 system_size = 620
num_particles = 6400 num_particles = 6400
## diffusion of liquid argon: https://www.sciencedirect.com/science/article/pii/0375960171901290 ## diffusion of liquid argon: https://www.sciencedirect.com/science/article/pii/0375960171901290
## units "60 * 3.41 AA * sqrt(119.5 k K / 40 dalton)" "AA**2/ns" ## units "60 * 3.41 AA * sqrt(119.5 k K / 40 dalton)" "AA**2/ns"
argon = ParticleType(name="Ar", argon = ParticleType(name="Ar",
mass=39.95, # dalton mass=39.95, # dalton
damping_coefficient=5000, # per ns damping_coefficient=5000, # per ns
custom_property="my_value", custom_property="my_value",
epsilon=0.177, # kcal_mol epsilon=0.177, # kcal_mol
radius=3.345/2) # Angstroms radius=3.345/2) # Angstroms
print(argon) print(argon)
``` ```
%% Output
<<class 'arbdmodel.ParticleType'> Ar>
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 2: Build a system ### Step 2: Build a system
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
## Generate Nx3 array of random cartesian coordinates ## Generate Nx3 array of random cartesian coordinates
positions = system_size*(np.random.random([num_particles,3])-0.5) positions = system_size*(np.random.random([num_particles,3])-0.5)
## Create a list of point particles located at those positions ## Create a list of point particles located at those positions
points = [PointParticle(type_=argon, name="Ar", position=positions[i]) points = [PointParticle(type_=argon, name="Ar", position=positions[i])
for i in range(num_particles)] for i in range(num_particles)]
model = ArbdModel(points, model = ArbdModel(points,
dimensions=[system_size for i in range(3)], # Ångstroms dimensions=[system_size for i in range(3)], # Ångstroms
timestep = 20e-6, # ns; can be specified below with engine instead timestep = 20e-6, # ns; can be specified below with engine instead
) )
print(model) print(model)
print(model.children[:2]) print(model.children[:2])
``` ```
%% Output
<arbdmodel.ArbdModel object at 0x7fe2d9889580>
[<arbdmodel.PointParticle object at 0x7fe2f051b850>, <arbdmodel.PointParticle object at 0x7fe2f052f100>]
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 3: Describe the interactions between the particles ### Step 3: Describe the interactions between the particles
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
from arbdmodel.interactions import LennardJones from arbdmodel.interactions import LennardJones
lj = LennardJones() lj = LennardJones()
model.add_nonbonded_interaction(lj) model.add_nonbonded_interaction(lj)
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 4: Run the simulation ### Step 4: Run the simulation
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
from arbdmodel import ArbdEngine from arbdmodel import ArbdEngine
engine = ArbdEngine(output_period = 1e2, num_steps = 1e5) engine = ArbdEngine(output_period = 1e2, num_steps = 1e5)
engine.simulate(model, output_name="1-lj", directory='sims/1-argon', engine.simulate(model, output_name="1-lj", directory='sims/1-argon',
num_steps=1e3, gpu=1 num_steps=1e3, gpu=1
) # simulation parameters here override those in constructor ) # simulation parameters here override those in constructor
``` ```
%% Output
arbdmodel: INFO: Running arbd with: /home/cmaffeo2/bin/arbd -g 1 1-lj.bd output/1-lj
–––––––––––––––––––––––––––––––––––––––––––––
| Atomic Resolution Brownian Dynamics |
–––––––––––––––––––––––––––––––––––––––––––––
Found 2 GPU(s)
[0] NVIDIA GeForce RTX 3080 (may timeout) | SM 8.6, 1.71GHz, 7.8GB RAM
[1] NVIDIA GeForce RTX 3080 | SM 8.6, 1.71GHz, 7.8GB RAM
Initializing devices... 1
READER: 0 seed 32472
READER: 1 timestep 2e-05
READER: 2 steps 1000
READER: 3 numberFluct 0 # deprecated
READER: 4 interparticleForce 1 # other values deprecated
READER: 5 fullLongRange 0 # deprecated
READER: 6 temperature 295
READER: 7 ParticleDynamicType Langevin
READER: 8 outputPeriod 100.0
READER: 9 outputEnergyPeriod 100.0
READER: 10 outputFormat dcd
READER: 11 decompPeriod 40
READER: 12 cutoff 50
READER: 13 pairlistDistance 10
READER: 14 origin -310.0 -310.0 -310.0
READER: 15 systemSize 620 620 620
READER: 16 particle Ar
READER: 17 num 6400
READER: 18 mass 39.95
READER: 19 transDamping 5000 5000 5000
READER: 20 gridFile potentials/null.dx
READER: 21 inputParticles 1-lj.particles.txt
READER: 22 tabulatedPotential 1
READER: 23 tabulatedFile 0@0@potentials/1-lj-nb.Ar-Ar.dat
Read config file 1-lj.bd
Reader::getValue(9) 100.0
Value 273a: 100.0
Reader::getValue(9) 100.0
Value 273b: 100.0
Reader::getValue(0) 32472
Reader::getValue(9) 100.0
Value 273.0: 100.0
Parsing seed: 32472
Reader::getValue(1) 2e-05
Reader::getValue(9) 100.0
Value 273.1: 100.0
Parsing timestep: 2e-05
Reader::getValue(2) 1000
Reader::getValue(9) 100.0
Value 273.2: 100.0
Parsing steps: 1000
Reader::getValue(3) 0 # deprecated
Reader::getValue(9) 100.0
Value 273.3: 100.0
Parsing numberFluct: 0 # deprecated
Reader::getValue(4) 1 # other values deprecated
Reader::getValue(9) 100.0
Value 273.4: 100.0
Parsing interparticleForce: 1 # other values deprecated
Reader::getValue(5) 0 # deprecated
Reader::getValue(9) 100.0
Value 273.5: 100.0
Parsing fullLongRange: 0 # deprecated
Reader::getValue(6) 295
Reader::getValue(9) 100.0
Value 273.6: 100.0
Parsing temperature: 295
Reader::getValue(7) Langevin
Reader::getValue(9) 100.0
Value 273.7: 100.0
Parsing ParticleDynamicType: Langevin
Reader::getValue(8) 100.0
Reader::getValue(9) 100.0
Value 273.8: 100.0
Parsing outputPeriod: 100.0
Reader::getValue(9) 100.0
Reader::getValue(9) 100.0
Value 273.9: 100.0
Parsing outputEnergyPeriod: 100.0
Reader::getValue(10) dcd
Reader::getValue(9) 100.0
Value 273.10: 100.0
Parsing outputFormat: dcd
Reader::getValue(11) 40
Reader::getValue(9) 100.0
Value 273.11: 100.0
Parsing decompPeriod: 40
Reader::getValue(12) 50
Reader::getValue(9) 100.0
Value 273.12: 100.0
Parsing cutoff: 50
Reader::getValue(13) 10
Reader::getValue(9) 100.0
Value 273.13: 100.0
Parsing pairlistDistance: 10
Reader::getValue(14) -310.0 -310.0 -310.0
Reader::getValue(9) 100.0
Value 273.14: 100.0
Parsing origin: -310.0 -310.0 -310.0
Reader::getValue(15) 620 620 620
Reader::getValue(9) 100.0
Value 273.15: 100.0
Parsing systemSize: 620 620 620
Reader::getValue(16) Ar
Reader::getValue(9) 100.0
Value 273.16: 100.0
Parsing particle: Ar
Reader::getValue(17) 6400
Reader::getValue(9) 100.0
Value 273.17: 100.0
Parsing num: 6400
Reader::getValue(18) 39.95
Reader::getValue(9) 100.0
Value 273.18: 100.0
Parsing mass: 39.95
Reader::getValue(19) 5000 5000 5000
Reader::getValue(9) 100.0
Value 273.19: 100.0
Parsing transDamping: 5000 5000 5000
Reader::getValue(20) potentials/null.dx
Reader::getValue(9) 100.0
Value 273.20: 100.0
Parsing gridFile: potentials/null.dx
The grid file name potentials/null.dx
potentials/null.dx Reader::getValue(21) 1-lj.particles.txt
Reader::getValue(9) 100.0
Value 273.21: 100.0
Parsing inputParticles: 1-lj.particles.txt
Reader::getValue(22) 1
Reader::getValue(9) 100.0
Value 273.22: 100.0
Parsing tabulatedPotential: 1
Reader::getValue(23) 0@0@potentials/1-lj-nb.Ar-Ar.dat
Reader::getValue(9) 100.0
Value 273.23: 100.0
Parsing tabulatedFile: 0@0@potentials/1-lj-nb.Ar-Ar.dat
6400 particles
0 particles attached to RBs
0 groups
Found 1 particle types.
Loading the potential grids...
Built 0 exclusions.
Cutting off the potential from 48 to 50.
Copying particle grids to GPU 1
Copying particle data to GPU 1
56
Setting up random number generator with seed 32472
Created Cell Decomposition (11, 11, 11)
Loading 1 tabulated non-bonded potentials...
Using 0 non-bonded exclusions
Configuration: 6400 particles | 1 replicas
Step 100 [10.00% complete | 0.474 ms/step | 3645.416 ns/day]
Step 200 [20.00% complete | 0.458 ms/step | 3770.621 ns/day]
Step 300 [30.00% complete | 0.441 ms/step | 3915.526 ns/day]
Step 400 [40.00% complete | 0.457 ms/step | 3778.370 ns/day]
Step 500 [50.00% complete | 0.441 ms/step | 3920.323 ns/day]
Step 600 [60.00% complete | 0.458 ms/step | 3772.020 ns/day]
Step 700 [70.00% complete | 0.441 ms/step | 3914.462 ns/day]
Step 800 [80.00% complete | 0.458 ms/step | 3774.904 ns/day]
Step 900 [90.00% complete | 0.441 ms/step | 3918.989 ns/day]
Final Step: 1000
Total Run Time: 0.74s
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 5: Visualize the results ### Step 5: Visualize the results
We recommend using [VMD](https://www.ks.uiuc.edu/Development/Download/download.cgi?PackageName=VMD) for visualizing the simulations trajectories with: We recommend using [VMD](https://www.ks.uiuc.edu/Development/Download/download.cgi?PackageName=VMD) for visualizing the simulations trajectories with:
```bash ```bash
vmd sims/1-argon/1-lj.psf sims/1-argon/output/1-lj.dcd vmd sims/1-argon/1-lj.psf sims/1-argon/output/1-lj.dcd
``` ```
However, for convenience we have provided commands to use the nglview module, which provides browser-based visualization. However, for convenience we have provided commands to use the nglview module, which provides browser-based visualization.
note: you must enable certain jupyter extensions to use nglview as follows before running jupyter-notebook note: you must enable certain jupyter extensions to use nglview as follows before running jupyter-notebook
```bash ```bash
jupyter-nbextension enable --py --sys-prefix widgetsnbextension jupyter-nbextension enable --py --sys-prefix widgetsnbextension
jupyter-nbextension enable nglview --py --sys-prefix jupyter-nbextension enable nglview --py --sys-prefix
``` ```
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
import MDAnalysis as mda import MDAnalysis as mda
import nglview as nv import nglview as nv
u = mda.Universe("sims/1-argon/1-lj.pdb", "sims/1-argon/output/1-lj.dcd") u = mda.Universe("sims/1-argon/1-lj.pdb", "sims/1-argon/output/1-lj.dcd")
w = nv.show_mdanalysis(u) w = nv.show_mdanalysis(u)
w.clear_representations() w.clear_representations()
w.add_ball_and_stick("all",radius=10) w.add_ball_and_stick("all",radius=10)
w w
``` ```
%% Output
/data/server1/cmaffeo2/miniconda3/lib/python3.8/site-packages/MDAnalysis/coordinates/chemfiles.py:108: DeprecationWarning: distutils Version classes are deprecated. Use packaging.version instead.
MIN_CHEMFILES_VERSION = LooseVersion("0.9")
/data/server1/cmaffeo2/miniconda3/lib/python3.8/site-packages/MDAnalysis/topology/PDBParser.py:317: UserWarning: Element information is missing, elements attribute will not be populated. If needed these can be guessed using MDAnalysis.topology.guessers.
warnings.warn("Element information is missing, elements attribute "
/data/server1/cmaffeo2/miniconda3/lib/python3.8/site-packages/MDAnalysis/topology/guessers.py:146: UserWarning: Failed to guess the mass for the following atom types: A
warnings.warn("Failed to guess the mass for the following atom types: {}".format(atom_type))
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 6: Customize the interactions ### Step 6: Customize the interactions
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
from arbdmodel.interactions import AbstractPotential from arbdmodel.interactions import AbstractPotential
## Clear nonbonded interactions ## Clear nonbonded interactions
model.nonbonded_interactions = [] model.nonbonded_interactions = []
class BuckinghamPotential(AbstractPotential): class BuckinghamPotential(AbstractPotential):
def potential(self, r, types): def potential(self, r, types):
## https://royalsocietypublishing.org/doi/10.1098/rspa.1938.0173 ## https://royalsocietypublishing.org/doi/10.1098/rspa.1938.0173
## optionally could pull information from typeA, typeB ## optionally could pull information from typeA, typeB
typeA, typeB = types typeA, typeB = types
return 143932.65 * ( 1.69 * np.exp( -r/0.273 ) - 102e-4 / r**6 ) return 143932.65 * ( 1.69 * np.exp( -r/0.273 ) - 102e-4 / r**6 )
pot = BuckinghamPotential() pot = BuckinghamPotential()
model.add_nonbonded_interaction( pot ) model.add_nonbonded_interaction( pot )
engine.simulate(model, output_name="2-buck", directory='sims/1-argon', engine.simulate(model, output_name="2-buck", directory='sims/1-argon',
num_steps = 1e3, gpu=1 num_steps = 1e3, gpu=1
) )
``` ```
%% Output
arbdmodel: INFO: Running arbd with: /home/cmaffeo2/bin/arbd -g 1 2-buck.bd output/2-buck
–––––––––––––––––––––––––––––––––––––––––––––
| Atomic Resolution Brownian Dynamics |
–––––––––––––––––––––––––––––––––––––––––––––
Found 2 GPU(s)
[0] NVIDIA GeForce RTX 3080 (may timeout) | SM 8.6, 1.71GHz, 7.8GB RAM
[1] NVIDIA GeForce RTX 3080 | SM 8.6, 1.71GHz, 7.8GB RAM
Initializing devices... 1
READER: 0 seed 63096
READER: 1 timestep 2e-05
READER: 2 steps 1000
READER: 3 numberFluct 0 # deprecated
READER: 4 interparticleForce 1 # other values deprecated
READER: 5 fullLongRange 0 # deprecated
READER: 6 temperature 295
READER: 7 ParticleDynamicType Langevin
READER: 8 outputPeriod 100.0
READER: 9 outputEnergyPeriod 100.0
READER: 10 outputFormat dcd
READER: 11 decompPeriod 40
READER: 12 cutoff 50
READER: 13 pairlistDistance 10
READER: 14 origin -310.0 -310.0 -310.0
READER: 15 systemSize 620 620 620
READER: 16 particle Ar
READER: 17 num 6400
READER: 18 mass 39.95
READER: 19 transDamping 5000 5000 5000
READER: 20 gridFile potentials/null.dx
READER: 21 inputParticles 2-buck.particles.txt
READER: 22 tabulatedPotential 1
READER: 23 tabulatedFile 0@0@potentials/1-lj-nb.Ar-Ar.dat
Read config file 2-buck.bd
Reader::getValue(9) 100.0
Value 273a: 100.0
Reader::getValue(9) 100.0
Value 273b: 100.0
Reader::getValue(0) 63096
Reader::getValue(9) 100.0
Value 273.0: 100.0
Parsing seed: 63096
Reader::getValue(1) 2e-05
Reader::getValue(9) 100.0
Value 273.1: 100.0
Parsing timestep: 2e-05
Reader::getValue(2) 1000
Reader::getValue(9) 100.0
Value 273.2: 100.0
Parsing steps: 1000
Reader::getValue(3) 0 # deprecated
Reader::getValue(9) 100.0
Value 273.3: 100.0
Parsing numberFluct: 0 # deprecated
Reader::getValue(4) 1 # other values deprecated
Reader::getValue(9) 100.0
Value 273.4: 100.0
Parsing interparticleForce: 1 # other values deprecated
Reader::getValue(5) 0 # deprecated
Reader::getValue(9) 100.0
Value 273.5: 100.0
Parsing fullLongRange: 0 # deprecated
Reader::getValue(6) 295
Reader::getValue(9) 100.0
Value 273.6: 100.0
Parsing temperature: 295
Reader::getValue(7) Langevin
Reader::getValue(9) 100.0
Value 273.7: 100.0
Parsing ParticleDynamicType: Langevin
Reader::getValue(8) 100.0
Reader::getValue(9) 100.0
Value 273.8: 100.0
Parsing outputPeriod: 100.0
Reader::getValue(9) 100.0
Reader::getValue(9) 100.0
Value 273.9: 100.0
Parsing outputEnergyPeriod: 100.0
Reader::getValue(10) dcd
Reader::getValue(9) 100.0
Value 273.10: 100.0
Parsing outputFormat: dcd
Reader::getValue(11) 40
Reader::getValue(9) 100.0
Value 273.11: 100.0
Parsing decompPeriod: 40
Reader::getValue(12) 50
Reader::getValue(9) 100.0
Value 273.12: 100.0
Parsing cutoff: 50
Reader::getValue(13) 10
Reader::getValue(9) 100.0
Value 273.13: 100.0
Parsing pairlistDistance: 10
Reader::getValue(14) -310.0 -310.0 -310.0
Reader::getValue(9) 100.0
Value 273.14: 100.0
Parsing origin: -310.0 -310.0 -310.0
Reader::getValue(15) 620 620 620
Reader::getValue(9) 100.0
Value 273.15: 100.0
Parsing systemSize: 620 620 620
Reader::getValue(16) Ar
Reader::getValue(9) 100.0
Value 273.16: 100.0
Parsing particle: Ar
Reader::getValue(17) 6400
Reader::getValue(9) 100.0
Value 273.17: 100.0
Parsing num: 6400
Reader::getValue(18) 39.95
Reader::getValue(9) 100.0
Value 273.18: 100.0
Parsing mass: 39.95
Reader::getValue(19) 5000 5000 5000
Reader::getValue(9) 100.0
Value 273.19: 100.0
Parsing transDamping: 5000 5000 5000
Reader::getValue(20) potentials/null.dx
Reader::getValue(9) 100.0
Value 273.20: 100.0
Parsing gridFile: potentials/null.dx
The grid file name potentials/null.dx
potentials/null.dx Reader::getValue(21) 2-buck.particles.txt
Reader::getValue(9) 100.0
Value 273.21: 100.0
Parsing inputParticles: 2-buck.particles.txt
Reader::getValue(22) 1
Reader::getValue(9) 100.0
Value 273.22: 100.0
Parsing tabulatedPotential: 1
Reader::getValue(23) 0@0@potentials/1-lj-nb.Ar-Ar.dat
Reader::getValue(9) 100.0
Value 273.23: 100.0
Parsing tabulatedFile: 0@0@potentials/1-lj-nb.Ar-Ar.dat
6400 particles
0 particles attached to RBs
0 groups
Found 1 particle types.
Loading the potential grids...
Built 0 exclusions.
Cutting off the potential from 48 to 50.
Copying particle grids to GPU 1
Copying particle data to GPU 1
56
Setting up random number generator with seed 63096
Created Cell Decomposition (11, 11, 11)
Loading 1 tabulated non-bonded potentials...
Using 0 non-bonded exclusions
Configuration: 6400 particles | 1 replicas
Step 100 [10.00% complete | 0.472 ms/step | 3664.511 ns/day]
Step 200 [20.00% complete | 0.458 ms/step | 3776.224 ns/day]
Step 300 [30.00% complete | 0.441 ms/step | 3920.234 ns/day]
Step 400 [40.00% complete | 0.457 ms/step | 3778.205 ns/day]
Step 500 [50.00% complete | 0.441 ms/step | 3921.302 ns/day]
Step 600 [60.00% complete | 0.457 ms/step | 3779.114 ns/day]
Step 700 [70.00% complete | 0.440 ms/step | 3927.451 ns/day]
Step 800 [80.00% complete | 0.614 ms/step | 2816.075 ns/day]
Step 900 [90.00% complete | 0.441 ms/step | 3921.835 ns/day]
Final Step: 1000
Total Run Time: 0.69s
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
import MDAnalysis as mda import MDAnalysis as mda
import nglview as nv import nglview as nv
u = mda.Universe("sims/1-argon/2-buck.pdb", "sims/1-argon/output/2-buck.dcd") u = mda.Universe("sims/1-argon/2-buck.pdb", "sims/1-argon/output/2-buck.dcd")
w = nv.show_mdanalysis(u) w = nv.show_mdanalysis(u)
w.clear_representations() w.clear_representations()
w.add_ball_and_stick("all",radius=10) w.add_ball_and_stick("all",radius=10)
w w
``` ```
%% Output
/data/server1/cmaffeo2/miniconda3/lib/python3.8/site-packages/MDAnalysis/topology/PDBParser.py:317: UserWarning: Element information is missing, elements attribute will not be populated. If needed these can be guessed using MDAnalysis.topology.guessers.
warnings.warn("Element information is missing, elements attribute "
/data/server1/cmaffeo2/miniconda3/lib/python3.8/site-packages/MDAnalysis/topology/guessers.py:146: UserWarning: Failed to guess the mass for the following atom types: A
warnings.warn("Failed to guess the mass for the following atom types: {}".format(atom_type))
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Step 7: Add bonds ### Step 7: Add bonds
Argon is, of course, a noble gas and will not form bonds; this is a demonstration only, and not representative of a physical system Argon is, of course, a noble gas and will not form bonds; this is a demonstration only, and not representative of a physical system
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
from arbdmodel.interactions import HarmonicBond from arbdmodel.interactions import HarmonicBond
bond = HarmonicBond(k=100, # kcal/mol AA**2 bond = HarmonicBond(k=100, # kcal/mol AA**2
r0 = 3 # AA r0 = 3 # AA
) )
## Go through every other index ## Go through every other index
for i in range(0,len(model.children),2): for i in range(0,len(model.children),2):
j = i+1 j = i+1
model.add_bond( model.children[i], model.add_bond( model.children[i],
model.children[j], model.children[j],
bond ) bond )
print("numbonds:", len(model.bonds)) print("numbonds:", len(model.bonds))
engine.simulate(model, output_name="3-bonds", directory='sims/1-argon', engine.simulate(model, output_name="3-bonds", directory='sims/1-argon',
num_steps = 1e3, gpu=1 num_steps = 1e3, gpu=1
) )
``` ```
%% Output
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
<ipython-input-5-c1d63cbeb26c> in <module>
5
6 ## Go through every other index
----> 7 for i in range(0,len(model.children),2):
8 j = i+1
9 model.add_bond( model.children[i],
NameError: name 'model' is not defined
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
import MDAnalysis as mda import MDAnalysis as mda
import nglview as nv import nglview as nv
u = mda.Universe("sims/1-argon/3-bonds.psf","sims/1-argon/3-bonds.pdb", "sims/1-argon/output/3-bonds.dcd") u = mda.Universe("sims/1-argon/3-bonds.psf","sims/1-argon/3-bonds.pdb", "sims/1-argon/output/3-bonds.dcd")
w = nv.show_mdanalysis(u) w = nv.show_mdanalysis(u)
w.clear_representations() w.clear_representations()
w.add_ball_and_stick("index < 10",radius=10) w.add_ball_and_stick("index < 10",radius=10)
w w
``` ```
%% Output
/data/server1/cmaffeo2/miniconda3/lib/python3.8/site-packages/MDAnalysis/coordinates/base.py:892: UserWarning: Reader has no dt information, set to 1.0 ps
warnings.warn("Reader has no dt information, set to 1.0 ps")
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` python ``` python
## Why are there so many long bonds? ## Why are there so many long bonds?
print(bond.range_) ## Distances are in Angstroms print(bond.range_) ## Distances are in Angstroms
r = np.linspace(*bond.range_, 20) r = np.linspace(*bond.range_, 20)
print(r) print(r)
print(bond.potential(r)) ## Energies are kcal/mol print(bond.potential(r)) ## Energies are kcal/mol
## Outside of this range, the bond applies zero force ## Outside of this range, the bond applies zero force
``` ```
%% Output
(0, 50)
[ 0. 2.63157895 5.26315789 7.89473684 10.52631579 13.15789474
15.78947368 18.42105263 21.05263158 23.68421053 26.31578947 28.94736842
31.57894737 34.21052632 36.84210526 39.47368421 42.10526316 44.73684211
47.36842105 50. ]
[4.50000000e+02 6.78670360e+00 2.56094183e+02 1.19792244e+03
2.83227147e+03 5.15914127e+03 8.17853186e+03 1.18904432e+04
1.62948753e+04 2.13918283e+04 2.71813019e+04 3.36632964e+04
4.08378116e+04 4.87048476e+04 5.72644044e+04 6.65164820e+04
7.64610803e+04 8.70981994e+04 9.84278393e+04 1.10450000e+05]
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### ARBD is configured through plain text files ### ARBD is configured through plain text files
Use your favorite editor to examine the files named `sims/1-argon/*.bd`. Use your favorite editor to examine the files named `sims/1-argon/*.bd`.
ARBD uses *key-value* pairs to configure most features. The configuration parser does not use an interpretter, so math and string operations may not be used inside .bd files. For a list of the available keywords, see the documentation that comes with ARBD. ARBD uses *key-value* pairs to configure most features. The configuration parser does not use an interpretter, so math and string operations may not be used inside .bd files. For a list of the available keywords, see the documentation that comes with ARBD.
Also, please note that ARBD is under development; some keywords have been deprecated, and the syntax is subject to evolve. The user-oriented interface of the arbdmodel package is less likely to change, which is one of the reasons we strongly recommend using arbdmodel to set up your simulation systems. Nevertheless, if you want to use ARBD, it's a good idea to familiarize yourself with the configuration file syntax. Also, please note that ARBD is under development; some keywords have been deprecated, and the syntax is subject to evolve. The user-oriented interface of the arbdmodel package is less likely to change, which is one of the reasons we strongly recommend using arbdmodel to set up your simulation systems. Nevertheless, if you want to use ARBD, it's a good idea to familiarize yourself with the configuration file syntax.
%% Cell type:markdown id: tags:
### More to come shortly...
First, we'll introduce some classes helpful for working with polymers, and we'll show the polymer models that ship with the class.
We'll show how beads can experience grid-based potentials.
Then we'll show how to construct simple rigid body models.
%% Cell type:code id: tags:
``` python
```
......
This diff is collapsed.
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment