Basic usage

This section provides a brief overview of how to use the adaptive_sampling package. For more detailed information, please refer to the tutorials and code documentation.

Molecular dynamics (MD) interfaces

The adaptive_sampling package provides interfaces to run molecular dynamics simulation on various levels of theory. These interfaces allow you to run simulations using the sampling and exploration tools provided by the package.

  • AseMD: For running ab initio molecular dynamics simulations using the calculators from the atomic simulation environment (ASE), which provides interfaces to various quantum chemistry packages.

  • AshMD: For running ab initio molecular dynamics simulations using the ASH package, which provides additional support for QM/MM simulations.

  • AdaptiveSamplingOpenMM: Allows performing molecular dynamics simulations using the OpenMM engine.

  • InterfaceMD_2D: Provides 2D potentials for testing of sampling methods.

The ASE interface

The AseMD class allows for running adaptive sampling algorithms using any ASE calculator. Refer to the ASE documentation for more information on available calculators and molecular modeling workflows in ASE. In the below example we will use the tblite package and the GFN2-xTB tight-binding method.

 1 # create an ASE fragment from an XYZ files
 2from ase.io import read
 3frag = read("structure.xyz")
 4
 5# Attach the tblite calculator to the fragment
 6from tblite.ase import TBLite
 7frag.calc = TBLite(method="GFN2-xTB")
 8
 9# Initialize the AseMD object
10from adaptive_sampling.interface import AseMD
11the_md = AseMD(
12    atoms=frag,               # ASE fragment
13    dt=1.0,                   # Timestep in fs
14    thermostat=True,          # Apply an langevin thermostat
15    target_temp=300.0,        # The thermostat temperature
16    scratch_dir="./scratch/", # Location where scratch files of the calculator should be written
17)
18
19# Initialize one or multiple sampling methods from `sampling_tools`
20from adaptive_sampling.sampling_tools import *
21the_bias = ...
22
23# Initialize the MD
24the_md.calc_init(
25    init_momenta="random",   # Method for momenta initialization
26    biaspots=[the_bias, ],   # Set the sampling algorithms
27    init_temp=300.0,         # The initial temperature, if `init_momenta="random"`
28    restart_file=None,       # Filename of restart file, if `init_momenta="read"`
29)
30
31# Ready for running the MD
32the_md.run(nsteps=10000)

The ASH interface

The AshMD class allows for running adaptive sampling algorithms using any ASH theory object. Refer to the ASH documentation for more information on available theory objects and molecular modeling workflows in ASH. In the below example we will perform a QM/MM MD simulation using the xTB method and an AMBER force field.

 1# Create the ASH fragment from a PDB
 2import ash
 3frag = ash.Fragment(
 4    pdbfile="structure.pdb",
 5    charge=0,                           # charge of the QM subsystem
 6    mult=1                              # multiplicity of the QM subsystem
 7)
 8
 9# Create the MM Theory
10mm_theory = ash.OpenMMTheory(
11    cluster_fragment=frag,              # ASH Fragment
12    Amberfiles=True,                    # Use Amber parameter file
13    amberprmtopfile=f"structure.parm7", # The Amber parameter file
14    hydrogenmass=1.0,                   # default: 1.5 (hydrogen mass repartitioning)
15    rigidwater=False,                   # constraints not compatible with QM/MM!
16    autoconstraints=None,               # constraints not compatible with QM/MM!
17    periodic=True,                      # Periodic boundary conditions or not.
18)
19
20# Create the QM Theory
21qm_theory = ash.xTBTheory(
22    xtbmethod="GFN2",                   # The xTB method
23    runmode="inputfile",                # Only "inputfile" supports QM/MM
24)
25
26# Create the QM/MM System
27qm_atoms = [i for i in range(0, 23)]    # Indices of QM atoms
28qmmm_theory = ash.QMMMTheory(
29    qm_theory=qm_theory,                # ASH QM Theory object
30    mm_theory=mm_theory,                # ASH MM Theory object (should be OpenMMTheory)
31    fragment=frag,                      # ASH Fragment
32    embedding="Elstat",                 # QM/MM embedding type
33    qmatoms=qm_atoms,                   # The QM atoms (list of atom indices)
34)
35
36# Initialize the AseMD interface
37from adaptive_sampling.interface.interfaceASH import AshMD
38the_md = AshMD(
39    fragment=frag,                      # ASH fragment
40    calculator=qmmm_theory,             # ASH calculator
41    dt=1.0,                             # Time step in fs
42    thermostat=True,                    # Apply Langevin thermostat
43    target_temp=300.0,                  # The target temperature in Kelvin
44    barostat=True,                      # Apply Monte-Carlo barostat
45    target_pressure=1.0,                # The target pressure in Bar
46    barostat_freq=25,                   # Frequency of updating the barostat
47)
48
49# Initialize one or multiple sampling methods from `sampling_tools`
50from adaptive_sampling.sampling_tools import *
51the_bias = ...
52
53# Initialize the MD
54the_md.calc_init(
55    init_momenta="random",              # Method for momenta initialization
56    biaspots=[the_bias, ],              # Set the sampling algorithms
57    init_temp=300.0,                    # The initial temperature, if `init_momenta="random"`
58    restart_file=None,                  # Filename of restart file, if `init_momenta="read"`
59)
60
61# Ready to run the MD
62the_md.run(nsteps=2)

