View Javadoc
1   // ******************************************************************************
2   //
3   // Title:       Force Field X.
4   // Description: Force Field X - Software for Molecular Biophysics.
5   // Copyright:   Copyright (c) Michael J. Schnieders 2001-2024.
6   //
7   // This file is part of Force Field X.
8   //
9   // Force Field X is free software; you can redistribute it and/or modify it
10  // under the terms of the GNU General Public License version 3 as published by
11  // the Free Software Foundation.
12  //
13  // Force Field X is distributed in the hope that it will be useful, but WITHOUT
14  // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
15  // FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
16  // details.
17  //
18  // You should have received a copy of the GNU General Public License along with
19  // Force Field X; if not, write to the Free Software Foundation, Inc., 59 Temple
20  // Place, Suite 330, Boston, MA 02111-1307 USA
21  //
22  // Linking this library statically or dynamically with other modules is making a
23  // combined work based on this library. Thus, the terms and conditions of the
24  // GNU General Public License cover the whole combination.
25  //
26  // As a special exception, the copyright holders of this library give you
27  // permission to link this library with independent modules to produce an
28  // executable, regardless of the license terms of these independent modules, and
29  // to copy and distribute the resulting executable under terms of your choice,
30  // provided that you also meet, for each linked independent module, the terms
31  // and conditions of the license of that module. An independent module is a
32  // module which is not derived from or based on this library. If you modify this
33  // library, you may extend this exception to your version of the library, but
34  // you are not obligated to do so. If you do not wish to do so, delete this
35  // exception statement from your version.
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   * OpenMM In-Plane Angle Force.
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     * Create an OpenMM Angle Force.
67     *
68     * @param openMMEnergy The OpenMM Energy instance that contains the angles.
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        // Clean up the Memory allocated by the OpenMMCustomCompoundBondForce constructor.
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          // Skip Normal angles unless this is ManyBody Titration.
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            // This is a place-holder Angle, in case the Normal Angle is switched to an
100           // In-Plane Angle during in the udpateInPlaneAngleForce.
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.log(Level.INFO, format("  In-Plane Angles \t%6d\t\t%1d", nAngles, forceGroup));
132     }
133   }
134 
135   /**
136    * Convenience method to construct an OpenMM In-Plane Angle Force.
137    *
138    * @param openMMEnergy The OpenMM Energy instance that contains the angles.
139    * @return An OpenMM Angle Force, or null if there are no angles.
140    */
141   public static Force constructForce(OpenMMEnergy openMMEnergy) {
142     Angle[] angles = openMMEnergy.getAngles();
143     if (angles == null || angles.length < 1) {
144       return null;
145     }
146     InPlaneAngleForce angleForce = new InPlaneAngleForce(openMMEnergy);
147     if (angleForce.nAngles > 0) {
148       return angleForce;
149     }
150     return null;
151   }
152 
153   /**
154    * Update an existing angle force for the OpenMM System.
155    *
156    * @param openMMEnergy The OpenMM Energy instance that contains the angles.
157    */
158   public void updateForce(OpenMMEnergy openMMEnergy) {
159     Angle[] angles = openMMEnergy.getAngles();
160     if (angles == null || angles.length < 1) {
161       return;
162     }
163     IntArray particles = new IntArray(0);
164     DoubleArray parameters = new DoubleArray(0);
165     int index = 0;
166     for (Angle angle : angles) {
167       AngleType.AngleMode angleMode = angle.angleType.angleMode;
168       if (!manyBodyTitration && angleMode == AngleType.AngleMode.NORMAL) {
169         // Skip Normal angles unless this is ManyBody Titration.
170       } else {
171         double theta0 = angle.angleType.angle[angle.nh];
172         double k = OpenMM_KJPerKcal * angle.angleType.angleUnit * angle.angleType.forceConstant;
173         int i1 = angle.getAtom(0).getXyzIndex() - 1;
174         int i2 = angle.getAtom(1).getXyzIndex() - 1;
175         int i3 = angle.getAtom(2).getXyzIndex() - 1;
176         // There is no 4th atom for normal angles, so set the index to first atom.
177         int i4 = 0;
178         if (angleMode == AngleType.AngleMode.NORMAL) {
179           // Zero the force constant for Normal Angles.
180           k = 0.0;
181           Atom fourthAtom = angle.getFourthAtomOfTrigonalCenter();
182           if (fourthAtom != null) {
183             i4 = fourthAtom.getXyzIndex() - 1;
184           } else {
185             while (i1 == i4 || i2 == i4 || i3 == i4) {
186               i4++;
187             }
188           }
189         } else {
190           i4 = angle.getAtom4().getXyzIndex() - 1;
191         }
192         particles.append(i1);
193         particles.append(i2);
194         particles.append(i3);
195         particles.append(i4);
196         parameters.append(theta0);
197         parameters.append(k);
198         setBondParameters(index++, particles, parameters);
199         particles.resize(0);
200         parameters.resize(0);
201       }
202     }
203     particles.destroy();
204     parameters.destroy();
205     updateParametersInContext(openMMEnergy.getContext());
206   }
207 }