A general method is proposed for smoothing torsional potentials when the dihedral angle becomes undefined.

Torsional potentials impose steric or geometric constraints between four bonded particles, e.g., the rotational barriers between conformations of the carbon backbone in atomistic models for hydrocarbons.1,2 Typically, torsional potentials are functions of the dihedral angle $ϕ$ between a plane A containing three particles at {r1, r2, r3} and a plane B containing particle positions {r2, r3, r4} (Fig. 1). The planes A and B have normal vectors $nA=r12×r23$ and $nB=r23×r34$, respectively, where rij = rjri. The dihedral angle between the planes is3

(1)

In most cases, this definition can be straightforwardly applied to compute the torsional energy and forces.4 However, a degenerate state exists when {r1, r2, r3} or {r2, r3, r4} are collinear, and it is no longer possible to uniquely construct the required planes. The dihedral angle and torsional potential are then undefined, and the torsional forces are singular. Near this degenerate state, particles can be subjected to large forces, resulting in numerical integration errors or a complete crash for a molecular dynamics simulation.3 The singular nature of the torsional forces may also lead to unphysical motion or undesired restriction of conformations.

FIG. 1.

Dihedral angle $ϕ$ and bond angles $𝜃ijk$ formed from particles {r1, r2, r3, r4}. The positive orientation of $ϕ$ is indicated according to the IUPAC convention.5 Image rendered using Visual Molecular Dynamics 1.9.3.6

FIG. 1.

Dihedral angle $ϕ$ and bond angles $𝜃ijk$ formed from particles {r1, r2, r3, r4}. The positive orientation of $ϕ$ is indicated according to the IUPAC convention.5 Image rendered using Visual Molecular Dynamics 1.9.3.6

Close modal

Fortunately, such degenerate conformations are rarely, if ever, observed in atomistic models for particles participating in torsional potentials because the bond angle $𝜃ijk$ between them (Fig. 1),3

(2)

is restricted by additional interactions to values far from the collinear conformations. However, some coarse-grained polymer models impose significantly softer interactions,7–9 leading to a finite probability of observing a degenerate conformation, $𝜃ijk = 180°$. It is necessary to somehow modify the torsional potential in such a case in order to obtain physically, well-defined behavior in molecular dynamics simulations.

The steric interactions approximated by a torsional potential are expected to weaken when the end particles are farther apart and should be negligible when $??ijk = 180°$. Singular torsional forces can then be removed by modifying the potential so that forces go smoothly to zero when $𝜃ijk$ approaches $180°$. Other authors have demonstrated the effectiveness of this approach for certain polymer models.10–13 However, the functional forms of the potentials chosen in their works offer little control over the range of bond angles for which the potential is modified, which may be problematic for some geometries.

In this note, we propose a simple, general approach for smoothing a torsional potential as a function of the bond angle. The proposed method can be applied to existing torsional potentials over a controllable smoothing range and gives well-defined behavior in the degenerate case. Motivated by a common method for smoothing truncated pair potentials in molecular simulations, we define a smoothed torsional potential,

(3)

where $V0(ϕ)$ is the unsmoothed potential, Vs is the value to which the potential is smoothed, and s is a smoothing function. The force on particle n is

(4)

where $F0,n(ϕ)$ is the force on particle n in the absence of smoothing and $∇n$ is the gradient with respect to rn. The smoothing function must satisfy the following conditions:

1. $s(cos𝜃s)=1$ at the onset of smoothing (angle $𝜃s$);

2. $s(cos𝜃c)=0$ at the end of smoothing (angle $𝜃c$);

3. the first derivatives of s are zero at $𝜃s$ and $𝜃c$ so that Fn = F0,n at $𝜃s$ and $Fn=0$ at $𝜃c$; and

4. s must decay sufficiently quickly so that V and Fn have well-defined limits for all bond angles.

A suitable polynomial for s is

(5)

where $xs=cos𝜃s$, $xc=cos𝜃c$, and $m≥2$. (In writing s piecewise, we assume without loss of generality that $𝜃c>𝜃s$.) The values of $𝜃s$, $𝜃c$, and m set the range and steepness of smoothing. The gradients of s can be computed by chain rule, $∇ns(cos𝜃ijk)=s′(cos𝜃ijk)∇ncos𝜃ijk$, where

(6)

and

(7)

Here $∇n⊗rij$ denotes the outer product and is equal to the identity matrix I when n = j, −I when n = i, and the zero matrix 0 otherwise. In the degenerate limit, the unsmoothed torsional forces4 diverge as $1∕sin𝜃ijk$. If $𝜃c=0°$ or $180°$, $s′ ∼ (1−cos2𝜃ijk)m−1 ∼ sin2(m−1)𝜃ijk$, and Fn has a well-defined limit for $m≥2$.

