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 import jdk.incubator.vector.DoubleVector;
41
42 import static ffx.numerics.multipole.GKMultipoleOrder.DIPOLE;
43 import static ffx.numerics.multipole.GKMultipoleOrder.MONOPOLE;
44 import static ffx.numerics.multipole.GKMultipoleOrder.QUADRUPOLE;
45 import static ffx.numerics.multipole.GKTensorMode.BORN;
46 import static ffx.numerics.multipole.GKTensorMode.POTENTIAL;
47
48 /**
49 * The GKEnergyQI class computes the Generalized Kirkwood energy and forces using a QI frame.
50 */
51 public class GKEnergyQISIMD {
52
53 private final GKSourceSIMD gkSource;
54 private final GKTensorQISIMD gkMonopole;
55 private final GKTensorQISIMD gkDipole;
56 private final GKTensorQISIMD gkQuadrupole;
57
58 private final DoubleVector one = DoubleVector.zero(DoubleVector.SPECIES_PREFERRED).add(1.0);
59
60 /**
61 * Compute the GK Energy using a QI frame.
62 *
63 * @param soluteDielectric Solute dielectric constant.
64 * @param solventDielectric Solvent dielectric constant.
65 * @param gkc The GK interaction parameter.
66 * @param gradient If true, the gradient will be computed.
67 */
68 public GKEnergyQISIMD(double soluteDielectric, double solventDielectric, double gkc, boolean gradient) {
69 int monopoleOrder = 2;
70 int dipoleOrder = 3;
71 int quadrupoleOrder = 4;
72 if (gradient) {
73 monopoleOrder = 3;
74 dipoleOrder = 4;
75 quadrupoleOrder = 5;
76 }
77 gkSource = new GKSourceSIMD(quadrupoleOrder, gkc);
78 gkMonopole = new GKTensorQISIMD(MONOPOLE, monopoleOrder, gkSource, soluteDielectric, solventDielectric);
79 gkDipole = new GKTensorQISIMD(DIPOLE, dipoleOrder, gkSource, soluteDielectric, solventDielectric);
80 gkQuadrupole = new GKTensorQISIMD(QUADRUPOLE, quadrupoleOrder, gkSource, soluteDielectric, solventDielectric);
81 }
82
83 /**
84 * Initialize the potential.
85 *
86 * @param r The separation. vector.
87 * @param r2 The squared separation.
88 * @param rbi The Born radius of atom i.
89 * @param rbk The Born radius of atom k.
90 */
91 public void initPotential(DoubleVector[] r, DoubleVector r2, DoubleVector rbi, DoubleVector rbk) {
92 gkSource.generateSource(POTENTIAL, QUADRUPOLE, r2, rbi, rbk);
93 gkMonopole.setR(r);
94 gkDipole.setR(r);
95 gkQuadrupole.setR(r);
96 gkMonopole.generateTensor();
97 gkDipole.generateTensor();
98 gkQuadrupole.generateTensor();
99 }
100
101 /**
102 * Initialize for computing Born chain-rule terms.
103 *
104 * @param r The separation vector.
105 * @param r2 The squared separation.
106 * @param rbi The Born radius of atom i.
107 * @param rbk The Born radius of atom k.
108 */
109 public void initBorn(DoubleVector[] r, DoubleVector r2, DoubleVector rbi, DoubleVector rbk) {
110 gkSource.generateSource(BORN, QUADRUPOLE, r2, rbi, rbk);
111 gkMonopole.setR(r);
112 gkDipole.setR(r);
113 gkQuadrupole.setR(r);
114 gkMonopole.generateTensor();
115 gkDipole.generateTensor();
116 gkQuadrupole.generateTensor();
117 }
118
119 /**
120 * Compute the multipole energy.
121 *
122 * @param mI The polarizable multipole of atom i.
123 * @param mK The polarizable multipole of atom k.
124 * @return The multipole energy.
125 */
126 public DoubleVector multipoleEnergy(PolarizableMultipoleSIMD mI, PolarizableMultipoleSIMD mK) {
127 DoubleVector em = gkMonopole.multipoleEnergy(mI, mK);
128 DoubleVector ed = gkDipole.multipoleEnergy(mI, mK);
129 DoubleVector eq = gkQuadrupole.multipoleEnergy(mI, mK);
130 return em.add(ed).add(eq);
131 }
132
133 /**
134 * Compute the polarization energy.
135 *
136 * @param mI The polarizable multipole of atom i.
137 * @param mK The polarizable multipole of atom k.
138 * @return The polarization energy.
139 */
140 public DoubleVector polarizationEnergy(PolarizableMultipoleSIMD mI, PolarizableMultipoleSIMD mK) {
141 DoubleVector emp = gkMonopole.polarizationEnergy(mI, mK);
142 DoubleVector edp = gkDipole.polarizationEnergy(mI, mK);
143 DoubleVector eqp = gkQuadrupole.polarizationEnergy(mI, mK);
144 return emp.add(edp).add(eqp);
145 }
146
147 /**
148 * Compute the multipole energy and gradient.
149 *
150 * @param mI The polarizable multipole of atom i.
151 * @param mK The polarizable multipole of atom k.
152 * @param gI The gradient for atom i.
153 * @param tI The torque on atom i.
154 * @param tK The torque on atom k.
155 * @return The multipole energy.
156 */
157 public DoubleVector multipoleEnergyAndGradient(PolarizableMultipoleSIMD mI, PolarizableMultipoleSIMD mK,
158 DoubleVector[] gI, DoubleVector[] tI, DoubleVector[] tK) {
159 DoubleVector[] gK = new DoubleVector[3];
160 DoubleVector em = gkMonopole.multipoleEnergyAndGradient(mI, mK, gI, gK, tI, tK);
161 DoubleVector ed = gkDipole.multipoleEnergyAndGradient(mI, mK, gI, gK, tI, tK);
162 DoubleVector eq = gkQuadrupole.multipoleEnergyAndGradient(mI, mK, gI, gK, tI, tK);
163 return em.add(ed).add(eq);
164 }
165
166 /**
167 * Compute the polarization energy and gradient.
168 *
169 * @param mI The polarizable multipole of atom i.
170 * @param mK The polarizable multipole of atom k.
171 * @param mutualMask The mutual polarization mask.
172 * @param gI The gradient for atom i.
173 * @param tI The torque on atom i.
174 * @param tK The torque on atom k.
175 * @return The polarization energy.
176 */
177 public DoubleVector polarizationEnergyAndGradient(PolarizableMultipoleSIMD mI, PolarizableMultipoleSIMD mK, DoubleVector mutualMask,
178 DoubleVector[] gI, DoubleVector[] tI, DoubleVector[] tK) {
179 DoubleVector emp = gkMonopole.polarizationEnergyAndGradient(mI, mK, one, one, mutualMask, gI, tI, tK);
180 DoubleVector edp = gkDipole.polarizationEnergyAndGradient(mI, mK, one, one, mutualMask, gI, tI, tK);
181 DoubleVector eqp = gkQuadrupole.polarizationEnergyAndGradient(mI, mK, one, one, mutualMask, gI, tI, tK);
182 // Sum the GK polarization interaction energy.
183 return emp.add(edp).add(eqp);
184 }
185
186 /**
187 * Compute the Born chain-rule term for the multipole energy.
188 *
189 * @param mI The polarizable multipole of atom i.
190 * @param mK The polarizable multipole of atom k.
191 * @return The Born chain-rule term.
192 */
193 public DoubleVector multipoleEnergyBornGrad(PolarizableMultipoleSIMD mI, PolarizableMultipoleSIMD mK) {
194 DoubleVector db = gkMonopole.multipoleEnergyBornGrad(mI, mK);
195 db = db.add(gkDipole.multipoleEnergyBornGrad(mI, mK));
196 db = db.add(gkQuadrupole.multipoleEnergyBornGrad(mI, mK));
197 return db;
198 }
199
200 /**
201 * Compute the Born chain-rule term for the polarization energy.
202 *
203 * @param mI The polarizable multipole of atom i.
204 * @param mK The polarizable multipole of atom k.
205 * @param mutual If true, compute the mutual polarization contribution.
206 * @return The Born chain-rule term.
207 */
208 public DoubleVector polarizationEnergyBornGrad(PolarizableMultipoleSIMD mI, PolarizableMultipoleSIMD mK, boolean mutual) {
209 // Compute the GK polarization Born chain-rule term.
210 DoubleVector db = gkMonopole.polarizationEnergyBornGrad(mI, mK);
211 db = db.add(gkDipole.polarizationEnergyBornGrad(mI, mK));
212 db = db.add(gkQuadrupole.polarizationEnergyBornGrad(mI, mK));
213 // Add the mutual polarization contribution to Born chain-rule term.
214 if (mutual) {
215 db = db.add(gkDipole.mutualPolarizationEnergyBornGrad(mI, mK));
216 }
217 return db;
218 }
219
220 }