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  import static org.junit.Assert.assertEquals;
41  
42  import java.util.Arrays;
43  import java.util.Collection;
44  
45  import ffx.utilities.FFXTest;
46  import org.junit.Test;
47  import org.junit.runner.RunWith;
48  import org.junit.runners.Parameterized;
49  import org.junit.runners.Parameterized.Parameters;
50  
51  /**
52   * Parameterized Test of the MultipoleTensor class.
53   *
54   * @author Michael J. Schnieders
55   * @since 1.0
56   */
57  @RunWith(Parameterized.class)
58  public class MultipoleTensorTest extends FFXTest {
59  
60    private final int order;
61    private final int tensorCount;
62    private final String info;
63  
64    // Water Dimer O-O interaction
65    protected static final double[] xyz = {2.97338529, 0.0, 0.03546452, 0.0};
66    // Rotated multipole at site I.
67    protected final static double[] Qi = {-0.51966,
68        0.06979198988239577, 0.0, 0.0289581620819011,
69        0.024871041044109393, -0.1170771231287098, 0.09220608208460039,
70        0.0, -0.03374891685535346, 0.0};
71    // Rotated multipole at site K.
72    protected final static double[] Qk = {-0.51966,
73        0.05872406108747119, 0.0, 0.047549780780788455,
74        0.048623413357888695, -0.1170771231287098, 0.06845370977082109,
75        0.0, -0.04662811558421081, 0.0};
76  
77    // Permanent multipole Coulomb energy and torque.
78    protected final static double permanentEnergy = 0.08861401274679;
79    protected final static double[] permTorqueI = {0.0, 0.00039200904382, 0.0};
80    protected final static double[] permTorqueK = {0.0, -0.00400122461628, 0.0};
81    // Permanent multipole Ewald energy and torque.
82    protected final static double permanentEnergyEwald = 0.001275693067120;
83    protected final static double[] permTorqueIEwald = {0.0, -0.000304652548761, 0.0};
84    protected final static double[] permTorqueKEwald = {0.0, -0.000745949268609, 0.0};
85  
86    // Water Dimer O-O Coulomb mutual polarization Energy, grad and torque in vacuum
87    protected final static double scaleMutual = 1.0;
88    protected final static double[] Ui = {0.04886563833303603, 0.0, -0.0018979726219775425};
89    protected final static double[] Uk = {-0.040839567654139396, 0.0, -5.982126263609587E-4};
90    protected final static double polarizationEnergyCoulomb = -0.002576831234958;
91    protected final static double[] polarGradICoulomb = {-0.003233364149133, 0.0, 0.000081374809911};
92    protected final static double[] polarTorqueICoulomb = {0.0, 0.000015759704950, 0.0};
93    protected final static double[] polarTorqueKCoulomb = {0.0, 0.000340863064066, 0.0};
94    // Water Dimer O-O Thole correction for mutual Polarization energy, grad and torque in vacuum.
95    protected final static double polarizationEnergyThole = 0.000000036672468;
96    protected final static double[] polarGradIThole = {0.000001005045968, 0.0, 0.000000090221488};
97    protected final static double[] polarTorqueIThole = {0.0, 0.000000072500602, 0.0};
98    protected final static double[] polarTorqueKThole = {0.0, 0.000000151677615, 0.0};
99  
100   // Water Dimer O-O Coulomb mutual polarization Energy, grad and torque under Ewald.
101   protected final static double[] UiEwald = {0.017212550663556914, 0.0, 0.0026849371538929636};
102   protected final static double[] UkEwald = {0.0020051189004858657, 0.0, 0.005476556978973183};
103   protected final static double polarizationEnergyEwald = -0.000078870483232;
104   protected final static double[] polarGradIEwald = {-0.000345793382509, 0.0, 0.000025382496373};
105   protected final static double[] polarTorqueIEwald = {0.0, -0.000003385814434, 0.0};
106   protected final static double[] polarTorqueKEwald = {0.0, 0.000080790360336, 0.0};
107 
108   // AMOEBA+ Water Dimer O-O interaction
109   protected final static double amoebaPlusMPoleEnergyOO = 8.6959247973876330E-002;
110   protected final static double amoebaPlusMPoleEnergyOH = -3.3755171537654061E-002;
111   protected final static double amoebaPlusIndDipoleEnergyOO = 1.4935109017088957E-003/2;
112   protected final static double amoebaPlusIndDipoleEnergyOH = -0.2283;
113   protected final static double amoebaPlusChargeTransferEnergy = -0.0427;
114   protected final static double chargePenAlphaOx = 4.0047;
115   protected final static double chargePenAlphaHyd = 3.2541;
116   protected final static double[] chrgTranParamsOx = new double[]{9.4879, 3.8982}; // alpha, beta?
117   // Multipole & Induced dipoles (undo what PolarizableMultipole Applies - 1/3
118   // (tinker doesn't double, so we need to))
119   protected final static double[] QiXYZ = {1.467089, 0.067722, 0.000000};
120   protected final static double Zi = 8.0;
121   // Oxygen 1
122   protected final static double[] QiAmoebaP = {
123           -0.50316456711796764, // charge
124           -4.8391389372990617E-002,  -9.9697768378871676E-002,  0.0000000000000000, // dipole
125           1.2064199052004655E-002*3,   1.7432698501819319E-002*3,  -2.9496897553823977E-002*3, // trace
126           3.4088807853371586E-003*3,   0.0000000000000000,        0.0000000000000000 // off diag
127   };
128   protected final static double[] QiInduced = {
129           -1.9757159038393233E-002,
130           -1.1052620863319413E-002,
131           -2.8861245105701660E-021,
132   };
133   protected final static double[] QkXYZ = {-1.405210, -0.064738, 0.000000};
134   protected final static double Zk = 8.0;
135   // Oxygen 2
136   protected final static double[] QkAmoebaP = {
137           -0.50645532267523641,
138           -5.6497094042130665E-002,   9.5338606817858315E-002,   0.0000000000000000,
139           -1.6869890623243792E-002*3,   6.4602970864044542E-003*3,   1.0409593536839330E-002*3,
140           -2.1308020694714599E-002*3,   -8.9034686819118132E-019*3,
141           3.2464171803135626E-018*3
142   };
143   protected final static double[] QkInduced = {
144           -5.3919378773639018E-002,
145            6.6658070142521734E-003,
146           -9.1447209902403356E-022,
147   };
148   // Hydrogen 2 (H-O interaction)
149   protected final static double[] QjAmoebaP = {
150           0.25152712323368620,
151           0.11140287316108205,       -4.7139571392261465E-002,   0.0000000000000000,
152           1.1431494801174508E-004*3,   3.6656033650649092E-003*3,  -3.7799183130766555E-003*3,
153           9.0058239922283824E-003*3,   0.0000000000000000,        0.0000000000000000
154   };
155   protected final static double[] QjInduced = {-0.0461,0.0348,-0.0000};
156   protected final static double[] QjXYZ = {1.887074, -0.783266, 0.000000};
157   protected final static double Zj = 1.0;
158 
159   // Oxy - 1
160   protected final static double[] QlAPlusEwal = {-0.51294009372635441,
161           -4.7449978254291744E-002,  -2.3671746096901102E-002,   9.7311456599758192E-002,
162           5.8029182615047855E-003*3,  -1.6704434050841255E-002*3,   1.0901515789336465E-002*3,
163           -1.2263684370259969E-002*3, -9.4608247580866937E-003*3, -1.4686497793089748E-002*3};
164   // Hydrogen - 536
165   protected final static double[] QmAplusEwal = {0.26162508784642219,
166           -0.10389679865169722,       -2.6854217213802947E-002,  -5.5830507675936697E-002,
167           3*-6.6549319715359873E-003,  9.8396726583062210E-003*3.0,  -3.1847406867702320E-003*3.0,
168           3*-7.8527872497808878E-004,  -1.5585386297907382E-003*3,   4.0376954110960461E-003*3};
169   protected final static double[] xyzEwaldAPlus = new double[]
170           {-3.2149529999999999,        2.1452700000000000,       -2.0970249999999999};
171 
172   public MultipoleTensorTest(String info, int order) {
173     this.info = info;
174     this.order = order;
175     this.tensorCount = MultipoleUtilities.tensorCount(order);
176   }
177 
178   @Parameters
179   public static Collection<Object[]> data() {
180     return Arrays.asList(
181         new Object[][] {
182             {
183                 "Order 6", 6
184             },
185             {
186                 "Order 5", 5
187             },
188             {
189                 "Order 4", 4
190             }
191         });
192   }
193 
194   /**
195    * Test of tensorCount method, of class MultipoleTensor.
196    *
197    * @since 1.0
198    */
199   @Test
200   public void tensorCountTest() {
201     int result = MultipoleUtilities.tensorCount(order);
202     assertEquals(info, tensorCount, result);
203   }
204 
205   /** Test of ti method, of class MultipoleTensor. */
206   @Test
207   public void tensorIndexTest() {
208     int dx = 1;
209     int dy = 0;
210     int dz = 0;
211     int expResult = 1;
212     int result = MultipoleUtilities.ti(dx, dy, dz, order);
213     assertEquals(info, expResult, result);
214   }
215 
216 }