Class CustomIntegrator
- Direct Known Subclasses:
CustomMTSIntegrator
,CustomMTSLangevinIntegrator
To create an integration algorithm, you first define a set of variables the integrator will compute. Variables come in two types: global variables have a single value, while per-DOF variables have a value for every degree of freedom (x, y, or z coordinate of a particle). You can define as many variables as you want of each type. The value of any variable can be computed by the integration algorithm, or set directly by calling a method on the CustomIntegrator. All variables are persistent between integration steps; once a value is set, it keeps that value until it is changed by the user or recomputed in a later integration step.
Next, you define the algorithm as a series of computations. To execute a time step, the integrator performs the list of computations in order. Each computation updates the value of one global or per-DOF value. There are several types of computations that can be done:
- Global: You provide a mathematical expression involving only global variables. It is evaluated and stored into a global variable.
- Per-DOF: You provide a mathematical expression involving both global and per-DOF variables. It is evaluated once for every degree of freedom, and the values are stored into a per-DOF variable.
- Sum: You provide a mathematical expression involving both global and per-DOF variables. It is evaluated once for every degree of freedom. All of those values are then added together, and the sum is stored into a global variable.
- Constrain Positions: The particle positions are updated so that all distance constraints are satisfied.
- Constrain Velocities: The particle velocities are updated so the net velocity along any constrained distance is 0.
Like all integrators, CustomIntegrator ignores any particle whose mass is 0. It is skipped when doing per-DOF computations, and is not included when computing sums over degrees of freedom.
In addition to the variables you define by calling addGlobalVariable() and addPerDofVariable(), the integrator provides the following pre-defined variables:
- dt: (global) This is the step size being used by the integrator.
- energy: (global, read-only) This is the current potential energy of the system.
- energy0, energy1, energy2, ...: (global, read-only) This is similar to energy, but includes only the contribution from forces in one force group. A single computation step may only depend on a single energy variable (energy, energy0, energy1, etc.).
- x: (per-DOF) This is the current value of the degree of freedom (the x, y, or z coordinate of a particle).
- v: (per-DOF) This is the current velocity associated with the degree of freedom (the x, y, or z component of a particle's velocity).
- f: (per-DOF, read-only) This is the current force acting on the degree of freedom (the x, y, or z component of the force on a particle).
- f0, f1, f2, ...: (per-DOF, read-only) This is similar to f, but includes only the contribution from forces in one force group. A single computation step may only depend on a single force variable (f, f0, f1, etc.).
- m: (per-DOF, read-only) This is the mass of the particle the degree of freedom is associated with.
- uniform: (either global or per-DOF, read-only) This is a uniformly distributed random number between 0 and 1. Every time an expression is evaluated, a different value will be used. When used in a per-DOF expression, a different value will be used for every degree of freedom. Note, however, that if this variable appears multiple times in a single expression, the same value is used everywhere it appears in that expression.
- gaussian: (either global or per-DOF, read-only) This is a Gaussian distributed random number with mean 0 and variance 1. Every time an expression is evaluated, a different value will be used. When used in a per-DOF expression, a different value will be used for every degree of freedom. Note, however, that if this variable appears multiple times in a single expression, the same value is used everywhere it appears in that expression.
- A global variable is created for every adjustable parameter defined in the integrator's Context.
The following example uses a CustomIntegrator to implement a velocity Verlet integrator:
CustomIntegrator integrator = new CustomIntegrator(0.001);
integrator.addComputePerDof("v", "v+0.5*dt*f/m");
integrator.addComputePerDof("x", "x+dt*v");
integrator.addComputePerDof("v", "v+0.5*dt*f/m");
The first step updates the velocities based on the current forces. The second step updates the positions based on the new velocities, and the third step updates the velocities again. Although the first and third steps look identical, the forces used in them are different. You do not need to tell the integrator that; it will recognize that the positions have changed and know to recompute the forces automatically.
-
Field Summary
Fields inherited from class ffx.openmm.Integrator
pointer
-
Constructor Summary
Constructors -
Method Summary
Modifier and TypeMethodDescriptionint
addComputeGlobal
(String variable, String expression) Add a computation that computes a global value.int
addComputePerDof
(String variable, String expression) Add a per-DOF computation to this Integrator.int
addComputeSum
(String variable, String expression) Add a computation that computes a sum over degrees of freedom.int
Add a position constraint to this Integrator.int
Add a velocity constraint to this Integrator.int
addGlobalVariable
(String name, double initialValue) Add a global variable to this Integrator.int
addPerDofVariable
(String name, double initialValue) Add a per-DOF variable to this Integrator.int
addTabulatedFunction
(String name, com.sun.jna.ptr.PointerByReference function) Add a tabulated function that may appear in expressions.int
Add an update context state to this Integrator.int
beginIfBlock
(String condition) Begin a new if block.int
beginWhileBlock
(String condition) Begin a new while block.void
destroy()
Destroy the integrator.int
endBlock()
End the most recently begun if or while block.void
getComputationStep
(int index, com.sun.jna.ptr.IntByReference type, com.sun.jna.ptr.PointerByReference variable, com.sun.jna.ptr.PointerByReference expression) Get information about a computation in the integrator.void
getComputationStep
(int index, IntBuffer type, com.sun.jna.ptr.PointerByReference variable, com.sun.jna.ptr.PointerByReference expression) Get information about a computation in the integrator.double
getGlobalVariable
(int index) Get the value of a global variable.double
Get the value of a global variable, specified by name.getGlobalVariableName
(int index) Get the name of a global variable.Get the expression used to compute the kinetic energy.int
Get the number of computation steps defined for this integrator.int
Get the number of global variables that have been defined.int
Get the number of per-DOF variables that have been defined.int
Get the number of tabulated functions that have been defined.void
getPerDofVariable
(int index, com.sun.jna.ptr.PointerByReference variable) Get the values of a per-DOF variable.void
getPerDofVariableByName
(String name, com.sun.jna.ptr.PointerByReference variable) Get the values of a per-DOF variable, specified by name.getPerDofVariableName
(int index) Get the name of a per-DOF variable.int
Get the random number seed.com.sun.jna.ptr.PointerByReference
getTabulatedFunction
(int index) Get a reference to a tabulated function that may appear in expressions.getTabulatedFunctionName
(int index) Get the name of a tabulated function that may appear in expressions.void
setGlobalVariable
(int index, double value) Set the value of a global variable.void
setGlobalVariableByName
(String name, double value) Set the value of a global variable, specified by name.void
setKineticEnergyExpression
(String expression) Set the expression used to compute the kinetic energy.void
setPerDofVariable
(int index, com.sun.jna.ptr.PointerByReference variable) Set the values of a per-DOF variable.void
setPerDofVariableByName
(String name, com.sun.jna.ptr.PointerByReference variable) Set the values of a per-DOF variable, specified by name.void
setRandomNumberSeed
(int seed) Set the random number seed.void
step
(int steps) Advance a simulation through time by taking a series of time steps.Methods inherited from class ffx.openmm.Integrator
getConstraintTolerance, getIntegrationForceGroups, getPointer, getStepSize, setConstraintTolerance, setIntegrationForceGroups, setPointer, setStepSize
-
Constructor Details
-
CustomIntegrator
public CustomIntegrator(double dt) Constructor.- Parameters:
dt
- The time step.
-
-
Method Details
-
addComputeGlobal
Add a computation that computes a global value.- Parameters:
variable
- The name of the variable to store the computed value in.expression
- The expression to evaluate.- Returns:
- The index of the computation that was added.
-
addComputePerDof
Add a per-DOF computation to this Integrator.- Parameters:
variable
- The name of the variable to store the computed value in.expression
- The expression to evaluate.- Returns:
- The index of the computation that was added.
-
addComputeSum
Add a computation that computes a sum over degrees of freedom.- Parameters:
variable
- The name of the variable to store the computed value in.expression
- The expression to evaluate and sum.- Returns:
- The index of the computation that was added.
-
addConstrainPositions
public int addConstrainPositions()Add a position constraint to this Integrator.- Returns:
- The index of the computation that was added.
-
addConstrainVelocities
public int addConstrainVelocities()Add a velocity constraint to this Integrator.- Returns:
- The index of the computation that was added.
-
addGlobalVariable
Add a global variable to this Integrator.- Parameters:
name
- The name of the variable to create.initialValue
- The initial value of the variable.- Returns:
- The index of the variable that was added.
-
addPerDofVariable
Add a per-DOF variable to this Integrator.- Parameters:
name
- The name of the variable to create.initialValue
- The initial value of the variable.- Returns:
- The index of the variable that was added.
-
addTabulatedFunction
Add a tabulated function that may appear in expressions.- Parameters:
name
- The name of the function as it appears in expressions.function
- A TabulatedFunction object defining the function.- Returns:
- The index of the function that was added.
-
addUpdateContextState
public int addUpdateContextState()Add an update context state to this Integrator.- Returns:
- The index of the computation that was added.
-
beginIfBlock
Begin a new if block.- Parameters:
condition
- The condition under which the block should be executed.- Returns:
- The index of the computation that was added.
-
beginWhileBlock
Begin a new while block.- Parameters:
condition
- The condition under which the block should be executed repeatedly.- Returns:
- The index of the computation that was added.
-
destroy
public void destroy()Destroy the integrator.- Specified by:
destroy
in classIntegrator
-
endBlock
public int endBlock()End the most recently begun if or while block.- Returns:
- The index of the computation that was added.
-
getComputationStep
public void getComputationStep(int index, com.sun.jna.ptr.IntByReference type, com.sun.jna.ptr.PointerByReference variable, com.sun.jna.ptr.PointerByReference expression) Get information about a computation in the integrator.- Parameters:
index
- The index of the computation to get.type
- The type of the computation (output).variable
- The variable the computation is stored in (output).expression
- The expression to evaluate (output).
-
getComputationStep
public void getComputationStep(int index, IntBuffer type, com.sun.jna.ptr.PointerByReference variable, com.sun.jna.ptr.PointerByReference expression) Get information about a computation in the integrator.- Parameters:
index
- The index of the computation to get.type
- The type of the computation (output).variable
- The variable the computation is stored in (output).expression
- The expression to evaluate (output).
-
getGlobalVariable
public double getGlobalVariable(int index) Get the value of a global variable.- Parameters:
index
- The index of the variable to get.- Returns:
- The value of the variable.
-
getGlobalVariableByName
Get the value of a global variable, specified by name.- Parameters:
name
- The name of the variable to get.- Returns:
- The value of the variable.
-
getGlobalVariableName
Get the name of a global variable.- Parameters:
index
- The index of the variable to get.- Returns:
- The name of the variable.
-
getKineticEnergyExpression
Get the expression used to compute the kinetic energy.- Returns:
- The expression used to compute the kinetic energy.
-
getNumComputations
public int getNumComputations()Get the number of computation steps defined for this integrator.- Returns:
- The number of computation steps.
-
getNumGlobalVariables
public int getNumGlobalVariables()Get the number of global variables that have been defined.- Returns:
- The number of global variables.
-
getNumPerDofVariables
public int getNumPerDofVariables()Get the number of per-DOF variables that have been defined.- Returns:
- The number of per-DOF variables.
-
getNumTabulatedFunctions
public int getNumTabulatedFunctions()Get the number of tabulated functions that have been defined.- Returns:
- The number of tabulated functions.
-
getPerDofVariable
public void getPerDofVariable(int index, com.sun.jna.ptr.PointerByReference variable) Get the values of a per-DOF variable.- Parameters:
index
- The index of the variable to get.variable
- The values of the variable (output).
-
getPerDofVariableByName
Get the values of a per-DOF variable, specified by name.- Parameters:
name
- The name of the variable to get.variable
- The values of the variable (output).
-
getPerDofVariableName
Get the name of a per-DOF variable.- Parameters:
index
- The index of the variable to get.- Returns:
- The name of the variable.
-
getRandomNumberSeed
public int getRandomNumberSeed()Get the random number seed. See setRandomNumberSeed() for details.- Returns:
- The random number seed.
-
getTabulatedFunction
public com.sun.jna.ptr.PointerByReference getTabulatedFunction(int index) Get a reference to a tabulated function that may appear in expressions.- Parameters:
index
- The index of the function to get.- Returns:
- The TabulatedFunction object defining the function.
-
getTabulatedFunctionName
Get the name of a tabulated function that may appear in expressions.- Parameters:
index
- The index of the function to get.- Returns:
- The name of the function as it appears in expressions.
-
setGlobalVariable
public void setGlobalVariable(int index, double value) Set the value of a global variable.- Parameters:
index
- The index of the variable to set.value
- The new value of the variable.
-
setGlobalVariableByName
Set the value of a global variable, specified by name.- Parameters:
name
- The name of the variable to set.value
- The new value of the variable.
-
setKineticEnergyExpression
Set the expression used to compute the kinetic energy.- Parameters:
expression
- The expression used to compute the kinetic energy.
-
setPerDofVariable
public void setPerDofVariable(int index, com.sun.jna.ptr.PointerByReference variable) Set the values of a per-DOF variable.- Parameters:
index
- The index of the variable to set.variable
- The new values of the variable.
-
setPerDofVariableByName
Set the values of a per-DOF variable, specified by name.- Parameters:
name
- The name of the variable to set.variable
- The new values of the variable.
-
setRandomNumberSeed
public void setRandomNumberSeed(int seed) Set the random number seed. The precise meaning of this parameter is undefined, and is left up to each Platform to interpret in an appropriate way. It is guaranteed that if two simulations are run with different random number seeds, the sequence of random numbers will be different. On the other hand, no guarantees are made about the behavior of simulations that use the same seed. In particular, Platforms are permitted to use non-deterministic algorithms which produce different results on successive runs, even if those runs were initialized identically.If seed is set to 0 (which is the default value assigned), a unique seed is chosen when a Context is created from this Force. This is done to ensure that each Context receives unique random seeds without you needing to set them explicitly.
- Parameters:
seed
- The random number seed.
-
step
public void step(int steps) Advance a simulation through time by taking a series of time steps.- Overrides:
step
in classIntegrator
- Parameters:
steps
- The number of time steps to take.
-