Skip to content

Commit 8cdc798

Browse files
momokchungmomokchung
committed
Added exchpol tests, fixed bug in pair_alterpol.h
1 parent b9c319f commit 8cdc798

File tree

14 files changed

+756
-7
lines changed

14 files changed

+756
-7
lines changed

‎include/seq/pair_alterpol.h‎

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -82,18 +82,18 @@ inline void pair_dexpol(ExpolScr scrtyp, real r, real r2, real pscale, real cut,
8282
zr = 0.;
8383
real dot = ai[0][2];
8484
real eps = 1. / REAL_SQRT(2);
85-
if (abs(dot) >= eps) {
85+
if (fabs(dot) > eps) {
8686
xr = 0.;
8787
yr = 1.;
8888
dot = ai[1][2];
8989
}
9090
xr = xr - dot * ai[0][2];
9191
yr = yr - dot * ai[1][2];
9292
zr = zr - dot * ai[2][2];
93-
r = REAL_SQRT(xr * xr + yr * yr + zr * zr);
94-
ai[0][0] = xr / r;
95-
ai[1][0] = yr / r;
96-
ai[2][0] = zr / r;
93+
real dr = REAL_SQRT(xr * xr + yr * yr + zr * zr);
94+
ai[0][0] = xr / dr;
95+
ai[1][0] = yr / dr;
96+
ai[2][0] = zr / dr;
9797
ai[0][1] = ai[2][0] * ai[1][2] - ai[1][0] * ai[2][2];
9898
ai[1][1] = ai[0][0] * ai[2][2] - ai[2][0] * ai[0][2];
9999
ai[2][1] = ai[1][0] * ai[0][2] - ai[0][0] * ai[1][2];

‎test/aplusgas.cpp‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -21,7 +21,7 @@ TEST_CASE("APlus-Gas-Alyz", "[ff][aplus]")
2121
auto ref_g = r.getGradient();
2222

2323
const double eps_e = testGetEps(0.0030, 0.0001);
24-
const double eps_g = testGetEps(0.0010, 0.0001);
24+
const double eps_g = testGetEps(0.0013, 0.0001);
2525
const double eps_v = testGetEps(0.002, 0.001);
2626

2727
const char* argv[] = {"dummy", xn};
@@ -77,7 +77,7 @@ TEST_CASE("APlus-Gas", "[ff][aplus]")
7777
auto ref_g = r.getGradient();
7878

7979
const double eps_e = testGetEps(0.0030, 0.0001);
80-
const double eps_g = testGetEps(0.0010, 0.0001);
80+
const double eps_g = testGetEps(0.0013, 0.0001);
8181
const double eps_v = testGetEps(0.002, 0.001);
8282

8383
const char* argv[] = {"dummy", xn};

‎test/cmakesrc.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@ chgtrn.cpp
1313
disp.cpp
1414
emhippo.cpp
1515
ephippo.cpp
16+
expol.cpp
1617
geom.cpp
1718
improp.cpp
1819
imptor.cpp

‎test/expol.cpp‎

Lines changed: 230 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,230 @@
1+
#include "ff/amoebamod.h"
2+
3+
#include "test.h"
4+
#include "testrt.h"
5+
6+
using namespace tinker;
7+
8+
TEST_CASE("ExchangePolarization-1", "[ff][hippo][expol]")
9+
{
10+
TestFile fx1(TINKER9_DIRSTR "/test/file/expol/NaCl.xyz");
11+
TestFile fk1(TINKER9_DIRSTR "/test/file/expol/expol.key");
12+
TestFile fp1(TINKER9_DIRSTR "/test/file/expol/expol.prm");
13+
14+
const char* xn = "NaCl.xyz";
15+
const char* kn = "expol.key";
16+
const char* argv[] = {"dummy", xn, "-k", kn};
17+
int argc = 4;
18+
19+
const double eps_e = testGetEps(0.0001, 0.0001);
20+
const double eps_g = testGetEps(0.0001, 0.0001);
21+
const double eps_v = testGetEps(0.001, 0.001);
22+
23+
TestReference r(TINKER9_DIRSTR "/test/ref/expol.1.txt");
24+
auto ref_c = r.getCount();
25+
auto ref_e = r.getEnergy();
26+
auto ref_v = r.getVirial();
27+
auto ref_g = r.getGradient();
28+
29+
rc_flag = calc::xyz | calc::vmask;
30+
testBeginWithArgs(argc, argv);
31+
initialize();
32+
33+
energy(calc::v0);
34+
COMPARE_REALS(esum, ref_e, eps_e);
35+
36+
energy(calc::v1);
37+
COMPARE_REALS(esum, ref_e, eps_e);
38+
COMPARE_GRADIENT(ref_g, eps_g);
39+
for (int i = 0; i < 3; ++i)
40+
for (int j = 0; j < 3; ++j)
41+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
42+
43+
energy(calc::v3);
44+
COMPARE_REALS(esum, ref_e, eps_e);
45+
COMPARE_INTS(countReduce(nep), ref_c);
46+
47+
energy(calc::v4);
48+
COMPARE_REALS(esum, ref_e, eps_e);
49+
COMPARE_GRADIENT(ref_g, eps_g);
50+
51+
energy(calc::v5);
52+
COMPARE_GRADIENT(ref_g, eps_g);
53+
54+
energy(calc::v6);
55+
COMPARE_GRADIENT(ref_g, eps_g);
56+
for (int i = 0; i < 3; ++i)
57+
for (int j = 0; j < 3; ++j)
58+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
59+
60+
finish();
61+
testEnd();
62+
}
63+
64+
TEST_CASE("ExchangePolarization-2", "[ff][hippo][expol]")
65+
{
66+
TestFile fx1(TINKER9_DIRSTR "/test/file/expol/Nawater.xyz");
67+
TestFile fk1(TINKER9_DIRSTR "/test/file/expol/expol.key");
68+
TestFile fp1(TINKER9_DIRSTR "/test/file/expol/expol.prm");
69+
70+
const char* xn = "Nawater.xyz";
71+
const char* kn = "expol.key";
72+
const char* argv[] = {"dummy", xn, "-k", kn};
73+
int argc = 4;
74+
75+
const double eps_e = testGetEps(0.0001, 0.0001);
76+
const double eps_g = testGetEps(0.0001, 0.0001);
77+
const double eps_v = testGetEps(0.001, 0.001);
78+
79+
TestReference r(TINKER9_DIRSTR "/test/ref/expol.2.txt");
80+
auto ref_c = r.getCount();
81+
auto ref_e = r.getEnergy();
82+
auto ref_v = r.getVirial();
83+
auto ref_g = r.getGradient();
84+
85+
rc_flag = calc::xyz | calc::vmask;
86+
testBeginWithArgs(argc, argv);
87+
initialize();
88+
89+
energy(calc::v0);
90+
COMPARE_REALS(esum, ref_e, eps_e);
91+
92+
energy(calc::v1);
93+
COMPARE_REALS(esum, ref_e, eps_e);
94+
COMPARE_GRADIENT(ref_g, eps_g);
95+
for (int i = 0; i < 3; ++i)
96+
for (int j = 0; j < 3; ++j)
97+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
98+
99+
energy(calc::v3);
100+
COMPARE_REALS(esum, ref_e, eps_e);
101+
COMPARE_INTS(countReduce(nep), ref_c);
102+
103+
energy(calc::v4);
104+
COMPARE_REALS(esum, ref_e, eps_e);
105+
COMPARE_GRADIENT(ref_g, eps_g);
106+
107+
energy(calc::v5);
108+
COMPARE_GRADIENT(ref_g, eps_g);
109+
110+
energy(calc::v6);
111+
COMPARE_GRADIENT(ref_g, eps_g);
112+
for (int i = 0; i < 3; ++i)
113+
for (int j = 0; j < 3; ++j)
114+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
115+
116+
finish();
117+
testEnd();
118+
}
119+
120+
TEST_CASE("ExchangePolarization-3", "[ff][hippo][expol]")
121+
{
122+
TestFile fx1(TINKER9_DIRSTR "/test/file/expol/Clwater.xyz");
123+
TestFile fk1(TINKER9_DIRSTR "/test/file/expol/expol.key");
124+
TestFile fp1(TINKER9_DIRSTR "/test/file/expol/expol.prm");
125+
126+
const char* xn = "Clwater.xyz";
127+
const char* kn = "expol.key";
128+
const char* argv[] = {"dummy", xn, "-k", kn};
129+
int argc = 4;
130+
131+
const double eps_e = testGetEps(0.0001, 0.0001);
132+
const double eps_g = testGetEps(0.0001, 0.0001);
133+
const double eps_v = testGetEps(0.001, 0.001);
134+
135+
TestReference r(TINKER9_DIRSTR "/test/ref/expol.3.txt");
136+
auto ref_c = r.getCount();
137+
auto ref_e = r.getEnergy();
138+
auto ref_v = r.getVirial();
139+
auto ref_g = r.getGradient();
140+
141+
rc_flag = calc::xyz | calc::vmask;
142+
testBeginWithArgs(argc, argv);
143+
initialize();
144+
145+
energy(calc::v0);
146+
COMPARE_REALS(esum, ref_e, eps_e);
147+
148+
energy(calc::v1);
149+
COMPARE_REALS(esum, ref_e, eps_e);
150+
COMPARE_GRADIENT(ref_g, eps_g);
151+
for (int i = 0; i < 3; ++i)
152+
for (int j = 0; j < 3; ++j)
153+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
154+
155+
energy(calc::v3);
156+
COMPARE_REALS(esum, ref_e, eps_e);
157+
COMPARE_INTS(countReduce(nep), ref_c);
158+
159+
energy(calc::v4);
160+
COMPARE_REALS(esum, ref_e, eps_e);
161+
COMPARE_GRADIENT(ref_g, eps_g);
162+
163+
energy(calc::v5);
164+
COMPARE_GRADIENT(ref_g, eps_g);
165+
166+
energy(calc::v6);
167+
COMPARE_GRADIENT(ref_g, eps_g);
168+
for (int i = 0; i < 3; ++i)
169+
for (int j = 0; j < 3; ++j)
170+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
171+
172+
finish();
173+
testEnd();
174+
}
175+
176+
TEST_CASE("ExchangePolarization-4", "[ff][hippo][expol]")
177+
{
178+
TestFile fx1(TINKER9_DIRSTR "/test/file/expol/water2Na2Clbox.xyz");
179+
TestFile fk1(TINKER9_DIRSTR "/test/file/expol/expol.key");
180+
TestFile fp1(TINKER9_DIRSTR "/test/file/expol/expol.prm");
181+
182+
const char* xn = "water2Na2Clbox.xyz";
183+
const char* kn = "expol.key";
184+
const char* argv[] = {"dummy", xn, "-k", kn};
185+
int argc = 4;
186+
187+
const double eps_e = testGetEps(0.0001, 0.0001);
188+
const double eps_g = testGetEps(0.0001, 0.0001);
189+
const double eps_v = testGetEps(0.001, 0.001);
190+
191+
TestReference r(TINKER9_DIRSTR "/test/ref/expol.4.txt");
192+
auto ref_c = r.getCount();
193+
auto ref_e = r.getEnergy();
194+
auto ref_v = r.getVirial();
195+
auto ref_g = r.getGradient();
196+
197+
rc_flag = calc::xyz | calc::vmask;
198+
testBeginWithArgs(argc, argv);
199+
initialize();
200+
201+
energy(calc::v0);
202+
COMPARE_REALS(esum, ref_e, eps_e);
203+
204+
energy(calc::v1);
205+
COMPARE_REALS(esum, ref_e, eps_e);
206+
COMPARE_GRADIENT(ref_g, eps_g);
207+
for (int i = 0; i < 3; ++i)
208+
for (int j = 0; j < 3; ++j)
209+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
210+
211+
energy(calc::v3);
212+
COMPARE_REALS(esum, ref_e, eps_e);
213+
COMPARE_INTS(countReduce(nep), ref_c);
214+
215+
energy(calc::v4);
216+
COMPARE_REALS(esum, ref_e, eps_e);
217+
COMPARE_GRADIENT(ref_g, eps_g);
218+
219+
energy(calc::v5);
220+
COMPARE_GRADIENT(ref_g, eps_g);
221+
222+
energy(calc::v6);
223+
COMPARE_GRADIENT(ref_g, eps_g);
224+
for (int i = 0; i < 3; ++i)
225+
for (int j = 0; j < 3; ++j)
226+
COMPARE_REALS(vir[i * 3 + j], ref_v[i][j], eps_v);
227+
228+
finish();
229+
testEnd();
230+
}

‎test/file/expol/Clwater.xyz‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,5 @@
1+
4 Cl-_water
2+
1 O 0.00000000 0.00000000 0.00000000 1 2 3
3+
2 H -0.05052100 0.95541000 0.00000000 2 1
4+
3 H 0.97780300 -0.13434900 0.00000000 2 1
5+
4 Cl- 3.08435800 0.00000000 0.00000000 15

‎test/file/expol/NaCl.xyz‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,3 @@
1+
2 Na+_Cl-
2+
1 Na+ 0.000000 0.000000 0.000000 7
3+
2 Cl- 0.000000 0.000000 2.370000 15

‎test/file/expol/Nawater.xyz‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,5 @@
1+
4 Na+_water
2+
1 O 0.00000000 0.00000000 0.00000000 1 2 3
3+
2 H -0.58876200 0.75914800 0.00000000 2 1
4+
3 H -0.58876200 -0.75914800 0.00000000 2 1
5+
4 Na+ 2.21989500 0.00000000 0.00000000 7

‎test/file/expol/expol.key‎

Lines changed: 25 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,25 @@
1+
parameters expol.prm
2+
3+
exchpol 2 7.57187611 1.00774942 22.92257722 0
4+
exchpol 1 8.52625781 0.26279239 2.34042498 0
5+
exchpol 7 1.28799890 14.03512371 7.21629192 1
6+
exchpol 15 0.80994925 2.03295455 1.28410603 1
7+
exchpol 8 1.09975765 6.89625706 4.51567248 1
8+
exchpol 16 0.90156301 2.17737556 1.28742270 1
9+
10+
chgpen 7 4.090412102728663 5.346221041987701
11+
chgpen 15 9.937217338945683 3.403540591081243
12+
chgpen 8 8.689356786155788 4.550459010432590
13+
chgpen 16 12.129167017481114 3.222907140342866
14+
chgpen 20 6.298463497273733 6.755370803213201
15+
chgpen 21 12.760026532535434 5.683167591839673
16+
17+
polarize 1 1.21889959 2
18+
polarize 2 0.14304110 1
19+
20+
screening s2u
21+
digits 10
22+
23+
polar-eps 1.0E-14
24+
25+
polarizeterm only

0 commit comments

Comments
 (0)