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-2025.
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.numerics.multipole;
39  
40  /**
41   * The PolarizableMultipole class defines a polarizable multipole.
42   *
43   * @author Michael J. Schnieders
44   * @since 1.0
45   */
46  public class PolarizableMultipole {
47  
48    private static final double oneThird = 1.0 / 3.0;
49    private static final double twoThirds = 2.0 / 3.0;
50  
51    /**
52     * Partial charge.
53     */
54    protected double q;
55    /**
56     * Dipole x-component.
57     */
58    protected double dx;
59    /**
60     * Dipole y-component.
61     */
62    protected double dy;
63    /**
64     * Dipole z-component.
65     */
66    protected double dz;
67    /**
68     * Quadrupole xx-component multiplied by 1/3.
69     */
70    protected double qxx;
71    /**
72     * Quadrupole yy-component multiplied by 1/3.
73     */
74    protected double qyy;
75    /**
76     * Quadrupole zz-component multiplied by 1/3.
77     */
78    protected double qzz;
79    /**
80     * Quadrupole xy-component multiplied by 2/3.
81     */
82    protected double qxy;
83    /**
84     * Quadrupole xz-component multiplied by 2/3.
85     */
86    protected double qxz;
87    /**
88     * Quadrupole xz-component multiplied by 2/3.
89     */
90    protected double qyz;
91  
92    /**
93     * Induced dipole x-component.
94     */
95    protected double ux;
96    /**
97     * Induced dipole y-component.
98     */
99    protected double uy;
100   /**
101    * Induced dipole z-component.
102    */
103   protected double uz;
104   /**
105    * Induced dipole chain rule x-component.
106    */
107   protected double px;
108   /**
109    * Induced dipole chain rule y-component.
110    */
111   protected double py;
112   /**
113    * Induced dipole chain rule z-component.
114    */
115   protected double pz;
116   /**
117    * Averaged induced dipole + induced dipole chain-rule x-component: sx = 0.5 * (ux + px).
118    */
119   protected double sx;
120   /**
121    * Averaged induced dipole + induced dipole chain-rule y-component: sy = 0.5 * (uy + py).
122    */
123   protected double sy;
124   /**
125    * Averaged induced dipole + induced dipole chain-rule z-component: sz = 0.5 * (uz + pz).
126    */
127   protected double sz;
128   protected double Z;
129 
130 
131   /**
132    * PolarizableMultipole constructor with zero moments.
133    */
134   public PolarizableMultipole() {
135   }
136 
137   /**
138    * PolarizableMultipole constructor.
139    *
140    * @param Q   Multipole Q[q, dx, dy, dz, qxx, qyy, qzz, qxy, qxz, qyz]
141    * @param u   Induced dipole u[ux, uy, uz]
142    * @param uCR Induced dipole chain-rule uCR[ux, uy, uz]
143    */
144   public PolarizableMultipole(double[] Q, double[] u, double[] uCR) {
145     setPermanentMultipole(Q);
146     setInducedDipole(u, uCR);
147   }
148 
149   /**
150    * PolarizableMultipole constructor.
151    *
152    * @param Q   Multipole Q[q, dx, dy, dz, qxx, qyy, qzz, qxy, qxz, qyz]
153    * @param u   Induced dipole u[ux, uy, uz]
154    * @param uCR Induced dipole chain-rule uCR[ux, uy, uz]
155    */
156   public PolarizableMultipole(double[] Q, double[] u, double[] uCR, double Z) {
157     setPermanentMultipole(Q);
158     setInducedDipole(u, uCR);
159     this.Z = Z;
160     this.q = -Z + q;
161   }
162 
163   /**
164    * Set the permanent multipole.
165    * <p>
166    * Note that the quadrupole trace components are multiplied by 1/3 and the
167    * off-diagonal components are multiplied by 2/3.
168    *
169    * @param Q   Multipole Q[q, dx, dy, dz, qxx, qyy, qzz, qxy, qxz, qyz]
170    * @param u   Induced dipole u[ux, uy, uz]
171    * @param uCR Induced dipole chain-rule uCR[ux, uy, uz]
172    */
173   public void set(double[] Q, double[] u, double[] uCR) {
174     setPermanentMultipole(Q);
175     setInducedDipole(u, uCR);
176   }
177 
178   /**
179    * Set the permanent multipole.
180    * <p>
181    * Note that the quadrupole trace components are multiplied by 1/3 and the
182    * off-diagonal components are multiplied by 2/3.
183    *
184    * @param Q Multipole Q[q, dx, dy, dz, qxx, qyy, qzz, qxy, qxz, qyz]
185    */
186   public final void setPermanentMultipole(double[] Q) {
187     q = Q[0];
188     dx =  Q[1];
189     dy =  Q[2];
190     dz =  Q[3];
191     qxx = Q[4] * oneThird;
192     qyy = Q[5] * oneThird;
193     qzz = Q[6] * oneThird;
194     qxy = Q[7] * twoThirds;
195     qxz = Q[8] * twoThirds;
196     qyz = Q[9] * twoThirds;
197   }
198 
199   /**
200    * Set the induced dipole.
201    *
202    * @param u   Induced dipole u[ux, uy, uz]
203    * @param uCR Induced dipole chain-rule uCR[ux, uy, uz]
204    */
205   public final void setInducedDipole(double[] u, double[] uCR) {
206     ux = u[0];
207     uy = u[1];
208     uz = u[2];
209     px = uCR[0];
210     py = uCR[1];
211     pz = uCR[2];
212     sx = 0.5 * (ux + px);
213     sy = 0.5 * (uy + py);
214     sz = 0.5 * (uz + pz);
215   }
216 
217   /**
218    * Clear the induced dipoles.
219    */
220   public void clearInducedDipoles() {
221     ux = 0.0;
222     uy = 0.0;
223     uz = 0.0;
224     px = 0.0;
225     py = 0.0;
226     pz = 0.0;
227     sx = 0.0;
228     sy = 0.0;
229     sz = 0.0;
230   }
231 
232   /**
233    * Compute the scaled and averaged induced dipole.
234    *
235    * @param scaleInduction Induction mask scale factor.
236    * @param scaleEnergy    Energy mask scale factor.
237    */
238   public final void applyMasks(double scaleInduction, double scaleEnergy) {
239     // [Ux, Uy, Uz] resulted from induction masking rules, and we now apply the energy mask.
240     // [Px, Py, Pz] resulted from energy masking rules, and we now apply the induction mask.
241     sx = 0.5 * (ux * scaleEnergy + px * scaleInduction);
242     sy = 0.5 * (uy * scaleEnergy + py * scaleInduction);
243     sz = 0.5 * (uz * scaleEnergy + pz * scaleInduction);
244   }
245 
246 }