The OpenMM interface

This section provides a minimal example for using the OpenMM interface together with AMBER style prmtop and crd files. Note that also other file formats are supported by OpenMM. You can initialize the AdaptiveSamplingOpenMM class as follows:

 1from sys import stdout
 2
 3from openmm import *
 4from openmm.app import *
 5from openmm.unit import *
 6
 7from adaptive_sampling.interface import AdaptiveSamplingOpenMM
 8
 9# Setup OpenMM
10prmtop = AmberPrmtopFile(f"structure.prmtop")
11crd = AmberInpcrdFile(f"structure.crd")
12system = prmtop.createSystem(
13    nonbondedMethod=NoCutoff,
14    constraints=HBonds,
15)
16
17# Initialize the `AdaptiveSamplingOpenMM` interface to couple the OpenMM simulation to a bias potential
18the_md = AdaptiveSamplingOpenMM(
19    crd.positions,
20    prmtop.topology,
21    system,
22    dt=2.0,               # timestep in fs
23    equil_temp=300.0,     # temperature of simulation in Kelvin
24    langevin_damping=1.0, # langevin damping in 1/ps
25)
26
27# The OpenMM `simulation` and `integrator` objects are set up internally, but can still be modified by calling `the_md.simulation` or `the_md.integrator`
28the_md.integrator.setConstraintTolerance(0.00001)
29the_md.simulation.reporters.append(DCDReporter('trajectory.dcd', 100))

Before running the MD any importance sampling algorithm from sampling_tools has to be attached to the OpenMM interface:

1from adaptive_sampling.sampling_tools import *
2the_bias = ... # init sampling algorithm
3the_md.set_sampling_algorithm(the_bias)

If you want to apply multiple sampling algorithms, you can specify those as a list:

1from adaptive_sampling.sampling_tools import *
2the_bias1 = ... # init first sampling algorithm
3the_bias2 = ... # init second sampling algorithm, e.g. additional harmonic constraint
4the_md.set_sampling_algorithm([the_bias1, the_bias2])

Now, the MD is ready to run:

1the_md.run(nsteps=500000) # 500000 * 2 fs = 1 ns

Sampling tools

Importance sampling algorithms facilitate the calculation of reaction and activation free energies by sampling molecular transitions. In the adaptive sampling package these are located in the sampling_tools subpackage.

1from adaptive_sampling.sampling_tools import *

Implemented are a wide range of sampling algorithms, including well-tempered metadynamics (WTM) and its successor on-the-fly probability enhanced sampling (OPES), the adaptive biasing force (ABF) method or extended-system based hybrid methods (WTM-eABF, OPES-eABF).

To apply simple harmonic constraints, like for example often used in Umbrella Sampling (US), you can use the Harmonic_Constraint class.

 1from adaptive_sampling.sampling_tools import Harmonic_Constraint
 2the_bias = Harmonic_Constraint(
 3    the_md,
 4    force_constants=100.0,             # In kJ/mol/(CV unit)^2, can also be a list for multiple harmonic constraints
 5    equil_positions=1.0,               # In the unit of the CV, can also be a list for multiple harmonic constraints
 6    colvars=[
 7        ['distance', [0,1]],           # Definition of the Collective Variable
 8        # Additional CVs can be added to the list
 9    ],
10    outputfile='constraints_traj.dat', # The output file
11    output_stride=1,                   # The output stride
12)

The OPES-eABF hybrid sampling method can be applied as follows:

1from adaptive_sampling.sampling_tools import OPESeABF
2the_cv = [
3    # TYPE       IDX    MIN  MAX  BIN WIDTH
4    ['distance', [0,1], 1.0, 5.0, 0.1]
5]
6the_bias = OPESeABF(
7    the_md,
8    the_cv,
9)

Here, default parameters of OPES-eABF are applied, estimating parameters like the coupling width to the extended system and the kernel standard deviation from the first 1000 MD steps, before applying any bias.

For more information on OPES-eABF parameters and other sampling methods visit the corresponding section of the documentation, as well as the Code Documentation.

Exploration tools