1#![allow(clippy::excessive_precision)]
2
3use approx::assert_abs_diff_eq;
4use dftd4::prelude::*;
5
6fn test_rational_damping_noargs() {
7 assert!(DFTD4RationalDampingParamBuilder::default().build().is_err());
9
10 let builder = DFTD4RationalDampingParamBuilder::default().a1(0.4).a2(5.0);
11 assert!(builder.build().is_err_and(|e| match &e {
12 DFTD4Error::BuilderError(f) => f.field_name() == "s8",
13 _ => false,
14 }));
15
16 let builder = DFTD4RationalDampingParamBuilder::default().s8(1.0).a2(5.0);
17 assert!(builder.build().is_err_and(|e| match &e {
18 DFTD4Error::BuilderError(f) => f.field_name() == "a1",
19 _ => false,
20 }));
21
22 let builder = DFTD4RationalDampingParamBuilder::default().s8(1.0).a1(0.4);
23 assert!(builder.build().is_err_and(|e| match &e {
24 DFTD4Error::BuilderError(f) => f.field_name() == "a2",
25 _ => false,
26 }));
27}
28
29fn test_structure() {
30 let numbers = vec![6, 7, 6, 7, 6, 6, 6, 8, 7, 6, 8, 7, 6, 6, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1];
32 #[rustfmt::skip]
33 let positions = vec![
34 2.02799738646442, 0.09231312124713, -0.14310895950963,
35 4.75011007621000, 0.02373496014051, -0.14324124033844,
36 6.33434307654413, 2.07098865582721, -0.14235306905930,
37 8.72860718071825, 1.38002919517619, -0.14265542523943,
38 8.65318821103610, -1.19324866489847, -0.14231527453678,
39 6.23857175648671, -2.08353643730276, -0.14218299370797,
40 5.63266886875962, -4.69950321056008, -0.13940509630299,
41 3.44931709749015, -5.48092386085491, -0.14318454855466,
42 7.77508917214346, -6.24427872938674, -0.13107140408805,
43 10.30229550927022, -5.39739796609292, -0.13672168520430,
44 12.07410272485492, -6.91573621641911, -0.13666499342053,
45 10.70038521493902, -2.79078533715849, -0.14148379504141,
46 13.24597858727017, -1.76969072232377, -0.14218299370797,
47 7.40891694074004, -8.95905928176407, -0.11636933482904,
48 1.38702118184179, 2.05575746325296, -0.14178615122154,
49 1.34622199478497, -0.86356704498496, 1.55590600570783,
50 1.34624089204623, -0.86133716815647, -1.84340893849267,
51 5.65596919189118, 4.00172183859480, -0.14131371969009,
52 14.67430918222276, -3.26230980007732, -0.14344911021228,
53 13.50897177220290, -0.60815166181684, 1.54898960808727,
54 13.50780014200488, -0.60614855212345, -1.83214617078268,
55 5.41408424778406, -9.49239668625902, -0.11022772492007,
56 8.31919801555568, -9.74947502841788, 1.56539243085954,
57 8.31511620712388, -9.76854236502758, -1.79108242206824,
58 ];
59
60 assert!(DFTD4Structure::new_f(&numbers, &vec![0.0; 72], None, None, None)
62 .is_err_and(|e| e.to_string().contains("Too close interatomic distances found")));
63
64 assert!(DFTD4Structure::new_f(&[1, 1, 1], &positions, None, None, None)
66 .is_err_and(|e| e.to_string().contains("Invalid dimension for positions")));
67
68 assert!(DFTD4Structure::new_f(&numbers, &[0.0; 7], None, None, None)
70 .is_err_and(|e| e.to_string().contains("Invalid dimension for positions")));
71
72 let mut mol = DFTD4Structure::new(&numbers, &positions, None, None, None);
74
75 assert!(mol
77 .update_f(&[0.0; 7], None)
78 .is_err_and(|e| e.to_string().contains("Invalid dimension for positions")));
79
80 assert!(mol
82 .update_f(&positions, Some(&[0.0; 7]))
83 .is_err_and(|e| e.to_string().contains("Invalid dimension for lattice")));
84
85 assert!(mol
87 .update_f(&[0.0; 72], None)
88 .is_err_and(|e| e.to_string().contains("Too close interatomic distances found"),));
89}
90
91fn test_blypd4() {
92 let thr = 1.0e-7;
94
95 let numbers = vec![1, 1, 6, 5, 1, 15, 8, 17, 13, 15, 5, 1, 9, 15, 1, 15];
96 #[rustfmt::skip]
97 let positions = vec![
98 2.79274810283778, 3.82998228828316, -2.79287054959216,
99 -1.43447454186833, 0.43418729987882, 5.53854345129809,
100 -3.26268343665218, -2.50644032426151, -1.56631149351046,
101 2.14548759959147, -0.88798018953965, -2.24592534506187,
102 -4.30233097423181, -3.93631518670031, -0.48930754109119,
103 0.06107643564880, -3.82467931731366, -2.22333344469482,
104 0.41168550401858, 0.58105573172764, 5.56854609916143,
105 4.41363836635653, 3.92515871809283, 2.57961724984000,
106 1.33707758998700, 1.40194471661647, 1.97530004949523,
107 3.08342709834868, 1.72520024666801, -4.42666116106828,
108 -3.02346932078505, 0.04438199934191, -0.27636197425010,
109 1.11508390868455, -0.97617412809198, 6.25462847718180,
110 0.61938955433011, 2.17903547389232, -6.21279842416963,
111 -2.67491681346835, 3.00175899761859, 1.05038813614845,
112 -4.13181080289514, -2.34226739863660, -3.44356159392859,
113 2.85007173009739, -2.64884892757600, 0.71010806424206,
114 ];
115
116 let model = DFTD4Model::new(&numbers, &positions, None, None, None);
117
118 let param = DFTD4Param::load_rational_damping("blyp", false);
119 let res = model.get_dispersion(¶m, false);
120
121 assert_abs_diff_eq!(res.energy, -0.06991716314879085, epsilon = thr);
122
123 let res = model.get_dispersion(¶m, true);
124
125 assert_abs_diff_eq!(res.energy, -0.06991716314879085, epsilon = thr);
126}
127
128#[cfg(feature = "d4s")]
129fn test_tpssd4s() {
130 let thr = 1.0e-7;
132
133 let numbers = vec![1, 1, 6, 5, 1, 15, 8, 17, 13, 15, 5, 1, 9, 15, 1, 15];
134 #[rustfmt::skip]
135 let positions = vec![
136 2.79274810283778, 3.82998228828316, -2.79287054959216,
137 -1.43447454186833, 0.43418729987882, 5.53854345129809,
138 -3.26268343665218, -2.50644032426151, -1.56631149351046,
139 2.14548759959147, -0.88798018953965, -2.24592534506187,
140 -4.30233097423181, -3.93631518670031, -0.48930754109119,
141 0.06107643564880, -3.82467931731366, -2.22333344469482,
142 0.41168550401858, 0.58105573172764, 5.56854609916143,
143 4.41363836635653, 3.92515871809283, 2.57961724984000,
144 1.33707758998700, 1.40194471661647, 1.97530004949523,
145 3.08342709834868, 1.72520024666801, -4.42666116106828,
146 -3.02346932078505, 0.04438199934191, -0.27636197425010,
147 1.11508390868455, -0.97617412809198, 6.25462847718180,
148 0.61938955433011, 2.17903547389232, -6.21279842416963,
149 -2.67491681346835, 3.00175899761859, 1.05038813614845,
150 -4.13181080289514, -2.34226739863660, -3.44356159392859,
151 2.85007173009739, -2.64884892757600, 0.71010806424206,
152 ];
153
154 let model = DFTD4Model::new_d4s(&numbers, &positions, None, None, None);
155
156 let param = DFTD4Param::load_rational_damping("tpss", false);
157 let res = model.get_dispersion(¶m, false);
158
159 assert_abs_diff_eq!(res.energy, -0.046233140236052253, epsilon = thr);
160
161 let res = model.get_dispersion(¶m, true);
162
163 assert_abs_diff_eq!(res.energy, -0.046233140236052253, epsilon = thr);
164}
165
166fn test_pbed4() {
167 let thr = 1.0e-7;
169
170 let numbers = vec![1, 9, 15, 13, 1, 1, 13, 5, 3, 15, 8, 1, 1, 5, 16, 1];
171 #[rustfmt::skip]
172 let positions = vec![
173 -2.14132037405479, -1.34402701877044, -2.32492500904728,
174 4.46671289205392, -2.04800110524830, 0.44422406067087,
175 -4.92212517643478, -1.73734240529793, 0.96890323821450,
176 -1.30966093045696, -0.52977363497805, 3.44453452239668,
177 -4.34208759006189, -4.30470270977329, 0.39887431726215,
178 0.61788392767516, 2.62484136683297, -3.28228926932647,
179 4.23562873444840, -1.68839322682951, -3.53824299552792,
180 2.23130060612446, 1.93579813100155, -1.80384647554323,
181 -2.32285463652832, 2.90603947535842, -1.39684847191937,
182 2.34557941578250, 2.86074312333371, 1.82827238641666,
183 -3.66431367659153, -0.42910188232667, -1.81957402856634,
184 -0.34927881505446, -1.75988134003940, 5.98017466326572,
185 0.29500802281217, -2.00226104143537, 0.53023447931897,
186 2.10449364205058, -0.56741404446633, 0.30975625014335,
187 -1.59355304432499, 3.69176153150419, 2.87878226787916,
188 4.34858700256050, 2.39171478113440, -2.61802993563738,
189 ];
190
191 let model = DFTD4Model::new(&numbers, &positions, None, None, None);
192
193 let param = DFTD4Param::load_rational_damping("pbe", false);
194 let res = model.get_dispersion(¶m, false);
195
196 assert_abs_diff_eq!(res.energy, -0.028415184156428127, epsilon = thr);
197
198 let res = model.get_dispersion(¶m, true);
199
200 assert_abs_diff_eq!(res.energy, -0.028415184156428127, epsilon = thr);
201}
202
203#[cfg(feature = "d4s")]
204fn test_r2scan3c() {
205 let thr = 1.0e-8;
207
208 let numbers = vec![1, 9, 15, 13, 1, 1, 13, 5, 3, 15, 8, 1, 1, 5, 16, 1];
209 #[rustfmt::skip]
210 let positions = vec![
211 -2.14132037405479, -1.34402701877044, -2.32492500904728,
212 4.46671289205392, -2.04800110524830, 0.44422406067087,
213 -4.92212517643478, -1.73734240529793, 0.96890323821450,
214 -1.30966093045696, -0.52977363497805, 3.44453452239668,
215 -4.34208759006189, -4.30470270977329, 0.39887431726215,
216 0.61788392767516, 2.62484136683297, -3.28228926932647,
217 4.23562873444840, -1.68839322682951, -3.53824299552792,
218 2.23130060612446, 1.93579813100155, -1.80384647554323,
219 -2.32285463652832, 2.90603947535842, -1.39684847191937,
220 2.34557941578250, 2.86074312333371, 1.82827238641666,
221 -3.66431367659153, -0.42910188232667, -1.81957402856634,
222 -0.34927881505446, -1.75988134003940, 5.98017466326572,
223 0.29500802281217, -2.00226104143537, 0.53023447931897,
224 2.10449364205058, -0.56741404446633, 0.30975625014335,
225 -1.59355304432499, 3.69176153150419, 2.87878226787916,
226 4.34858700256050, 2.39171478113440, -2.61802993563738,
227 ];
228
229 let model = DFTD4Model::custom_d4s(&numbers, &positions, 2.0, 1.0, None, None, None);
230
231 let param = DFTD4RationalDampingParamBuilder::default()
232 .s8(0.0)
233 .a1(0.42)
234 .a2(5.65)
235 .s9(2.0)
236 .build()
237 .unwrap()
238 .new_param();
239 let res = model.get_dispersion(¶m, false);
240
241 assert_abs_diff_eq!(res.energy, -0.008016697276824889, epsilon = thr);
242
243 let res = model.get_dispersion(¶m, true);
244
245 assert_abs_diff_eq!(res.energy, -0.008016697276824889, epsilon = thr);
246}
247
248fn test_pair_resolved() {
249 let thr = 1.0e-8;
251
252 let numbers = vec![16, 16, 16, 16, 16, 16, 16, 16];
253 #[rustfmt::skip]
254 let positions = vec![
255 -4.15128787379191, 1.71951973863958, -0.93066267097296,
256 -4.15128787379191, -1.71951973863958, 0.93066267097296,
257 -1.71951973863958, -4.15128787379191, -0.93066267097296,
258 1.71951973863958, -4.15128787379191, 0.93066267097296,
259 4.15128787379191, -1.71951973863958, -0.93066267097296,
260 4.15128787379191, 1.71951973863958, 0.93066267097296,
261 1.71951973863958, 4.15128787379191, -0.93066267097296,
262 -1.71951973863958, 4.15128787379191, 0.93066267097296,
263 ];
264 #[rustfmt::skip]
265 let pair_disp2 = vec![
266 -0.00000000e-0, -5.80599854e-4, -4.74689854e-4, -2.11149449e-4, -1.63163128e-4, -2.11149449e-4, -4.74689854e-4, -5.80599854e-4,
267 -5.80599854e-4, -0.00000000e-0, -5.80599854e-4, -4.74689854e-4, -2.11149449e-4, -1.63163128e-4, -2.11149449e-4, -4.74689854e-4,
268 -4.74689854e-4, -5.80599854e-4, -0.00000000e-0, -5.80599854e-4, -4.74689854e-4, -2.11149449e-4, -1.63163128e-4, -2.11149449e-4,
269 -2.11149449e-4, -4.74689854e-4, -5.80599854e-4, -0.00000000e-0, -5.80599854e-4, -4.74689854e-4, -2.11149449e-4, -1.63163128e-4,
270 -1.63163128e-4, -2.11149449e-4, -4.74689854e-4, -5.80599854e-4, -0.00000000e-0, -5.80599854e-4, -4.74689854e-4, -2.11149449e-4,
271 -2.11149449e-4, -1.63163128e-4, -2.11149449e-4, -4.74689854e-4, -5.80599854e-4, -0.00000000e-0, -5.80599854e-4, -4.74689854e-4,
272 -4.74689854e-4, -2.11149449e-4, -1.63163128e-4, -2.11149449e-4, -4.74689854e-4, -5.80599854e-4, -0.00000000e-0, -5.80599854e-4,
273 -5.80599854e-4, -4.74689854e-4, -2.11149449e-4, -1.63163128e-4, -2.11149449e-4, -4.74689854e-4, -5.80599854e-4, -0.00000000e-0,
274 ];
275 #[rustfmt::skip]
276 let pair_disp3 = vec![
277 0.00000000e-0, 3.39353850e-7, 8.74462839e-7, 1.17634100e-6, 9.86937310e-7, 1.17634100e-6, 8.74462839e-7, 3.39353850e-7,
278 3.39353850e-7, 0.00000000e-0, 3.39353850e-7, 8.74462839e-7, 1.17634100e-6, 9.86937310e-7, 1.17634100e-6, 8.74462839e-7,
279 8.74462839e-7, 3.39353850e-7, 0.00000000e-0, 3.39353850e-7, 8.74462839e-7, 1.17634100e-6, 9.86937310e-7, 1.17634100e-6,
280 1.17634100e-6, 8.74462839e-7, 3.39353850e-7, 0.00000000e-0, 3.39353850e-7, 8.74462839e-7, 1.17634100e-6, 9.86937310e-7,
281 9.86937310e-7, 1.17634100e-6, 8.74462839e-7, 3.39353850e-7, 0.00000000e-0, 3.39353850e-7, 8.74462839e-7, 1.17634100e-6,
282 1.17634100e-6, 9.86937310e-7, 1.17634100e-6, 8.74462839e-7, 3.39353850e-7, 0.00000000e-0, 3.39353850e-7, 8.74462839e-7,
283 8.74462839e-7, 1.17634100e-6, 9.86937310e-7, 1.17634100e-6, 8.74462839e-7, 3.39353850e-7, 0.00000000e-0, 3.39353850e-7,
284 3.39353850e-7, 8.74462839e-7, 1.17634100e-6, 9.86937310e-7, 1.17634100e-6, 8.74462839e-7, 3.39353850e-7, 0.00000000e-0,
285 ];
286
287 let model = DFTD4Model::new(&numbers, &positions, None, None, None);
288
289 let param = DFTD4RationalDampingParamBuilder::default()
290 .s8(1.20065498)
291 .a1(0.40085597)
292 .a2(5.02928789)
293 .build()
294 .unwrap()
295 .new_param();
296 let res = model.get_pairwise_dispersion(¶m);
297
298 for (a, b) in res.pair_energy2.iter().zip(pair_disp2.iter()) {
299 assert_abs_diff_eq!(a, b, epsilon = thr);
300 }
301
302 for (a, b) in res.pair_energy3.iter().zip(pair_disp3.iter()) {
303 assert_abs_diff_eq!(a, b, epsilon = thr);
304 }
305}
306
307fn test_properties() {
308 let thr = 1.0e-7;
310
311 let numbers = vec![7, 7, 7, 7, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1];
312 #[rustfmt::skip]
313 let positions = vec![
314 1.97420621099560, 1.97415497783241, 1.97424596974304,
315 6.82182427659395, 2.87346383480995, 7.72099517560089,
316 7.72104957181201, 6.82177051521773, 2.87336561318016,
317 2.87343220660781, 7.72108897828386, 6.82187093171878,
318 3.51863272100286, 2.63865333484548, 1.00652979981286,
319 2.63877594964754, 1.00647313885594, 3.51882748086447,
320 1.00639728563189, 3.51850454450845, 2.63869202592387,
321 8.36624975982697, 2.20896711017229, 8.68870955681018,
322 7.48639684558259, 3.84114715917956, 6.17640982573725,
323 5.85401675167715, 1.32911569888797, 7.05654606696031,
324 7.05646299938990, 5.85409590282274, 1.32879923864813,
325 8.68882633853582, 8.36611541129785, 2.20894120662207,
326 3.84121223226912, 6.17673669892998, 7.48629723649480,
327 1.32897854262127, 7.05658604099926, 5.85414031368096,
328 2.20884896069885, 8.68875820985799, 8.36643568423387,
329 6.17659142004652, 7.48627051643848, 3.84109594690835,
330 ];
331 #[rustfmt::skip]
332 let lattice = vec![
333 9.69523775911749, 0.0, 0.0,
334 0.0, 9.69523775911749, 0.0,
335 0.0, 0.0, 9.69523775911749,
336 ];
337 #[rustfmt::skip]
338 let cn = [
339 2.5775156287150218,
340 2.5775155620078536,
341 2.5775157938667150,
342 2.5775157704485387,
343 0.8591731475439074,
344 0.8591680526841657,
345 0.8591744284869478,
346 0.8591732359038715,
347 0.8591678667769283,
348 0.8591744593270527,
349 0.8591684383407867,
350 0.8591754863625011,
351 0.8591751690053771,
352 0.8591719636497469,
353 0.8591686377934131,
354 0.8591718691634261,
355 ];
356 #[rustfmt::skip]
357 let charges = [
358 -0.86974543285813199,
359 -0.86974501326130316,
360 -0.86974658178316333,
361 -0.86974643466661739,
362 0.28992010784471012,
363 0.28989847135893759,
364 0.28992738637932841,
365 0.28992042841052362,
366 0.28989746859382759,
367 0.28992753735979965,
368 0.28990032127437559,
369 0.28993196014690176,
370 0.28993044356403513,
371 0.28991421277784613,
372 0.28990119160907674,
373 0.28991393324985348,
374 ];
375 #[rustfmt::skip]
376 let alpha = [
377 9.9853045768095097,
378 9.9853025181287016,
379 9.9853102162086920,
380 9.9853094944097140,
381 1.3513315505023076,
382 1.3513939255952379,
383 1.3513106323935196,
384 1.3513306244831862,
385 1.3513968091133417,
386 1.3513101978580895,
387 1.3513885997538553,
388 1.3512974505020905,
389 1.3513018164322617,
390 1.3513485145877058,
391 1.3513860914395939,
392 1.3513493246534483,
393 ];
394
395 let model = DFTD4Model::new(&numbers, &positions, None, Some(&lattice), None);
396
397 let res = model.get_properties();
398
399 for (a, b) in res.cn.iter().zip(cn.iter()) {
400 assert_abs_diff_eq!(a, b, epsilon = thr);
401 }
402
403 for (a, b) in res.charges.iter().zip(charges.iter()) {
404 assert_abs_diff_eq!(a, b, epsilon = thr);
405 }
406
407 for (a, b) in res.alpha.iter().zip(alpha.iter()) {
408 assert_abs_diff_eq!(a, b, epsilon = thr);
409 }
410}
411
412#[test]
413fn test() {
414 test_rational_damping_noargs();
415 test_structure();
416 test_blypd4();
417 #[cfg(feature = "d4s")]
418 test_tpssd4s();
419 test_pbed4();
420 #[cfg(feature = "d4s")]
421 test_r2scan3c();
422 test_pair_resolved();
423 test_properties();
424}
425
426fn main() {
427 test_rational_damping_noargs();
428 test_structure();
429 test_blypd4();
430 #[cfg(feature = "d4s")]
431 test_tpssd4s();
432 test_pbed4();
433 #[cfg(feature = "d4s")]
434 test_r2scan3c();
435 test_pair_resolved();
436 test_properties();
437}