add_image_charges

mdhelper.openmm.system.add_image_charges(system: System, topology: Topology, positions: ndarray[float] | Quantity, temp: float | Quantity, fric: float | Quantity, dt: float | Quantity, *, gamma: float = -1, n_cells: int = 2, nbforce: NonbondedForce = None, cnbforces: dict[CustomNonbondedForce, dict[str, Any]] = None, wall_indices: ndarray[int] = None, exclude: bool = False) → tuple[Quantity, Integrator][source]

Implements the method of image charges for dielectric boundaries. For more information about the method, see Refs. [1], [2], and [3] for perfectly conducting boundaries, and Refs. [4], [5], and [6] otherwise.

Note

The boundaries must be in the \(xy\)-plane and along the \(z\)-axis.

Parameters:
systemopenmm.System

OpenMM molecular system.

topologyopenmm.app.Topology

Topological information about an OpenMM system.

positionsnumpy.ndarray

Positions of the \(N\) particles in the real system.

Shape: \((N,\,3)\).

Reference unit: \(\mathrm{nm}\).

tempfloat or openmm.unit.Quantity

System temperature \(T\).

Reference unit: \(\mathrm{K}\).

fricfloat or openmm.unit.Quantity

Friction coefficient \(\gamma\) that couples the system to the heat bath.

Reference unit: \(\mathrm{ps}^{-1}\).

dtfloat or openmm.unit.Quantity

Integration step size \(\Delta t\).

Reference unit: \(\mathrm{ps}\).

gammafloat, keyword-only, default: -1

Scaled image charge magnitude \(\gamma\).

n_cellsint, keyword-only, default: 2

Total number \(n_\mathrm{cell}\) of cells, with \(1\) real and \(n_\mathrm{cell}-1\) image cells.

Tip

Use \(n_\mathrm{cell}>2\) to instantiate a highly-confined simulation system with periodic dimensions large enough to satisfy the constraint that the force field cutoff be less than half the box size.

Valid values: Even integers greater than \(1\) when \(\gamma=-1\).

nbforceopenmm.NonbondedForce, keyword-only, optional

Standard OpenMM pair potential object implementing the nonbonded Lennard-Jones and Coulomb potentials. For the image charges, the charges \(q\) have their signs flipped, and the LJ parameters \(\sigma\) and \(\epsilon\) are both set to \(0\).

cnbforcesdict, keyword-only, optional

Custom pair potential objects implementing other non-standard pair potentials. The keys are the potentials and the values are the corresponding per-particle parameters to use, which can be provided as follows:

  • None: All per-particle parameters set to \(0\).

  • dict: Information on how to transform the real particle parameters into those for the image charges. The valid (optional) key-value pairs are:

    • "charge": Index (int) where the charge \(q\) information is stored.

    • "zero": Index (int) or indices (list, slice, or numpy.ndarray) of parameters to zero out.

    • "replace": dict containing key-value pairs of indices (int) of parameters to change and their replacement value(s). If the value is an int, all particles receive that value for the current pair potential. If the value is another dict, the new parameter value (value of inner dict) is determined by the current parameter value (key of inner dict).

To prevent unintended behavior, ensure that each parameter index for a custom pair potential is used only once across the keys listed above. If one does appear more than once, the order in which parameter values are updated follows that of the list above.

See the Examples section below for a simple illustration of the possibilities outlined above.

wall_indicesarray-like, keyword-only, optional

Atom indices corresponding to wall particles. If not provided, the wall particles are guessed using the system dimensions.

excludebool, keyword-only, default: False

Specifies whether interactions between a wall particle and all image wall particles should be excluded. If False, only the interaction between a wall particle and its own image is removed.

Tip

If you want accurate forces acting on wall particles, set exclude=True.

Returns:
positionsnumpy.ndarray

Positions of the \(N\) particles in the real system and their images.

Shape: \((2N,\,3)\).

Reference unit: \(\mathrm{nm}\).

integratoropenmm.Integrator

Integrator that simulates a system using Langevin dynamics, with the LFMiddle discretization.

References

[1]

Hautman, J.; Halley, J. W.; Rhee, Y. ‐J. Molecular Dynamics Simulation of Water Beween Two Ideal Classical Metal Walls. J. Chem. Phys. 1989, 91 (1), 467–472. https://doi.org/10.1063/1.457481.

[2]

Dwelle, K. A.; Willard, A. P. Constant Potential, Electrochemically Active Boundary Conditions for Electrochemical Simulation. J. Phys. Chem. C 2019, 123 (39), 24095–24103. https://doi.org/10.1021/acs.jpcc.9b06635.

[3]

Son, C. Y.; Wang, Z.-G. Image-Charge Effects on Ion Adsorption near Aqueous Interfaces. Proc. Natl. Acad. Sci. U.S.A. 2021, 118 (19), e2020615118. https://doi.org/10.1073/pnas.2020615118.

[4]

Dos Santos, A. P.; Girotto, M.; Levin, Y. Simulations of Coulomb Systems with Slab Geometry Using an Efficient 3D Ewald Summation Method. The Journal of Chemical Physics 2016, 144 (14), 144103. https://doi.org/10.1063/1.4945560.

[5]

Dos Santos, A. P.; Girotto, M.; Levin, Y. Simulations of Coulomb Systems Confined by Polarizable Surfaces Using Periodic Green Functions. The Journal of Chemical Physics 2017, 147 (18), 184105. https://doi.org/10.1063/1.4997420.

[6]

Son, C. Y.; Wang, Z.-G. Manuscript in preparation.

Examples

Prepare the cnbforces dict argument for a system with the following custom pair potentials:

  • pair_gauss: Gaussian potential with per-particle parameters type and beta.

    • type is the particle type, and is used to look up a tabulated function for parameter values. As such, the new values will vary and depend on the old ones. For example, for a system with \(2\) types of particles, \(0\) becomes \(2\) and \(1\) becomes \(3\).

    • beta is assumed to be constant, with the value for the real particles differing from that for image charges. For this example, the new value is \(\beta = 42\) (arbitrary units).

  • pair_wca: Weeks–Chander–Andersen potential with per- particle parameters sigma and epsilon, both of which should be set to \(0\) to turn off this potential for image charges.

  • pair_coul_recip: Smeared Coulomb potential with per- particle parameter charge, which should be flipped for image charges, and dummy, a value that should be set to \(0\).

NEW_BETA_VALUE = 42
cnbforces = {
    pair_gauss: {"replace": {0: {0: 2, 1: 3}, 1: NEW_BETA_VALUE}},
    pair_wca: None,
    pair_coul_recip: {"charge": 0, "zero": 1}
}