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.ForceFieldEnergy;
43 import ffx.potential.bonded.Torsion;
44 import ffx.potential.parameters.TorsionType;
45 import ffx.potential.terms.RestrainTorsionPotentialEnergy;
46
47 import java.util.logging.Level;
48 import java.util.logging.Logger;
49
50 import static edu.uiowa.jopenmm.OpenMMAmoebaLibrary.OpenMM_KJPerKcal;
51 import static edu.uiowa.jopenmm.OpenMMAmoebaLibrary.OpenMM_RadiansPerDegree;
52 import static java.lang.String.format;
53
54
55
56
57 public class RestrainTorsionsForce extends PeriodicTorsionForce {
58
59 private static final Logger logger = Logger.getLogger(RestrainTorsionsForce.class.getName());
60
61
62
63
64
65
66 public RestrainTorsionsForce(RestrainTorsionPotentialEnergy restrainTorsionPotentialEnergy) {
67 Torsion[] restrainTorsions = restrainTorsionPotentialEnergy.getRestrainTorsionArray();
68 for (Torsion restrainTorsion : restrainTorsions) {
69 int a1 = restrainTorsion.getAtom(0).getArrayIndex();
70 int a2 = restrainTorsion.getAtom(1).getArrayIndex();
71 int a3 = restrainTorsion.getAtom(2).getArrayIndex();
72 int a4 = restrainTorsion.getAtom(3).getArrayIndex();
73 TorsionType torsionType = restrainTorsion.torsionType;
74 int nTerms = torsionType.phase.length;
75 for (int j = 0; j < nTerms; j++) {
76 addTorsion(a1, a2, a3, a4, j + 1,
77 torsionType.phase[j] * OpenMM_RadiansPerDegree,
78 OpenMM_KJPerKcal * torsionType.torsionUnit * torsionType.amplitude[j]);
79 }
80 }
81 int forceGroup = restrainTorsionPotentialEnergy.getForceGroup();
82 setForceGroup(forceGroup);
83 logger.log(Level.INFO, format(" Restrain-Torsions \t%6d\t\t%1d", restrainTorsions.length, forceGroup));
84 }
85
86
87
88
89
90
91
92 public RestrainTorsionsForce(RestrainTorsionPotentialEnergy restrainTorsionPotentialEnergy, int topology, OpenMMDualTopologyEnergy openMMDualTopologyEnergy) {
93 ForceFieldEnergy forceFieldEnergy = openMMDualTopologyEnergy.getForceFieldEnergy(topology);
94 Torsion[] restrainTorsions = restrainTorsionPotentialEnergy.getRestrainTorsionArray();
95
96 double scale = openMMDualTopologyEnergy.getTopologyScale(topology);
97
98 for (Torsion restrainTorsion : restrainTorsions) {
99 int a1 = restrainTorsion.getAtom(0).getArrayIndex();
100 int a2 = restrainTorsion.getAtom(1).getArrayIndex();
101 int a3 = restrainTorsion.getAtom(2).getArrayIndex();
102 int a4 = restrainTorsion.getAtom(3).getArrayIndex();
103 a1 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a1);
104 a2 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a2);
105 a3 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a3);
106 a4 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a4);
107 TorsionType torsionType = restrainTorsion.torsionType;
108 int nTerms = torsionType.phase.length;
109 for (int j = 0; j < nTerms; j++) {
110 double k = restrainTorsion.getTorsionScale() * torsionType.torsionUnit * torsionType.amplitude[j];
111 k = k * (1 - scale);
112 addTorsion(a1, a2, a3, a4, j + 1,
113 torsionType.phase[j] * OpenMM_RadiansPerDegree,
114 OpenMM_KJPerKcal * k);
115 }
116 }
117
118 int forceGroup = forceFieldEnergy.getMolecularAssembly().getForceField().getInteger("RESTRAIN_TORSION_FORCE_GROUP", 0);
119 setForceGroup(forceGroup);
120 logger.info(format(" Restrain-Torsions: %10d", restrainTorsions.length));
121 logger.fine(format(" Force Group: %10d", forceGroup));
122 }
123
124
125
126
127
128
129
130 public static Force constructForce(OpenMMEnergy openMMEnergy) {
131 RestrainTorsionPotentialEnergy restrainTorsionPotentialEnergy =
132 openMMEnergy.getRestrainTorsionPotentialEnergy();
133 if (restrainTorsionPotentialEnergy == null) {
134 return null;
135 }
136 return new RestrainTorsionsForce(restrainTorsionPotentialEnergy);
137 }
138
139
140
141
142
143
144
145
146 public static Force constructForce(int topology, OpenMMDualTopologyEnergy openMMDualTopologyEnergy) {
147 ForceFieldEnergy forceFieldEnergy = openMMDualTopologyEnergy.getForceFieldEnergy(topology);
148 RestrainTorsionPotentialEnergy restrainTorsionPotentialEnergy = forceFieldEnergy.getRestrainTorsionPotentialEnergy();
149 if (restrainTorsionPotentialEnergy == null) {
150 return null;
151 }
152 return new RestrainTorsionsForce(restrainTorsionPotentialEnergy, topology, openMMDualTopologyEnergy);
153 }
154
155
156
157
158
159
160 public void updateForce(OpenMMEnergy openMMEnergy) {
161
162 RestrainTorsionPotentialEnergy restrainTorsionPotentialEnergy =
163 openMMEnergy.getRestrainTorsionPotentialEnergy();
164 if (restrainTorsionPotentialEnergy == null) {
165 return;
166 }
167 Torsion[] restrainTorsions = restrainTorsionPotentialEnergy.getRestrainTorsionArray();
168 int index = 0;
169 for (Torsion restrainTorsion : restrainTorsions) {
170 TorsionType torsionType = restrainTorsion.torsionType;
171 int nTerms = torsionType.phase.length;
172 int a1 = restrainTorsion.getAtom(0).getArrayIndex();
173 int a2 = restrainTorsion.getAtom(1).getArrayIndex();
174 int a3 = restrainTorsion.getAtom(2).getArrayIndex();
175 int a4 = restrainTorsion.getAtom(3).getArrayIndex();
176 for (int j = 0; j < nTerms; j++) {
177 double forceConstant = OpenMM_KJPerKcal * torsionType.torsionUnit * torsionType.amplitude[j];
178 setTorsionParameters(index++, a1, a2, a3, a4, j + 1,
179 torsionType.phase[j] * OpenMM_RadiansPerDegree, forceConstant);
180 }
181 }
182 updateParametersInContext(openMMEnergy.getContext());
183 }
184
185
186
187
188
189
190
191 public void updateForce(int topology, OpenMMDualTopologyEnergy openMMDualTopologyEnergy) {
192 ForceFieldEnergy forceFieldEnergy = openMMDualTopologyEnergy.getForceFieldEnergy(topology);
193 RestrainTorsionPotentialEnergy restrainTorsionPotentialEnergy = forceFieldEnergy.getRestrainTorsionPotentialEnergy();
194 if (restrainTorsionPotentialEnergy == null) {
195 return;
196 }
197
198 Torsion[] restrainTorsions = restrainTorsionPotentialEnergy.getRestrainTorsionArray();
199
200 double scale = openMMDualTopologyEnergy.getTopologyScale(topology);
201
202 int index = 0;
203 for (Torsion restrainTorsion : restrainTorsions) {
204 TorsionType torsionType = restrainTorsion.torsionType;
205 int nTerms = torsionType.phase.length;
206 int a1 = restrainTorsion.getAtom(0).getArrayIndex();
207 int a2 = restrainTorsion.getAtom(1).getArrayIndex();
208 int a3 = restrainTorsion.getAtom(2).getArrayIndex();
209 int a4 = restrainTorsion.getAtom(3).getArrayIndex();
210 a1 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a1);
211 a2 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a2);
212 a3 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a3);
213 a4 = openMMDualTopologyEnergy.mapToDualTopologyIndex(topology, a4);
214 for (int j = 0; j < nTerms; j++) {
215 double forceConstant = restrainTorsion.getTorsionScale() * torsionType.torsionUnit * torsionType.amplitude[j];
216 forceConstant = forceConstant * (1 - scale);
217 setTorsionParameters(index++, a1, a2, a3, a4, j + 1,
218 torsionType.phase[j] * OpenMM_RadiansPerDegree, OpenMM_KJPerKcal * forceConstant);
219 }
220 }
221 updateParametersInContext(openMMDualTopologyEnergy.getContext());
222 }
223
224 }