Coverage for cvpack/torsion_similarity.py: 100%

1""" 

2.. class:: TorsionSimilarity 

3 :platform: Linux, MacOS, Windows 

4 :synopsis: The degree of similarity between pairs of torsion angles. 

5 

6.. classauthor:: Charlles Abreu <craabreu@gmail.com> 

7 

8""" 

9 

10from __future__ import annotations 

11 

12import numpy as np 

13import openmm 

14from numpy.typing import ArrayLike 

15from openmm import unit as mmunit 

16 

17from .collective_variable import CollectiveVariable 

18 

19 

20class TorsionSimilarity(CollectiveVariable, openmm.CustomCompoundBondForce): 

21 r""" 

22 The degree of similarity between `n` pairs of torsion angles: 

23 

24 .. math:: 

25 

26 s({\bf r}) = \frac{n}{2} + \frac{1}{2} \sum_{i=1}^n \cos\Big( 

27 \phi^{\rm 1st}_i({\bf r}) - \phi^{\rm 2nd}_i({\bf r}) 

28 \Big) 

29 

30 where :math:`\phi^{\rm kth}_i` is the torsion angle at position :math:`i` in the 

31 :math:`k`-th list. 

32 

33 .. note:: 

34 

35 In `PLUMED <https://www.plumed.org/doc-v2.8/user-doc/html/_colvar.html>`_, this 

36 collective variable is called ``DIHCOR``. 

37 

38 Parameters 

39 ---------- 

40 firstList 

41 A list of :math:`n` tuples of four atom indices defining the first torsion 

42 angle in each pair. 

43 secondList 

44 A list of :math:`n` tuples of four atom indices defining the second torsion 

45 angle in each pair. 

46 pbc 

47 Whether to use periodic boundary conditions in distance calculations. 

48 name 

49 The name of the collective variable. 

50 

51 Raises 

52 ------ 

53 ValueError 

54 If firstList and secondList have different lengths 

55 

56 Example 

57 ------- 

58 >>> import cvpack 

59 >>> import mdtraj 

60 >>> import openmm 

61 >>> import openmm.unit as mmunit 

62 >>> from openmmtools import testsystems 

63 >>> model = testsystems.LysozymeImplicit() 

64 >>> traj = mdtraj.Trajectory( 

65 ... model.positions, mdtraj.Topology.from_openmm(model.topology) 

66 ... ) 

67 >>> phi_atoms, _ = mdtraj.compute_phi(traj) 

68 >>> valid_atoms = traj.top.select("resid 59 to 79 and backbone") 

69 >>> phi_atoms = [ 

70 ... phi 

71 ... for phi in phi_atoms 

72 ... if all(atom in valid_atoms for atom in phi) 

73 ... ] 

74 >>> torsion_similarity = cvpack.TorsionSimilarity( 

75 ... phi_atoms[1:], phi_atoms[:-1] 

76 ... ) 

77 >>> torsion_similarity.addToSystem(model.system) 

78 >>> integrator = openmm.VerletIntegrator(0) 

79 >>> platform = openmm.Platform.getPlatformByName("Reference") 

80 >>> context = openmm.Context(model.system, integrator, platform) 

81 >>> context.setPositions(model.positions) 

82 >>> torsion_similarity.getValue(context) 

83 18.659... dimensionless 

84 """ 

85 

86 def __init__( 

87 self, 

88 firstList: ArrayLike, 

89 secondList: ArrayLike, 

90 pbc: bool = False, 

91 name: str = "torsion_similarity", 

92 ) -> None: 

93 firstList = [list(map(int, first)) for first in firstList] 

94 secondList = [list(map(int, second)) for second in secondList] 

95 if len(firstList) != len(secondList): 

96 raise ValueError("firstList and secondList must have the same length") 

97 function = f"0.5*(1 + cos(min(delta, {2 * np.pi} - delta)))" 

98 definition = "delta = dihedral(p1, p2, p3, p4) - dihedral(p5, p6, p7, p8)" 

99 super().__init__(8, f"{function}; {definition}") 

100 for first, second in zip(firstList, secondList): 

101 self.addBond([*first, *second], []) 

102 self.setUsesPeriodicBoundaryConditions(pbc) 

103 self._registerCV( 

104 name, 

105 mmunit.dimensionless, 

106 firstList=firstList, 

107 secondList=secondList, 

108 pbc=pbc, 

109 ) 

110 

111 

112TorsionSimilarity.registerTag("!cvpack.TorsionSimilarity")