The proposed smoothing function was tested for a simple torsional potential, $V0(ϕ) = − cosϕ∕2 + cos2ϕ$. The dihedral angle was initially fixed at $ϕ=120°$ with r1 = (−1, 0, 1), r2 = (−1, 0, 0), r3 = (0, 0, 0), and $r4=(0,−sinϕ,cosϕ)$. The angle $𝜃234$ was then varied from its initial value of $90°$ by rotating r4 around r3 in plane B, causing the dihedral angle to jump to $−60°$ as $𝜃234$ passed through $180°$. Figure 2 shows the z-component of the torsional force on particle 4, F4,z, as a function of $𝜃234$ for different smoothing ranges. As expected, the unsmoothed force diverges at $𝜃234=180°$. The smoothed force is equal to the unsmoothed force until $𝜃s$, at which point F4,z goes continuously to zero at $𝜃234=180°$ before returning to the unsmoothed value at $𝜃s$ for $ϕ=−60°$.

FIG. 2.

Smoothed torsional force in the z-direction on particle 4 as a function of bond angle $𝜃234$ for Vs = 0.5, m = 2, $𝜃c=180°$, and varied $𝜃s$. Inset: unsmoothed potential with initial (circle) and final (triangle) dihedral angles marked.

FIG. 2.

Smoothed torsional force in the z-direction on particle 4 as a function of bond angle $𝜃234$ for Vs = 0.5, m = 2, $𝜃c=180°$, and varied $𝜃s$. Inset: unsmoothed potential with initial (circle) and final (triangle) dihedral angles marked.

Close modal

The chosen form of s is particularly computationally convenient because (1) s and $∇ns$ can be calculated from quantities that are already required for the unsmoothed potential and (2) no additional square-root evaluations are necessary if smoothing is applied symmetrically around $90°$ to eliminate both possible degeneracies (i.e., $𝜃c=180°$ also implies smoothing toward $𝜃c=0°$). We encourage the use of the smoothed torsional potential especially during the development of coarse-grained models, using s as an additional adjustable parameter, in order to obtain correct simulation behavior and physics.

See supplementary material for an example implementation of the smoothed torsional potential.

This work was supported by the Princeton Center for Complex Materials (PCCM), a U.S. National Science Foundation Materials Research Science and Engineering Center (Award No. DMR-1420541). This research is part of the Blue Waters sustained-petascale computing project, which is supported by the National Science Foundation (Award Nos. OCI-0725070 and ACI-1238993) and the state of Illinois. Blue Waters is a joint effort of the University of Illinois at Urbana-Champaign and its National Center for Supercomputing Applications.

1.
W. L.
Jorgensen
,
J. D.
, and
C. J.
Swenson
,
J. Am. Chem. Soc.
106
,
6638
(
1984
).
2.
M. G.
Martin
and
J. I.
Siepmann
,
J. Phys. Chem. B
102
,
2569
(
1998
).
3.
M. P.
Allen
and
D. J.
Tildesley
,
Computer Simulation of Liquids
(
Oxford University Press
,
New York
,
1987
).
4.
A.
Blondel
and
M.
Karplus
,
J. Comput. Chem.
17
,
1132
(
1996
).
5.
IUPAC
, in
Compendium of Chemical Terminology
, 2nd ed., edited by
A. D.
McNaught
and
A.
Wilkinson
(
Blackwell Scientific Publications
,
Oxford
,
1997
).
6.
W.
Humphrey
,
A.
Dalke
, and
K.
Schulten
,
J. Mol. Graphics
14
,
33
(
1996
).
7.
D.
Brown
and
J. H. R.
Clarke
,
Macromolecules
24
,
2075
(
1991
).
8.
P.
Depa
,
C.
Chen
, and
J. K.
Maranas
,
J. Chem. Phys.
134
,
014903
(
2011
).
9.
D. D.
Hsu
,
W.
Xia
,
S. G.
Arturo
, and
S.
Keten
,
J. Chem. Theory Comput.
10
,
2514
(
2014
).
10.
R. E.
Tuzun
,
D. W.
Noid
, and
B. G.
Sumpter
,
Macromol. Theory Simul.
4
,
909
(
1995
).
11.
M.
Bulacu
and
E.
van der Giessen
,
J. Chem. Phys.
123
,
114901
(
2005
).
12.
M.
Bulacu
,
N.
Goga
,
W.
Zhao
,
G.
Rossi
,
L.
Monticelli
,
X.
Periole
,
D. P.
Tieleman
, and
S. J.
Marrink
,
J. Chem. Theory Comput.
9
,
3282
(
2013
).
13.
M.
Bernabei
,
A. J.
Moreno
, and
J.
Colmenero
,
Phys. Rev. Lett.
101
,
255701
(
2008
).