Download Notebook View in GitHub Open in Google Colab

Applying a Fixed External Force

There are many situations where it is useful to apply an explicitly specified, external force to particles. For example,

• Applying opposite forces to two molecules to pull them apart

• Applying a force to ions to represent an electic field across a membrane

• Applying a time varying force that is calculated by a separate code, such as in QM/MM simulations

This is easily done with a CustomExternalForce. Let’s start by loading a PDB file and creating a System from it.

from openmm.app import *
from openmm import *
from openmm.unit import *

pdb = PDBFile('ala_ala_ala.pdb')
forcefield = ForceField('amber14-all.xml')
system = forcefield.createSystem(pdb.topology, constraints=HBonds)

Now we can create a CustomExternalForce. We want it to apply a specified force (fx, fy, fz) to each particle. Recall that when creating a custom force we specify the energy, and that the force is the negative gradient of the energy. We therefore want the potential energy as a function of the position (x, y, z) to be -fx*x-fy*y-fz*z.

external = CustomExternalForce('-fx*x-fy*y-fz*z')
system.addForce(external)
external.addPerParticleParameter('fx')
external.addPerParticleParameter('fy')
external.addPerParticleParameter('fz')

2

You can apply whatever forces you want to whichever particles you want. For this example, let’s apply opposite forces to the two particles at the ends of the chain to stretch it out.

external.addParticle(0, (-10, 0, 0)*kilojoules_per_mole/nanometer)
external.addParticle(32, (10, 0, 0)*kilojoules_per_mole/nanometer)

1

Now we can run a simulation of it.

integrator = LangevinMiddleIntegrator(300*kelvin, 1/picosecond, 0.004*picoseconds)
simulation = Simulation(pdb.topology, system, integrator)
simulation.context.setPositions(pdb.positions)
simulation.step(10)

Changing the forces at any time is easy. Just update the values of the per-particle parameters, then call updateParametersInContext() to copy the new values over to the running simulation.

external.setParticleParameters(0, 0, (-20, 0, 0)*kilojoules_per_mole/nanometer)
external.updateParametersInContext(simulation.context)