From 2cb74fea71ff5baa27f3c9f12c06ab3d2aec1bba Mon Sep 17 00:00:00 2001 From: Felix Pultar Date: Thu, 17 Nov 2022 20:08:26 +0100 Subject: [PATCH 1/2] Increased significant digits for point charge gradients to 12 Signed-off-by: Felix Pultar --- src/xtb/calculator.f90 | 2 +- test/api/c_api_example.c | 29 +++++++++++++++++++++++++++-- 2 files changed, 28 insertions(+), 3 deletions(-) diff --git a/src/xtb/calculator.f90 b/src/xtb/calculator.f90 index 799cd9ac5..b2da1cf2b 100644 --- a/src/xtb/calculator.f90 +++ b/src/xtb/calculator.f90 @@ -282,7 +282,7 @@ subroutine singlepoint(self, env, mol, chk, printlevel, restart, & if (allocated(set%pcem_grad) .and. self%pcem%n > 0) then call open_file(ich,set%pcem_grad,'w') do i=1,self%pcem%n - write(ich,'(3f12.8)')self%pcem%grd(1:3,i) + write(ich,'(3f16.12)')self%pcem%grd(1:3,i) enddo call close_file(ich) endif diff --git a/test/api/c_api_example.c b/test/api/c_api_example.c index a4c04ee49..353d7b463 100644 --- a/test/api/c_api_example.c +++ b/test/api/c_api_example.c @@ -227,6 +227,15 @@ int testSecond() { 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7}; // results will live in here + double energy = 0.0; + double grad[3 * 31] = { + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; double pcgrad[3 * 32] = { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, @@ -260,6 +269,14 @@ int testSecond() { if (xtb_checkEnvironment(env)) goto error; + xtb_getEnergy(env, res, &energy); + if (xtb_checkEnvironment(env)) + goto error; + + xtb_getGradient(env, res, grad); + if (xtb_checkEnvironment(env)) + goto error; + xtb_getPCGradient(env, res, pcgrad); if (xtb_checkEnvironment(env)) goto error; @@ -273,9 +290,17 @@ int testSecond() { xtb_delete(mol); xtb_delete(env); - if (!check(pcgrad[0], 0.00000755, 1.0e-6, "pcgrad[0] does not match")) + if (!check(energy, -36.669892958754, 1.0e-9, "Energy does not match")) + goto error; + + if (!check(grad[0], -0.001992161613, 1.0e-8, "grad[0] does not match")) + goto error; + if (!check(grad[92], 0.000556642950, 1.0e-8, "grad[92] does not match")) + goto error; + + if (!check(pcgrad[0], 0.000007545625, 1.0e-8, "pcgrad[0] does not match")) goto error; - if (!check(pcgrad[95], 0.00001312, 1.0e-6, "pcgrad[95] does not match")) + if (!check(pcgrad[95], 0.000013117708, 1.0e-8, "pcgrad[95] does not match")) goto error; return 0; From 9fee0290143be0921856fe9f91b5e358a6572834 Mon Sep 17 00:00:00 2001 From: Felix Pultar Date: Thu, 17 Nov 2022 20:51:43 +0100 Subject: [PATCH 2/2] Second unit test in C API is water tetramers: https://xtb-docs.readthedocs.io/en/latest/pcem.html. Be more consistent with Fortran unit tests Signed-off-by: Felix Pultar --- test/api/c_api_example.c | 126 +++++++++++++++------------------------ 1 file changed, 48 insertions(+), 78 deletions(-) diff --git a/test/api/c_api_example.c b/test/api/c_api_example.c index 353d7b463..bc5e8b61a 100644 --- a/test/api/c_api_example.c +++ b/test/api/c_api_example.c @@ -162,88 +162,50 @@ int testFirst() { } int testSecond() { - // second molecule + // water tetramer : https://xtb-docs.readthedocs.io/en/latest/pcem.html - int const natoms = 31; - int const attyp[31] = {1, 6, 1, 1, 6, 1, 6, 1, 1, 1, 6, 1, 6, 1, 1, 6, - 1, 1, 6, 1, 6, 1, 1, 1, 6, 1, 1, 6, 1, 8, 1}; + int const natoms = 6; + int const attyp[6] = {8, 1, 1, 8, 1, 1}; double const charge = 0.0; int const uhf = 0; - double const coord[3 * 31] = { - -0.442947496, 11.001947210, 23.53843018, -2.236268461, 11.818985980, - 22.93889444, -1.841717705, 13.792091510, 22.49830981, -2.896722482, - 10.768954200, 21.29444454, -4.133665133, 11.695821690, 25.11713090, - -4.513283393, 9.713842463, 25.52960654, -6.643232696, 12.675979330, - 24.06753249, -6.497152995, 14.672778630, 23.58352390, -8.148942688, - 12.705927690, 25.47277389, -7.353363183, 11.592564650, 22.46591770, - -3.249632431, 12.876986410, 27.66935330, -4.846956553, 12.734322190, - 28.99201199, -2.467110848, 15.660343880, 27.33890614, -4.009123687, - 16.787220770, 26.56744596, -0.818133300, 15.913359350, 26.13075930, - -1.719833938, 16.708905520, 29.92782456, -1.175383185, 18.674510890, - 29.63963875, -3.270175176, 16.843998070, 31.30603239, 0.393734165, - 15.240980550, 31.24575662, 2.098110773, 15.276181590, 30.05626854, - 1.103721260, 16.300985820, 33.84032285, -0.510552700, 16.227748650, - 35.11752635, 1.786206150, 18.232492740, 33.62584759, 2.541275947, - 15.177545230, 34.79620160, -0.431661103, 12.490133300, 31.57724434, - 1.201728308, 11.377478740, 32.22095010, -1.982252711, 12.263731930, - 32.94292201, -1.094745099, 11.448025710, 28.96323675, 0.563579412, - 11.458509150, 27.70991388, -1.947354387, 8.933606299, 29.46637609, - -0.489290309, 7.918137207, 29.92393411}; - - int npc = 32; + double const coord[3 * 6] = {-2.75237178376284, 2.43247309226225, -0.01392519847964, + -0.93157260886974, 2.79621404458590, -0.01863384029005, + -3.43820531288547, 3.30583608421060, 1.42134539425148, + -2.43247309226225, -2.75237178376284, 0.01392519847964, + -2.79621404458590, -0.93157260886974, 0.01863384029005, + -3.30583608421060, -3.43820531288547, -1.42134539425148}; + + int npc = 6; // external point charges - double pc[32] = {0.431000, -0.664000, 0.173000, 0.020000, 0.020000, 0.020000, - 0.431000, -0.664000, 0.173000, 0.020000, 0.020000, 0.020000, - 0.431000, -0.664000, 0.173000, 0.020000, 0.020000, 0.020000, - 0.431000, -0.664000, 0.173000, 0.020000, 0.020000, 0.020000, - 0.431000, -0.664000, 0.173000, 0.020000, 0.020000, 0.020000, - 0.431000, -0.664000}; + double pc[6] = {-0.69645733, 0.36031084, 0.33614649, + -0.69645733, 0.36031084, 0.33614649}; // coordinates of point charges - double pccoord[3 * 32] = { - -1.696669514, 28.11745897, 55.50995136, -0.967547429, 28.88443423, - 54.00850230, -0.950868672, 31.58534217, 53.92332839, 0.439341438, - 32.17158366, 52.52085582, -0.346867500, 32.42710162, 55.70375073, - -2.886638164, 32.14280874, 53.49308978, 26.383000520, 21.74765157, - 24.17099786, 27.036716710, 22.30632379, 22.54790876, 26.261114830, - 20.54528062, 20.65041197, 27.011658890, 18.66770536, 21.04304684, - 24.203790340, 20.57032899, 20.55244459, 26.920273440, 21.10742805, - 18.78164663, 25.713072340, 18.66022959, 28.70604561, 26.111998120, - 18.26958272, 26.95615185, 27.664033370, 16.09612865, 26.54280365, - 29.614523860, 16.44665278, 27.10468085, 26.659753370, 14.55538454, - 27.47032282, 27.860595530, 15.78003352, 24.51695451, 12.692343200, - -11.99272547, 35.30000333, 11.005978380, -11.84086488, 36.01217832, - 9.538990432, -10.94052003, 33.92895508, 10.179104340, -11.53047861, - 32.06199633, 9.248569167, -8.901943552, 33.87774529, 7.688508235, - -11.83241211, 34.08195518, -11.936904730, 22.39399025, 52.97048306, - -10.467911780, 21.31627534, 53.20385904, -8.951231490, 21.02794412, - 50.98626089, -8.320084347, 22.92510643, 50.49100234, -10.224293460, - 20.31250935, 49.53331363, -7.443136243, 19.66216255, 51.30727310, - 15.373924870, 30.20492464, 53.87317117, 15.097532960, 31.53390939, - 52.63558511}; - - int numbers[32] = {7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, - 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7}; + double pccoord[3 * 6] = {2.75237178376284, -2.43247309226225, -0.01392519847964, + 0.93157260886974, -2.79621404458590, -0.01863384029005, + 3.43820531288547, -3.30583608421060, 1.42134539425148, + 2.43247309226225, 2.75237178376284, 0.01392519847964, + 2.79621404458590, 0.93157260886974, 0.01863384029005, + 3.30583608421060, 3.43820531288547, -1.42134539425148}; + + int numbers[6] = {8, 1, 1, 8, 1, 1}; // results will live in here double energy = 0.0; - double grad[3 * 31] = { - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; - double pcgrad[3 * 32] = { - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, - 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; + double grad[3 * 6] = {0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0}; + + double pcgrad[3 * 6] = {0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0}; xtb_TEnvironment env = xtb_newEnvironment(); xtb_TCalculator calc = xtb_newCalculator(); @@ -257,7 +219,7 @@ int testSecond() { if (xtb_checkEnvironment(env)) goto error; - xtb_loadGFN1xTB(env, mol, calc, NULL); + xtb_loadGFN2xTB(env, mol, calc, NULL); if (xtb_checkEnvironment(env)) goto error; @@ -290,17 +252,25 @@ int testSecond() { xtb_delete(mol); xtb_delete(env); - if (!check(energy, -36.669892958754, 1.0e-9, "Energy does not match")) + if (!check(energy, -10.160927754235, 1.0e-9, "Energy does not match")) goto error; - if (!check(grad[0], -0.001992161613, 1.0e-8, "grad[0] does not match")) + if (!check(grad[0 + 4 * 3], -2.0326750053811E-03, 1.0e-8, "grad[1,5] does not match")) + goto error; + if (!check(grad[1 + 1 * 3], 3.3125460324302E-03, 1.0e-8, "grad[2,2] does not match")) goto error; - if (!check(grad[92], 0.000556642950, 1.0e-8, "grad[92] does not match")) + if (!check(grad[0 + 3 * 3], -1.1929405354823E-03, 1.0e-8, "grad[1,4] does not match")) + goto error; + if (!check(grad[2 + 5 * 3], -1.6607683158587E-03, 1.0e-8, "grad[3,6] does not match")) goto error; - if (!check(pcgrad[0], 0.000007545625, 1.0e-8, "pcgrad[0] does not match")) + if (!check(pcgrad[0 + 4 * 3], -0.000248319582, 1.0e-8, "pcgrad[1,5] does not match")) + goto error; + if (!check(pcgrad[1 + 1 * 3], 0.001420844475, 1.0e-8, "pcgrad[2,2] does not match")) + goto error; + if (!check(pcgrad[0 + 3 * 3], 0.003746685275, 1.0e-8, "pcgrad[1,4] does not match")) goto error; - if (!check(pcgrad[95], 0.000013117708, 1.0e-8, "pcgrad[95] does not match")) + if (!check(pcgrad[2 + 5 * 3], 0.000651613444, 1.0e-8, "pcgrad[3,6] does not match")) goto error; return 0;