1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38 package ffx.potential.openmm;
39
40 import ffx.openmm.Force;
41 import ffx.openmm.PeriodicTorsionForce;
42 import ffx.potential.bonded.Torsion;
43 import ffx.potential.parameters.ForceField;
44 import ffx.potential.parameters.TorsionType;
45
46 import java.util.logging.Level;
47 import java.util.logging.Logger;
48
49 import static edu.uiowa.jopenmm.OpenMMAmoebaLibrary.OpenMM_KJPerKcal;
50 import static edu.uiowa.jopenmm.OpenMMAmoebaLibrary.OpenMM_RadiansPerDegree;
51 import static java.lang.String.format;
52
53 public class TorsionForce extends PeriodicTorsionForce {
54
55 private static final Logger logger = Logger.getLogger(TorsionForce.class.getName());
56
57 private final boolean manyBodyTitration;
58
59 private double lambdaTorsion = 1.0;
60
61 public TorsionForce(OpenMMEnergy openMMEnergy) {
62 ForceField forceField = openMMEnergy.getMolecularAssembly().getForceField();
63 manyBodyTitration = forceField.getBoolean("MANYBODY_TITRATION", false);
64 Torsion[] torsions = openMMEnergy.getTorsions();
65 if (torsions == null || torsions.length < 1) {
66 return;
67 }
68
69 for (Torsion torsion : torsions) {
70 int a1 = torsion.getAtom(0).getXyzIndex() - 1;
71 int a2 = torsion.getAtom(1).getXyzIndex() - 1;
72 int a3 = torsion.getAtom(2).getXyzIndex() - 1;
73 int a4 = torsion.getAtom(3).getXyzIndex() - 1;
74 TorsionType torsionType = torsion.torsionType;
75 int nTerms = torsionType.phase.length;
76 for (int j = 0; j < nTerms; j++) {
77 double k = torsion.getTorsionScale() * torsionType.torsionUnit * torsionType.amplitude[j];
78 addTorsion(a1, a2, a3, a4, j + 1,
79 torsionType.phase[j] * OpenMM_RadiansPerDegree, OpenMM_KJPerKcal * k);
80 }
81
82
83 if (manyBodyTitration) {
84 for (int j = nTerms; j < 6; j++) {
85 addTorsion(a1, a2, a3, a4, j + 1, 0.0, 0.0);
86 }
87 }
88 }
89
90 int forceGroup = forceField.getInteger("TORSION_FORCE_GROUP", 0);
91 setForceGroup(forceGroup);
92 logger.info(format(" Torsions: %10d", torsions.length));
93 logger.fine(format(" Force Group: %10d", forceGroup));
94 }
95
96
97
98
99
100
101 public void setLambdaTorsion(double lambdaTorsion) {
102 this.lambdaTorsion = lambdaTorsion;
103 }
104
105
106
107
108
109
110
111 public static Force constructForce(OpenMMEnergy openMMEnergy) {
112 Torsion[] torsions = openMMEnergy.getTorsions();
113 if (torsions == null || torsions.length < 1) {
114 return null;
115 }
116 return new TorsionForce(openMMEnergy);
117 }
118
119
120
121
122
123
124 public void updateForce(OpenMMEnergy openMMEnergy) {
125
126 Torsion[] torsions = openMMEnergy.getTorsions();
127 if (torsions == null || torsions.length < 1) {
128 return;
129 }
130
131 int index = 0;
132 for (Torsion torsion : torsions) {
133 TorsionType torsionType = torsion.torsionType;
134 int nTerms = torsionType.phase.length;
135 int a1 = torsion.getAtom(0).getXyzIndex() - 1;
136 int a2 = torsion.getAtom(1).getXyzIndex() - 1;
137 int a3 = torsion.getAtom(2).getXyzIndex() - 1;
138 int a4 = torsion.getAtom(3).getXyzIndex() - 1;
139 for (int j = 0; j < nTerms; j++) {
140 double k = torsion.getTorsionScale() * torsionType.torsionUnit * torsionType.amplitude[j] * lambdaTorsion;
141 setTorsionParameters(index++, a1, a2, a3, a4, j + 1,
142 torsionType.phase[j] * OpenMM_RadiansPerDegree, OpenMM_KJPerKcal * k);
143 }
144
145
146 if (manyBodyTitration) {
147 for (int j = nTerms; j < 6; j++) {
148 setTorsionParameters(index++, a1, a2, a3, a4, j + 1, 0.0, 0.0);
149 }
150 }
151 }
152 updateParametersInContext(openMMEnergy.getContext());
153 }
154
155 }