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