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.parameters;
39
40 import ffx.utilities.FFXProperty;
41 import org.w3c.dom.Document;
42 import org.w3c.dom.Element;
43
44 import java.util.Arrays;
45 import java.util.Comparator;
46 import java.util.HashMap;
47 import java.util.Map;
48 import java.util.logging.Level;
49 import java.util.logging.Logger;
50
51 import static ffx.potential.parameters.ForceField.ForceFieldType.OPBEND;
52 import static ffx.utilities.Constants.DEGREES_PER_RADIAN;
53 import static ffx.utilities.Constants.KCAL_TO_KJ;
54 import static ffx.utilities.PropertyGroup.EnergyUnitConversion;
55 import static ffx.utilities.PropertyGroup.LocalGeometryFunctionalForm;
56 import static ffx.utilities.PropertyGroup.PotentialFunctionParameter;
57 import static java.lang.Double.parseDouble;
58 import static java.lang.Integer.parseInt;
59 import static java.lang.String.format;
60 import static java.lang.String.valueOf;
61 import static java.util.Arrays.copyOf;
62 import static org.apache.commons.math3.util.FastMath.PI;
63 import static org.apache.commons.math3.util.FastMath.pow;
64
65
66
67
68
69
70
71 @FFXProperty(name = "opbend", clazz = String.class, propertyGroup = PotentialFunctionParameter, description = """
72 [4 integers and 1 real]
73 Provides the values for a single out-of-plane bending potential parameter.
74 The first integer modifier is the atom class of the out-of-plane atom and the second integer is the atom class of the central trigonal atom.
75 The third and fourth integers give the atom classes of the two remaining atoms attached to the trigonal atom.
76 Values of zero for the third and fourth integers are treated as wildcards, and can represent any atom type.
77 The real number modifier gives the force constant value for the out-of-plane angle.
78 The default units for the force constant are kcal/mole/radian^2, but this can be controlled via the opbendunit keyword.
79 """)
80 public final class OutOfPlaneBendType extends BaseType implements Comparator<String> {
81
82
83
84
85 public static final double DEFAULT_OPBEND_CUBIC = 0.0;
86
87
88
89 public static final double DEFAULT_OPBEND_QUARTIC = 0.0;
90
91
92
93 public static final double DEFAULT_OPBEND_PENTIC = 0.0;
94
95
96
97 public static final double DEFAULT_OPBEND_SEXTIC = 0.0;
98
99
100
101
102 @FFXProperty(name = "opbend-cubic", propertyGroup = LocalGeometryFunctionalForm, defaultValue = "0.0", description = """
103 Sets the value of the cubic term in the Taylor series expansion form of the out-of-plane bending potential energy.
104 The real number modifier gives the value of the coefficient as a multiple of the quadratic coefficient.
105 This term multiplied by the out-of-plane bending energy unit conversion factor, the force constant,
106 and the cube of the deviation of the out-of-plane angle from zero gives the cubic contribution to the out-of-plane bending energy.
107 The default value in the absence of the opbend-cubic keyword is zero; i.e., the cubic out-of-plane bending term is omitted.
108 """)
109 public double cubic = DEFAULT_OPBEND_CUBIC;
110
111
112
113
114 @FFXProperty(name = "opbend-quartic", propertyGroup = LocalGeometryFunctionalForm, defaultValue = "0.0", description = """
115 Sets the value of the quartic term in the Taylor series expansion form of the out-of-plane bending potential energy.
116 The real number modifier gives the value of the coefficient as a multiple of the quadratic coefficient.
117 This term multiplied by the out-of-plane bending energy unit conversion factor, the force constant,
118 and the forth power of the deviation of the out-of-plane angle from zero gives the quartic contribution to the out-of-plane bending energy.
119 The default value in the absence of the opbend-quartic keyword is zero; i.e., the quartic out-of-plane bending term is omitted.
120 """)
121 public double quartic = DEFAULT_OPBEND_QUARTIC;
122
123
124
125
126 @FFXProperty(name = "opbend-pentic", propertyGroup = LocalGeometryFunctionalForm, defaultValue = "0.0", description = """
127 Sets the value of the fifth power term in the Taylor series expansion form of the out-of-plane bending potential energy.
128 The real number modifier gives the value of the coefficient as a multiple of the quadratic coefficient.
129 This term multiplied by the out-of-plane bending energy unit conversion factor, the force constant,
130 and the fifth power of the deviation of the out-of-plane angle from zero gives the pentic contribution to the out-of-plane bending energy.
131 The default value in the absence of the opbend-pentic keyword is zero; i.e., the pentic out-of-plane bending term is omitted.
132 """)
133 public double pentic = DEFAULT_OPBEND_PENTIC;
134
135
136
137
138 @FFXProperty(name = "opbend-sextic", propertyGroup = LocalGeometryFunctionalForm, defaultValue = "0.0", description = """
139 Sets the value of the sixth power term in the Taylor series expansion form of the out-of-plane bending potential energy.
140 The real number modifier gives the value of the coefficient as a multiple of the quadratic coefficient.
141 This term multiplied by the out-of-plane bending energy unit conversion factor, the force constant,
142 and the sixth power of the deviation of the out-of-plane angle from zero gives the sextic contribution to the out-of-plane bending energy.
143 The default value in the absence of the opbend-sextic keyword is zero; i.e., the sextic out-of-plane bending term is omitted.
144 """)
145 public double sextic = DEFAULT_OPBEND_SEXTIC;
146
147
148
149
150
151
152
153
154
155
156 @FFXProperty(name = "opbendunit", propertyGroup = EnergyUnitConversion, defaultValue = "(Pi/180)^2", description = """
157 Sets the scale factor needed to convert the energy value computed by the out-of-plane bending potential into units of kcal/mole. "
158 The correct value is force field dependent and typically provided in the header of the master force field parameter file.
159 """)
160 public double opBendUnit = DEFAULT_OPBEND_UNIT;
161
162 public static final double DEFAULT_OPBEND_UNIT = pow(PI / 180.0, 2);
163
164
165
166 private static final Logger logger = Logger.getLogger(OutOfPlaneBendType.class.getName());
167
168
169
170 public final int[] atomClasses;
171
172
173
174 public final double forceConstant;
175
176
177
178
179
180
181
182 public OutOfPlaneBendType(int[] atomClasses, double forceConstant) {
183 super(OPBEND, sortKey(atomClasses));
184 this.atomClasses = atomClasses;
185 this.forceConstant = forceConstant;
186 }
187
188
189
190
191
192
193
194
195
196
197 public static OutOfPlaneBendType average(OutOfPlaneBendType outOfPlaneBendType1,
198 OutOfPlaneBendType outOfPlaneBendType2, int[] atomClasses) {
199 if (outOfPlaneBendType1 == null || outOfPlaneBendType2 == null || atomClasses == null) {
200 return null;
201 }
202
203 double forceConstant =
204 (outOfPlaneBendType1.forceConstant + outOfPlaneBendType2.forceConstant) / 2.0;
205
206 return new OutOfPlaneBendType(atomClasses, forceConstant);
207 }
208
209
210
211
212
213
214
215
216 public static OutOfPlaneBendType parse(String input, String[] tokens) {
217 if (tokens.length < 6) {
218 logger.log(Level.WARNING, "Invalid OPBEND type:\n{0}", input);
219 } else {
220 try {
221 int[] atomClasses = new int[4];
222 atomClasses[0] = parseInt(tokens[1]);
223 atomClasses[1] = parseInt(tokens[2]);
224 atomClasses[2] = parseInt(tokens[3]);
225 atomClasses[3] = parseInt(tokens[4]);
226 double forceConstant = parseDouble(tokens[5]);
227 return new OutOfPlaneBendType(atomClasses, forceConstant);
228 } catch (NumberFormatException e) {
229 String message = "Exception parsing OPBEND type:\n" + input + "\n";
230 logger.log(Level.SEVERE, message, e);
231 }
232 }
233 return null;
234 }
235
236
237
238
239
240
241
242 public static String sortKey(int[] c) {
243 if (c == null || c.length != 4) {
244 return null;
245 }
246 return c[0] + " " + c[1] + " " + c[2] + " " + c[3];
247 }
248
249
250
251
252 @Override
253 public int compare(String s1, String s2) {
254 String[] keys1 = s1.split(" ");
255 String[] keys2 = s2.split(" ");
256
257 for (int i = 0; i < 4; i++) {
258 int c1 = parseInt(keys1[i]);
259 int c2 = parseInt(keys2[i]);
260 if (c1 < c2) {
261 return -1;
262 } else if (c1 > c2) {
263 return 1;
264 }
265 }
266 return 0;
267 }
268
269
270
271
272 @Override
273 public boolean equals(Object o) {
274 if (this == o) {
275 return true;
276 }
277 if (o == null || getClass() != o.getClass()) {
278 return false;
279 }
280 OutOfPlaneBendType outOfPlaneBendType = (OutOfPlaneBendType) o;
281 return Arrays.equals(atomClasses, outOfPlaneBendType.atomClasses);
282 }
283
284
285
286
287 @Override
288 public int hashCode() {
289 return Arrays.hashCode(atomClasses);
290 }
291
292
293
294
295
296
297 public void incrementClasses(int increment) {
298 for (int i = 0; i < atomClasses.length; i++) {
299 if (atomClasses[i] > 0) {
300 atomClasses[i] += increment;
301 }
302 }
303 setKey(sortKey(atomClasses));
304 }
305
306
307
308
309
310
311
312 public OutOfPlaneBendType patchClasses(HashMap<AtomType, AtomType> typeMap) {
313 int count = 0;
314 int len = atomClasses.length;
315
316
317 for (AtomType newType : typeMap.keySet()) {
318
319 for (int atomClass : atomClasses) {
320 if (atomClass == newType.atomClass) {
321 count++;
322 }
323 }
324 }
325
326
327 if (count == 1) {
328 int[] newClasses = copyOf(atomClasses, len);
329 for (AtomType newType : typeMap.keySet()) {
330 for (int i = 0; i < len; i++) {
331 if (atomClasses[i] == newType.atomClass) {
332 AtomType knownType = typeMap.get(newType);
333 newClasses[i] = knownType.atomClass;
334 }
335 }
336 }
337 return new OutOfPlaneBendType(newClasses, forceConstant);
338 }
339 return null;
340 }
341
342
343
344
345
346
347 @Override
348 public String toString() {
349 return String.format("opbend %5d %5d %5d %5d %6.2f", atomClasses[0], atomClasses[1],
350 atomClasses[2], atomClasses[3], forceConstant);
351 }
352
353
354
355
356
357
358
359 public static Element getXMLForce(Document doc, ForceField forceField) {
360 Map<String, OutOfPlaneBendType> types = forceField.getOutOfPlaneBendTypes();
361 if (!types.values().isEmpty()) {
362 Element node = doc.createElement("AmoebaOutOfPlaneBendForce");
363 node.setAttribute("type", forceField.getString("opbendtype", "ALLINGER"));
364 node.setAttribute("opbend-cubic", valueOf(forceField.getDouble("opbend-cubic", DEFAULT_OPBEND_CUBIC)));
365 node.setAttribute("opbend-quartic", valueOf(forceField.getDouble("opbend-quartic", DEFAULT_OPBEND_QUARTIC)));
366 node.setAttribute("opbend-pentic", valueOf(forceField.getDouble("opbend-pentic", DEFAULT_OPBEND_PENTIC)));
367 node.setAttribute("opbend-sextic", valueOf(forceField.getDouble("opbend-sextic", DEFAULT_OPBEND_SEXTIC)));
368 for (OutOfPlaneBendType outOfPlaneBendType : types.values()) {
369 node.appendChild(outOfPlaneBendType.toXML(doc));
370 }
371 return node;
372 }
373 return null;
374 }
375
376
377
378
379 public Element toXML(Document doc) {
380 Element node = doc.createElement("Angle");
381 int i = 1;
382 for (int ac : atomClasses) {
383 if (ac == 0) {
384 node.setAttribute(format("class%d", i), "");
385 } else {
386 node.setAttribute(format("class%d", i), format("%d", ac));
387 }
388 i++;
389 }
390
391 node.setAttribute("k", format("%f", forceConstant * KCAL_TO_KJ / (DEGREES_PER_RADIAN * DEGREES_PER_RADIAN)));
392 return node;
393 }
394 }