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.DoubleArray;
41 import ffx.openmm.Force;
42 import ffx.openmm.IntArray;
43 import ffx.openmm.CustomCompoundBondForce;
44 import ffx.potential.bonded.Angle;
45 import ffx.potential.bonded.Atom;
46 import ffx.potential.parameters.AngleType;
47 import ffx.potential.parameters.ForceField;
48
49 import java.util.logging.Level;
50 import java.util.logging.Logger;
51
52 import static edu.uiowa.jopenmm.OpenMMAmoebaLibrary.OpenMM_KJPerKcal;
53 import static java.lang.String.format;
54
55
56
57
58 public class InPlaneAngleForce extends CustomCompoundBondForce {
59
60 private static final Logger logger = Logger.getLogger(InPlaneAngleForce.class.getName());
61
62 private int nAngles = 0;
63 private final boolean manyBodyTitration;
64
65
66
67
68
69
70 public InPlaneAngleForce(OpenMMEnergy openMMEnergy) {
71 super(4, openMMEnergy.getInPlaneAngleEnergyString());
72 ForceField forceField = openMMEnergy.getMolecularAssembly().getForceField();
73 manyBodyTitration = forceField.getBoolean("MANYBODY_TITRATION", false);
74 Angle[] angles = openMMEnergy.getAngles();
75 if (angles == null || angles.length < 1) {
76
77 destroy();
78 return;
79 }
80 addPerBondParameter("theta0");
81 addPerBondParameter("k");
82 setName("InPlaneAngle");
83
84 IntArray particles = new IntArray(0);
85 DoubleArray parameters = new DoubleArray(0);
86 for (Angle angle : angles) {
87 AngleType.AngleMode angleMode = angle.angleType.angleMode;
88
89 if (!manyBodyTitration && angleMode == AngleType.AngleMode.NORMAL) {
90
91 } else {
92 double theta0 = angle.angleType.angle[angle.nh];
93 double k = OpenMM_KJPerKcal * angle.angleType.angleUnit * angle.angleType.forceConstant;
94 int i1 = angle.getAtom(0).getXyzIndex() - 1;
95 int i2 = angle.getAtom(1).getXyzIndex() - 1;
96 int i3 = angle.getAtom(2).getXyzIndex() - 1;
97 int i4 = 0;
98 if (angleMode == AngleType.AngleMode.NORMAL) {
99
100
101 k = 0.0;
102 Atom fourthAtom = angle.getFourthAtomOfTrigonalCenter();
103 if (fourthAtom != null) {
104 i4 = fourthAtom.getXyzIndex() - 1;
105 } else {
106 while (i1 == i4 || i2 == i4 || i3 == i4) {
107 i4++;
108 }
109 }
110 } else {
111 i4 = angle.getAtom4().getXyzIndex() - 1;
112 }
113 particles.append(i1);
114 particles.append(i2);
115 particles.append(i3);
116 particles.append(i4);
117 parameters.append(theta0);
118 parameters.append(k);
119 addBond(particles, parameters);
120 nAngles++;
121 particles.resize(0);
122 parameters.resize(0);
123 }
124 }
125 particles.destroy();
126 parameters.destroy();
127
128 if (nAngles > 0) {
129 int forceGroup = forceField.getInteger("IN_PLANE_ANGLE_FORCE_GROUP", 0);
130 setForceGroup(forceGroup);
131 logger.info(format(" In-Plane Angles: %10d", nAngles));
132 logger.fine(format(" Force Group: %10d", forceGroup));
133 }
134 }
135
136
137
138
139
140
141
142 public static Force constructForce(OpenMMEnergy openMMEnergy) {
143 Angle[] angles = openMMEnergy.getAngles();
144 if (angles == null || angles.length < 1) {
145 return null;
146 }
147 InPlaneAngleForce angleForce = new InPlaneAngleForce(openMMEnergy);
148 if (angleForce.nAngles > 0) {
149 return angleForce;
150 }
151 return null;
152 }
153
154
155
156
157
158
159 public void updateForce(OpenMMEnergy openMMEnergy) {
160 Angle[] angles = openMMEnergy.getAngles();
161 if (angles == null || angles.length < 1) {
162 return;
163 }
164 IntArray particles = new IntArray(0);
165 DoubleArray parameters = new DoubleArray(0);
166 int index = 0;
167 for (Angle angle : angles) {
168 AngleType.AngleMode angleMode = angle.angleType.angleMode;
169 if (!manyBodyTitration && angleMode == AngleType.AngleMode.NORMAL) {
170
171 } else {
172 double theta0 = angle.angleType.angle[angle.nh];
173 double k = OpenMM_KJPerKcal * angle.angleType.angleUnit * angle.angleType.forceConstant;
174 int i1 = angle.getAtom(0).getXyzIndex() - 1;
175 int i2 = angle.getAtom(1).getXyzIndex() - 1;
176 int i3 = angle.getAtom(2).getXyzIndex() - 1;
177
178 int i4 = 0;
179 if (angleMode == AngleType.AngleMode.NORMAL) {
180
181 k = 0.0;
182 Atom fourthAtom = angle.getFourthAtomOfTrigonalCenter();
183 if (fourthAtom != null) {
184 i4 = fourthAtom.getXyzIndex() - 1;
185 } else {
186 while (i1 == i4 || i2 == i4 || i3 == i4) {
187 i4++;
188 }
189 }
190 } else {
191 i4 = angle.getAtom4().getXyzIndex() - 1;
192 }
193 particles.append(i1);
194 particles.append(i2);
195 particles.append(i3);
196 particles.append(i4);
197 parameters.append(theta0);
198 parameters.append(k);
199 setBondParameters(index++, particles, parameters);
200 particles.resize(0);
201 parameters.resize(0);
202 }
203 }
204 particles.destroy();
205 parameters.destroy();
206 updateParametersInContext(openMMEnergy.getContext());
207 }
208 }