diff --git a/Makefile.am b/Makefile.am index 1ff940b..7f735dd 100644 --- a/Makefile.am +++ b/Makefile.am @@ -52,6 +52,7 @@ test_qmckl_f = tests/qmckl_f.f90 test_qmckl_fo = tests/qmckl_f.o src_qmckl_f = src/qmckl_f.f90 src_qmckl_fo = src/qmckl_f.o +header_tests = tests/chbrclf.h fortrandir = $(datadir)/$(PACKAGE_NAME)/fortran/ dist_fortran_DATA = $(qmckl_f) @@ -59,12 +60,11 @@ dist_fortran_DATA = $(qmckl_f) AM_CPPFLAGS = -I$(srcdir)/src -I$(srcdir)/include lib_LTLIBRARIES = src/libqmckl.la -src_libqmckl_la_SOURCES = $(qmckl_h) $(src_qmckl_f) $(C_FILES) $(F_FILES) $(H_PRIVATE_FUNC_FILES) $(H_PRIVATE_TYPE_FILES) +src_libqmckl_la_SOURCES = $(qmckl_h) $(src_qmckl_f) $(C_FILES) $(F_FILES) $(H_PRIVATE_FUNC_FILES) $(H_PRIVATE_TYPE_FILES) $(header_tests) export qmckl_f qmckl_h srcdir -CLEANFILES+=$(test_qmckl_f) $(src_qmckl_f) $(test_qmckl_o) $(src_qmckl_o) \ - $(qmckl_h) $(qmckl_f) +CLEANFILES+=$(test_qmckl_f) $(src_qmckl_f) $(test_qmckl_o) $(src_qmckl_o) htmlize_el=share/doc/qmckl/html/htmlize.el @@ -107,9 +107,9 @@ if QMCKL_DEVEL dist_src_DATA = $(ORG_FILES) $(TANGLED_FILES) $(EXPORTED_FILES) -BUILT_SOURCES = $(C_FILES) $(F_FILES) $(FH_FUNC_FILES) $(FH_TYPE_FILES) $(H_FUNC_FILES) $(H_TYPE_FILES) $(H_PRIVATE_FUNC_FILES) $(H_PRIVATE_TYPE_FILES) $(qmckl_f) $(qmckl_h) +BUILT_SOURCES = $(C_FILES) $(F_FILES) $(FH_FUNC_FILES) $(FH_TYPE_FILES) $(H_FUNC_FILES) $(H_TYPE_FILES) $(H_PRIVATE_FUNC_FILES) $(H_PRIVATE_TYPE_FILES) $(qmckl_f) $(qmckl_h) $(header_tests) -CLEANFILES += $(BUILT_SOURCES) $(C_TEST_FILES) $(F_TEST_FILES) $(TANGLED_FILES) $(C_TEST_FILES) $(F_TEST_FILES) $(qmckl_f) $(qmckl_h) $(HTML_FILES) $(TEXT_FILES) share/doc/qmckl/html/index.html $(EXPORTED_FILES) +CLEANFILES += $(BUILT_SOURCES) $(C_TEST_FILES) $(F_TEST_FILES) $(TANGLED_FILES) $(C_TEST_FILES) $(F_TEST_FILES) $(qmckl_f) $(qmckl_h) $(HTML_FILES) $(TEXT_FILES) share/doc/qmckl/html/index.html $(EXPORTED_FILES) $(header_tests) EXTRA_DIST += \ tools/build_doc.sh \ diff --git a/org/chbrclf.png b/org/chbrclf.png new file mode 100644 index 0000000..ab8d3fb Binary files /dev/null and b/org/chbrclf.png differ diff --git a/org/qmckl_distance.org b/org/qmckl_distance.org index e565c94..4ea0e2b 100644 --- a/org/qmckl_distance.org +++ b/org/qmckl_distance.org @@ -217,7 +217,7 @@ end function qmckl_distance_sq_f *** Performance - This function might be more efficient when ~A~ and ~B~ are + This function is more efficient when ~A~ and ~B~ are transposed. ** C interface :noexport: @@ -326,7 +326,7 @@ integer(qmckl_exit_code) function test_qmckl_distance_sq(context) bind(C) if (test_qmckl_distance_sq /= 0) return - test_qmckl_distance_sq = -1 + test_qmckl_distance_sq = QMCKL_FAILURE do j=1,n do i=1,m @@ -342,7 +342,7 @@ integer(qmckl_exit_code) function test_qmckl_distance_sq(context) bind(C) if (test_qmckl_distance_sq /= 0) return - test_qmckl_distance_sq = -1 + test_qmckl_distance_sq = QMCKL_FAILURE do j=1,n do i=1,m @@ -358,7 +358,7 @@ integer(qmckl_exit_code) function test_qmckl_distance_sq(context) bind(C) if (test_qmckl_distance_sq /= 0) return - test_qmckl_distance_sq = -1 + test_qmckl_distance_sq = QMCKL_FAILURE do j=1,n do i=1,m @@ -374,7 +374,7 @@ integer(qmckl_exit_code) function test_qmckl_distance_sq(context) bind(C) if (test_qmckl_distance_sq /= 0) return - test_qmckl_distance_sq = -1 + test_qmckl_distance_sq = QMCKL_FAILURE do j=1,n do i=1,m @@ -385,15 +385,15 @@ integer(qmckl_exit_code) function test_qmckl_distance_sq(context) bind(C) end do end do - test_qmckl_distance_sq = 0 + test_qmckl_distance_sq = QMCKL_SUCCESS deallocate(A,B,C) end function test_qmckl_distance_sq #+end_src #+begin_src c :comments link :tangle (eval c_test) -int test_qmckl_distance_sq(qmckl_context context); -assert(0 == test_qmckl_distance_sq(context)); +qmckl_exit_code test_qmckl_distance_sq(qmckl_context context); +assert(test_qmckl_distance_sq(context) == QMCKL_SUCCESS); #+end_src * Distance @@ -411,6 +411,9 @@ assert(0 == test_qmckl_distance_sq(context)); C_{ij} = \sqrt{\sum_{k=1}^3 (A_{k,i}-B_{k,j})^2} \] + If the input array is normal (~'N'~), the xyz coordinates are in + the leading dimension: ~[n][3]~ in C and ~(3,n)~ in Fortran. + #+NAME: qmckl_distance_args | qmckl_context | context | in | Global state | | char | transa | in | Array ~A~ is ~'N'~: Normal, ~'T'~: Transposed | @@ -595,8 +598,7 @@ end function qmckl_distance_f *** Performance - This function might be more efficient when ~A~ and ~B~ are - transposed. + This function is more efficient when ~A~ and ~B~ are transposed. ** C interface :noexport: @@ -704,7 +706,7 @@ integer(qmckl_exit_code) function test_qmckl_dist(context) bind(C) if (test_qmckl_dist /= 0) return - test_qmckl_dist = -1 + test_qmckl_dist = QMCKL_FAILURE do j=1,n do i=1,m @@ -720,7 +722,7 @@ integer(qmckl_exit_code) function test_qmckl_dist(context) bind(C) if (test_qmckl_dist /= 0) return - test_qmckl_dist = -1 + test_qmckl_dist = QMCKL_FAILURE do j=1,n do i=1,m @@ -736,7 +738,7 @@ integer(qmckl_exit_code) function test_qmckl_dist(context) bind(C) if (test_qmckl_dist /= 0) return - test_qmckl_dist = -1 + test_qmckl_dist = QMCKL_FAILURE do j=1,n do i=1,m @@ -752,7 +754,7 @@ integer(qmckl_exit_code) function test_qmckl_dist(context) bind(C) if (test_qmckl_dist /= 0) return - test_qmckl_dist = -1 + test_qmckl_dist = QMCKL_FAILURE do j=1,n do i=1,m @@ -763,15 +765,15 @@ integer(qmckl_exit_code) function test_qmckl_dist(context) bind(C) end do end do - test_qmckl_dist = 0 + test_qmckl_dist = QMCKL_SUCCESS deallocate(A,B,C) end function test_qmckl_dist #+end_src #+begin_src c :comments link :tangle (eval c_test) -int test_qmckl_dist(qmckl_context context); -assert(0 == test_qmckl_dist(context)); +qmckl_exit_code test_qmckl_dist(qmckl_context context); +assert(test_qmckl_dist(context) == QMCKL_SUCCESS); #+end_src * End of files :noexport: diff --git a/org/qmckl_electron.org b/org/qmckl_electron.org index dded8ff..b54aad6 100644 --- a/org/qmckl_electron.org +++ b/org/qmckl_electron.org @@ -25,6 +25,10 @@ up-spin and down-spin electrons, and the electron coordinates. #ifdef HAVE_CONFIG_H #include "config.h" #endif + +#include "CHBrClF.h" + + int main() { qmckl_context context; context = qmckl_context_create(); @@ -128,7 +132,7 @@ bool qmckl_electron_provided(const qmckl_context context) { the data has not been provided. If the function returns successfully, the variable pointed by the pointer given in argument contains the requested data. Otherwise, this variable is untouched. - + #+NAME:post #+begin_src c :exports none if ( (ctx->electron.uninitialized & mask) != 0) { @@ -137,7 +141,7 @@ if ( (ctx->electron.uninitialized & mask) != 0) { #+end_src *** Number of electrons - + #+begin_src c :comments org :tangle (eval h_func) :exports none qmckl_exit_code qmckl_get_electron_num (const qmckl_context context, int64_t* num); qmckl_exit_code qmckl_get_electron_up_num (const qmckl_context context, int64_t* up_num); @@ -211,10 +215,10 @@ qmckl_get_electron_down_num (const qmckl_context context, int64_t* down_num) { #+end_src *** Number of walkers - + A walker is a set of electron coordinates that are arguments of the wave function. ~walk_num~ is the number of walkers. - + #+begin_src c :comments org :tangle (eval h_func) :exports none qmckl_exit_code qmckl_get_electron_walk_num (const qmckl_context context, int64_t* walk_num); #+end_src @@ -242,12 +246,12 @@ qmckl_get_electron_walk_num (const qmckl_context context, int64_t* walk_num) { #+end_src *** Electron coordinates - + Returns the current electron coordinates. The pointer is assumed to point on a memory block of size ~3 * elec_num * walk_num~. In C the order of the indices is ~[walk_num][3][elec_num]~ and in Fortran it is ~(elec_num,3,walk_num)~. - + #+begin_src c :comments org :tangle (eval h_func) :exports none qmckl_exit_code qmckl_get_electron_coord (const qmckl_context context, double* coord); #+end_src @@ -273,7 +277,6 @@ qmckl_get_electron_coord (const qmckl_context context, double* elec_coord) { } #+end_src - ** Initialization functions @@ -451,25 +454,34 @@ qmckl_set_electron_coord(qmckl_context context, const double* coord) { ** Test - #+begin_src c :tangle (eval c_test) + #+begin_src python :results output :exports none :tangle none +import numpy as np + + #+end_src + + #+begin_src c :tangle (eval c_test) /* Reference input data */ +const int64_t walk_num = chbrclf_walk_num; +const int64_t elec_num = chbrclf_elec_num; +const int64_t elec_up_num = chbrclf_elec_up_num; +const int64_t elec_dn_num = chbrclf_elec_dn_num; +const double*** elec_coord = chbrclf_elec_coord; -#define up_num ((int64_t) 3) -#define down_num ((int64_t) 2) -#define walk_num ((int64_t) 2) -#define num (up_num+down_num) +const int64_t nucl_num = chbrclf_nucl_num; +const double* charge = chbrclf_charge; +const double** nucl_coord = chbrclf_nucl_coord; -double coord[walk_num*3*num] = - { 7.303633091022677881e+00, 1.375868694453235719e+01, 1.167371490471771217e-01, - 4.547755371567960836e+00, 3.245907105524011182e+00, 2.410764357550297110e-01, - 5.932816068137344523e+00, 1.491671465549257469e+01, 3.825374039119375236e-01, - 7.347336142660052083e+00, 1.341946976062362129e+00, 1.648917914228352322e+00, - 5.735221530102248444e+00, 1.064667491680036271e+01, 4.227201772236627297e-01, - 8.099550978782254163e+00, 6.861498941099086757e+00, 4.015884841159429036e-02, - 1.014757367558326173e+01, 5.219335322173662917e+00, 5.037004126899931322e-02, - 1.484094322159507051e+01, 9.777903829455864226e+00, 5.243007994024882767e-02, - 9.081723054990456845e+00, 5.499568496038920173e+00, 2.910446438899221347e-02, - 2.583154239492383653e+00, 1.442282811294904432e+00, 6.387191629878670451e-02 }; +double* coord = (double*) malloc(walk_num*num*3*sizeof(double)); + +double* x = coord; +for (int i=0 ; inucleus.uninitialized & mask) != 0) { #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code -qmckl_get_nucleus_num (const qmckl_context context, int64_t* num) { +qmckl_get_nucleus_num (const qmckl_context context, int64_t* const num) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_INVALID_CONTEXT; } + if (num == NULL) { + return qmckl_failwith( context, + QMCKL_INVALID_ARG_2, + "qmckl_get_nucleus_num", + "num is a null pointer"); + } + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); int32_t mask = 1 << 0; if ( (ctx->nucleus.uninitialized & mask) != 0) { - return QMCKL_NOT_PROVIDED; + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_get_nucleus_num", + "nucleus data is not provided"); } assert (ctx->nucleus.num >= (int64_t) 0); @@ -130,19 +143,29 @@ qmckl_get_nucleus_num (const qmckl_context context, int64_t* num) { qmckl_exit_code -qmckl_get_nucleus_charge (const qmckl_context context, double* charge) { +qmckl_get_nucleus_charge (const qmckl_context context, double* const charge) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_INVALID_CONTEXT; } + if (charge == NULL) { + return qmckl_failwith( context, + QMCKL_INVALID_ARG_2, + "qmckl_get_nucleus_charge", + "charge is a null pointer"); + } + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); int32_t mask = 1 << 1; if ( (ctx->nucleus.uninitialized & mask) != 0) { - return QMCKL_NOT_PROVIDED; + return qmckl_failwith( context, + QMCKL_NOT_PROVIDED, + "qmckl_get_nucleus_charge", + "nucleus data is not provided"); } assert (ctx->nucleus.charge != NULL); @@ -151,24 +174,37 @@ qmckl_get_nucleus_charge (const qmckl_context context, double* charge) { qmckl_exit_code rc; rc = qmckl_get_nucleus_num(context, &nucl_num); if (rc != QMCKL_SUCCESS) return rc; - - double* result = memcpy(charge, ctx->nucleus.charge, nucl_num*sizeof(double)); - if (result == NULL) return QMCKL_FAILURE; - + + memcpy(charge, ctx->nucleus.charge, nucl_num*sizeof(double)); + return QMCKL_SUCCESS; } qmckl_exit_code -qmckl_get_nucleus_coord (const qmckl_context context, double* coord) { - +qmckl_get_nucleus_coord (const qmckl_context context, const char transp, double* const coord) { + if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { return QMCKL_INVALID_CONTEXT; } + if (transp != 'N' && transp != 'T') { + return qmckl_failwith( context, + QMCKL_INVALID_ARG_2, + "qmckl_get_nucleus_coord", + "transp should be 'N' or 'T'"); + } + + if (coord == NULL) { + return qmckl_failwith( context, + QMCKL_INVALID_ARG_3, + "qmckl_get_nucleus_coord", + "coord is a null pointer"); + } + qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); - + int32_t mask = 1 << 2; if ( (ctx->nucleus.uninitialized & mask) != 0) { @@ -179,12 +215,18 @@ qmckl_get_nucleus_coord (const qmckl_context context, double* coord) { qmckl_exit_code rc; rc = qmckl_get_nucleus_num(context, &nucl_num); if (rc != QMCKL_SUCCESS) return rc; - + assert (ctx->nucleus.coord != NULL); - double* result = memcpy(coord, ctx->nucleus.coord, 3*nucl_num*sizeof(double)); - if (result == NULL) return QMCKL_FAILURE; - + if (transp == 'N') { + rc = qmckl_transpose(context, nucl_num, 3, + ctx->nucleus.coord, nucl_num, + coord, 3); + if (rc != QMCKL_SUCCESS) return rc; + } else { + memcpy(coord, ctx->nucleus.coord, 3*nucl_num*sizeof(double)); + } + return QMCKL_SUCCESS; } #+end_src @@ -218,7 +260,7 @@ bool qmckl_nucleus_provided(const qmckl_context context) { #+begin_src c :comments org :tangle (eval h_func) qmckl_exit_code qmckl_set_nucleus_num (qmckl_context context, const int64_t num); qmckl_exit_code qmckl_set_nucleus_charge (qmckl_context context, const double* charge); -qmckl_exit_code qmckl_set_nucleus_coord (qmckl_context context, const double* coord); +qmckl_exit_code qmckl_set_nucleus_coord (qmckl_context context, const char transp, const double* coord); #+end_src #+NAME:pre2 @@ -240,7 +282,7 @@ return QMCKL_SUCCESS; To set the number of nuclei, use - + #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_set_nucleus_num(qmckl_context context, const int64_t num) { @@ -261,13 +303,20 @@ qmckl_set_nucleus_num(qmckl_context context, const int64_t num) { } #+end_src - The following function sets the nuclear charges of all the atoms. + The following function sets the nuclear charges of all the atoms. #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code qmckl_set_nucleus_charge(qmckl_context context, const double* charge) { <> + if (charge == NULL) { + return qmckl_failwith( context, + QMCKL_INVALID_ARG_2, + "qmckl_set_nucleus_charge", + "charge is a null pointer"); + } + int64_t num; qmckl_exit_code rc; @@ -275,7 +324,7 @@ qmckl_set_nucleus_charge(qmckl_context context, const double* charge) { rc = qmckl_get_nucleus_num(context, &num); if (rc != QMCKL_SUCCESS) return rc; - + if (ctx->nucleus.charge != NULL) { qmckl_free(context, ctx->nucleus.charge); ctx->nucleus.charge= NULL; @@ -294,7 +343,7 @@ qmckl_set_nucleus_charge(qmckl_context context, const double* charge) { } ctx->nucleus.charge= memcpy(ctx->nucleus.charge, charge, num*sizeof(double)); assert (ctx->nucleus.charge != NULL); - + <> } #+end_src @@ -304,24 +353,24 @@ qmckl_set_nucleus_charge(qmckl_context context, const double* charge) { #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none qmckl_exit_code -qmckl_set_nucleus_coord(qmckl_context context, const double* coord) { +qmckl_set_nucleus_coord(qmckl_context context, const char transp, const double* coord) { <> - int64_t num; + int64_t nucl_num; qmckl_exit_code rc; int32_t mask = 1 << 2; - rc = qmckl_get_nucleus_num(context, &num); + rc = qmckl_get_nucleus_num(context, &nucl_num); if (rc != QMCKL_SUCCESS) return rc; - + if (ctx->nucleus.coord != NULL) { qmckl_free(context, ctx->nucleus.coord); ctx->nucleus.coord = NULL; } qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = 3*num*sizeof(double); + mem_info.size = 3*nucl_num*sizeof(double); assert(ctx->nucleus.coord == NULL); @@ -332,9 +381,15 @@ qmckl_set_nucleus_coord(qmckl_context context, const double* coord) { "qmckl_set_nucleus_coord", NULL); } - ctx->nucleus.coord = memcpy(ctx->nucleus.coord, coord, 3*num*sizeof(double)); - assert (ctx->nucleus.coord != NULL); - + if (transp == 'N') { + rc = qmckl_transpose(context, 3, nucl_num, + coord, 3, + ctx->nucleus.coord, nucl_num); + if (rc != QMCKL_SUCCESS) return rc; + } else { + memcpy(ctx->nucleus.coord, coord, 3*nucl_num*sizeof(double)); + } + <> } #+end_src @@ -342,23 +397,10 @@ qmckl_set_nucleus_coord(qmckl_context context, const double* coord) { ** Test #+begin_src c :tangle (eval c_test) -/* Reference input data */ +const int64_t nucl_num = chbrclf_nucl_num; +const double* nucl_charge = chbrclf_charge; +const double* nucl_coord = &(chbrclf_nucl_coord[0][0]); -#define num ((int64_t) 9) - -double charge[num] = { 6., 6., 6., 7., 7., 1., 1., 1., 1. }; - -double coord[3*num] = - { 4.166279566732572e-01, -1.526183863767697e+00, 1.041604719335635e+00, - -1.903457631371503e+00, 2.242154435363994e+00, 6.550163404813796e-01, - -3.575005445908036e+00, -3.063638942318878e+00, 2.086739409279095e+00, - 2.060062599100338e+00, -1.623431626827498e+00, -1.930074272670425e+00, - 9.491495662916423e-01, 3.808343139803397e-01, 4.077482772289367e+00, - 1.841031662652821e+00, -2.945591662994877e+00, -3.670011011125464e+00, - 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00, - 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00, - 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00}; - /* --- */ qmckl_exit_code rc; @@ -370,42 +412,50 @@ rc = qmckl_get_nucleus_num (context, &n); assert(rc == QMCKL_NOT_PROVIDED); -rc = qmckl_set_nucleus_num (context, num); +rc = qmckl_set_nucleus_num (context, nucl_num); assert(rc == QMCKL_SUCCESS); assert(!qmckl_nucleus_provided(context)); rc = qmckl_get_nucleus_num (context, &n); assert(rc == QMCKL_SUCCESS); -assert(n == num); +assert(n == nucl_num); -double coord2[3*num]; +double nucl_coord2[3*nucl_num]; -rc = qmckl_get_nucleus_coord (context, coord2); +rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2); assert(rc == QMCKL_NOT_PROVIDED); -rc = qmckl_set_nucleus_coord (context, coord); +rc = qmckl_set_nucleus_coord (context, 'T', &(nucl_coord[0])); assert(rc == QMCKL_SUCCESS); -rc = qmckl_get_nucleus_coord (context, coord2); -assert(rc == QMCKL_SUCCESS); -for (size_t i=0 ; i<3*num ; ++i) { - assert( coord[i] == coord2[i] ); - } - assert(!qmckl_nucleus_provided(context)); -double charge2[num]; +rc = qmckl_get_nucleus_coord (context, 'N', nucl_coord2); +assert(rc == QMCKL_SUCCESS); +for (size_t k=0 ; k<3 ; ++k) { + for (size_t i=0 ; inucleus.provided) return QMCKL_NOT_PROVIDED; - + /* Allocate array */ if (ctx->nucleus.nn_distance == NULL) { - + qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; mem_info.size = ctx->nucleus.num * ctx->nucleus.num * sizeof(double); double* nn_distance = (double*) qmckl_malloc(context, mem_info); - + if (nn_distance == NULL) { return qmckl_failwith( context, QMCKL_ALLOCATION_FAILED, @@ -547,7 +597,7 @@ qmckl_exit_code qmckl_compute_nn_distance ( double* const nn_distance ); #+end_src - + #+CALL: generate_c_interface(table=qmckl_nn_distance_args,rettyp="qmckl_exit_code",fname="qmckl_compute_nn_distance") #+RESULTS: @@ -578,12 +628,11 @@ qmckl_exit_code qmckl_compute_nn_distance ( assert(qmckl_nucleus_provided(context)); -double distance[num*num]; -rc = qmckl_get_nucleus_nn_distance(context, distance); +double distance[nucl_num*nucl_num]; rc = qmckl_get_nucleus_nn_distance(context, distance); assert(distance[0] == 0.); -assert(distance[1] == distance[num]); -assert(fabs(distance[1]-4.164450441785663) < 1.e-12); +assert(distance[1] == distance[nucl_num]); +assert(fabs(distance[1]-2.070304721365169) < 1.e-12); #+end_src @@ -642,7 +691,7 @@ qmckl_exit_code qmckl_provide_nucleus_repulsion(qmckl_context context) rc = qmckl_provide_nn_distance(context); if (rc != QMCKL_SUCCESS) return rc; - + rc = qmckl_compute_nucleus_repulsion(context, ctx->nucleus.num, ctx->nucleus.charge, @@ -651,7 +700,7 @@ qmckl_exit_code qmckl_provide_nucleus_repulsion(qmckl_context context) if (rc != QMCKL_SUCCESS) { return rc; } - + ctx->nucleus.repulsion_date = ctx->date; return QMCKL_SUCCESS; @@ -695,10 +744,9 @@ integer function qmckl_compute_nucleus_repulsion_f(context, nucl_num, charge, nn energy = 0.d0 do j=2, nucl_num do i=1, j-1 - energy = energy + charge(i) * charge(j) / nn_distance(i,j) + energy = energy + charge(i) * charge(j) / nn_distance(i,j) end do end do - print *, energy end function qmckl_compute_nucleus_repulsion_f #+end_src @@ -746,7 +794,7 @@ assert(qmckl_nucleus_provided(context)); double rep; rc = qmckl_get_nucleus_repulsion(context, &rep); -assert(rep - 163.50434957121263 < 1.e-10); +assert(rep - 318.2309879436158 < 1.e-10); #+end_src @@ -787,8 +835,6 @@ assert(rep - 163.50434957121263 < 1.e-10); #+end_src #+RESULTS: - | | color | - | | listings | # -*- mode: org -*- diff --git a/org/qmckl_tests.org b/org/qmckl_tests.org new file mode 100644 index 0000000..c9e48d3 --- /dev/null +++ b/org/qmckl_tests.org @@ -0,0 +1,677 @@ +#+TITLE: Data for Tests + +# -*- org-image-actual-width: 300 -*- + +* CHBrClF + + This test is the all-electron Hartree-Fock wave function of CHClBr, + in the aug-cc-pVTZ basis set. This is a non-symmetric molecule made + of 5 atoms, heavy and light one. The aug-cc-pVTZ basis set has both + diffuse and compact Gaussians, with large and small contractions, + and with a high maximum angular momentum. + + [[./chbrclf.png]] + + | Number of atoms | 5 | + | Number of alpha electrons | 34 | + | Number of beta electrons | 34 | + | Max number of primitives | 15 | + | Highest angular momentum | F | + | Atomic basis set | aug-cc-pVTZ | + | Nuclear repulsion energy | 318.2309879436158 | + | Number of primitives | 502 | + | Number of cartesian basis functions | 263 | + | Number of molecular orbitals | 224 | + | Hartree-Fock energy | -3169.90467157 Ha | + + +** XYZ coordinates + +#+BEGIN_example + 5 +CHBrClF + C 0.580107 0.471341 0.411546 + H 0.618322 0.595674 1.499355 + F 0.786938 1.650849 -0.204021 +Cl 1.850884 -0.689476 -0.067323 +Br -1.218470 -0.187436 -0.028227 +#+END_example + + Nuclear coordinates are stored in atomic units in transposed format. + +#+begin_src c :tangle ../tests/chbrclf.h +#define chbrclf_nucl_num ((int64_t) 5) + +const double chbrclf_charge[chbrclf_nucl_num] = { 6., 1., 9., 17., 35. }; + +const double chbrclf_nucl_coord[3][chbrclf_nucl_num] = +{ { 1.096243353458458e+00, 1.168459237342663e+00, 1.487097297712132e+00, 3.497663849983889e+00, -2.302574592081335e+00 }, + { 8.907054016973815e-01, 1.125660720053393e+00, 3.119652484478797e+00, -1.302920810073182e+00, -3.542027060505035e-01 }, + { 7.777092280258892e-01, 2.833370314829343e+00, -3.855438138411500e-01, -1.272220319439064e-01, -5.334129934317614e-02 } }; +#+end_src + +** Atomic basis set + +#+BEGIN_example +HYDROGEN +S 5 +1 3.387000E+01 6.068000E-03 +2 5.095000E+00 4.530800E-02 +3 1.159000E+00 2.028220E-01 +4 3.258000E-01 5.039030E-01 +5 1.027000E-01 3.834210E-01 +S 1 +1 3.258000E-01 1.000000E+00 +S 1 +1 1.027000E-01 1.000000E+00 +S 1 +1 0.0252600 1.0000000 +P 1 +1 1.407000E+00 1.000000E+00 +P 1 +1 3.880000E-01 1.000000E+00 +P 1 +1 0.1020000 1.0000000 +D 1 +1 1.057000E+00 1.0000000 +D 1 +1 0.2470000 1.0000000 + +CARBON +S 10 +1 8.236000E+03 5.310000E-04 +2 1.235000E+03 4.108000E-03 +3 2.808000E+02 2.108700E-02 +4 7.927000E+01 8.185300E-02 +5 2.559000E+01 2.348170E-01 +6 8.997000E+00 4.344010E-01 +7 3.319000E+00 3.461290E-01 +8 9.059000E-01 3.937800E-02 +9 3.643000E-01 -8.983000E-03 +10 1.285000E-01 2.385000E-03 +S 10 +1 8.236000E+03 -1.130000E-04 +2 1.235000E+03 -8.780000E-04 +3 2.808000E+02 -4.540000E-03 +4 7.927000E+01 -1.813300E-02 +5 2.559000E+01 -5.576000E-02 +6 8.997000E+00 -1.268950E-01 +7 3.319000E+00 -1.703520E-01 +8 9.059000E-01 1.403820E-01 +9 3.643000E-01 5.986840E-01 +10 1.285000E-01 3.953890E-01 +S 1 +1 9.059000E-01 1.000000E+00 +S 1 +1 1.285000E-01 1.000000E+00 +S 1 +1 0.0440200 1.0000000 +P 5 +1 1.871000E+01 1.403100E-02 +2 4.133000E+00 8.686600E-02 +3 1.200000E+00 2.902160E-01 +4 3.827000E-01 5.010080E-01 +5 1.209000E-01 3.434060E-01 +P 1 +1 3.827000E-01 1.000000E+00 +P 1 +1 1.209000E-01 1.000000E+00 +P 1 +1 0.0356900 1.0000000 +D 1 +1 1.097000E+00 1.000000E+00 +D 1 +1 3.180000E-01 1.000000E+00 +D 1 +1 0.1000000 1.0000000 +F 1 +1 7.610000E-01 1.0000000 +F 1 +1 0.2680000 1.0000000 + +FLUORINE +S 10 +1 1.950000E+04 5.070000E-04 +2 2.923000E+03 3.923000E-03 +3 6.645000E+02 2.020000E-02 +4 1.875000E+02 7.901000E-02 +5 6.062000E+01 2.304390E-01 +6 2.142000E+01 4.328720E-01 +7 7.950000E+00 3.499640E-01 +8 2.257000E+00 4.323300E-02 +9 8.815000E-01 -7.892000E-03 +10 3.041000E-01 2.384000E-03 +S 10 +1 1.950000E+04 -1.170000E-04 +2 2.923000E+03 -9.120000E-04 +3 6.645000E+02 -4.717000E-03 +4 1.875000E+02 -1.908600E-02 +5 6.062000E+01 -5.965500E-02 +6 2.142000E+01 -1.400100E-01 +7 7.950000E+00 -1.767820E-01 +8 2.257000E+00 1.716250E-01 +9 8.815000E-01 6.050430E-01 +10 3.041000E-01 3.695120E-01 +S 1 +1 2.257000E+00 1.000000E+00 +S 1 +1 3.041000E-01 1.000000E+00 +S 1 +1 0.0915800 1.0000000 +P 5 +1 4.388000E+01 1.666500E-02 +2 9.926000E+00 1.044720E-01 +3 2.930000E+00 3.172600E-01 +4 9.132000E-01 4.873430E-01 +5 2.672000E-01 3.346040E-01 +P 1 +1 9.132000E-01 1.000000E+00 +P 1 +1 2.672000E-01 1.000000E+00 +P 1 +1 0.0736100 1.0000000 +D 1 +1 3.107000E+00 1.000000E+00 +D 1 +1 8.550000E-01 1.000000E+00 +D 1 +1 0.2920000 1.0000000 +F 1 +1 1.917000E+00 1.0000000 +F 1 +1 0.7240000 1.0000000 + +S 20 +1 1.063900E+07 7.000000E-07 +2 1.593400E+06 5.700000E-06 +3 3.626100E+05 3.030000E-05 +4 1.027000E+05 1.275000E-04 +5 3.350100E+04 4.659000E-04 +6 1.209300E+04 1.509600E-03 +7 4.715900E+03 4.485200E-03 +8 1.955600E+03 1.198350E-02 +9 8.526100E+02 2.895710E-02 +10 3.876700E+02 5.815660E-02 +11 1.826800E+02 8.881330E-02 +12 8.824500E+01 4.452440E-02 +13 3.926300E+01 -2.060387E-01 +14 1.923400E+01 -5.127017E-01 +15 9.405700E+00 -1.509349E-01 +16 4.160100E+00 6.789203E-01 +17 1.899500E+00 5.817697E-01 +18 6.047200E-01 4.675550E-02 +19 3.011400E-01 -1.118250E-02 +20 1.251500E-01 2.440200E-03 +S 20 +1 1.063900E+07 -2.000000E-07 +2 1.593400E+06 -1.800000E-06 +3 3.626100E+05 -9.300000E-06 +4 1.027000E+05 -3.910000E-05 +5 3.350100E+04 -1.428000E-04 +6 1.209300E+04 -4.628000E-04 +7 4.715900E+03 -1.375000E-03 +8 1.955600E+03 -3.678400E-03 +9 8.526100E+02 -8.898100E-03 +10 3.876700E+02 -1.795290E-02 +11 1.826800E+02 -2.757320E-02 +12 8.824500E+01 -1.409530E-02 +13 3.926300E+01 6.725610E-02 +14 1.923400E+01 1.766928E-01 +15 9.405700E+00 5.288610E-02 +16 4.160100E+00 -3.075955E-01 +17 1.899500E+00 -4.700658E-01 +18 6.047200E-01 2.558761E-01 +19 3.011400E-01 6.980341E-01 +20 1.251500E-01 2.967256E-01 +S 1 +1 6.047200E-01 1.000000E+00 +S 1 +1 1.251500E-01 1.000000E+00 +S 1 +1 0.0455930 1.0000000 +P 13 +1 8.676500E+03 4.357000E-04 +2 2.055900E+03 3.781500E-03 +3 6.662300E+02 2.047820E-02 +4 2.531000E+02 7.928340E-02 +5 1.061200E+02 2.178473E-01 +6 4.724200E+01 3.878585E-01 +7 2.182500E+01 3.594350E-01 +8 9.968400E+00 1.121995E-01 +9 4.517100E+00 4.387400E-03 +10 1.998200E+00 1.780900E-03 +11 7.098800E-01 -4.576000E-04 +12 2.814500E-01 2.122000E-04 +13 1.020400E-01 -7.340000E-05 +P 9 +1 6.633000E+02 -6.521450E-04 +2 1.568000E+02 -5.194450E-03 +3 4.998000E+01 -2.469380E-02 +4 1.842000E+01 -7.281670E-02 +5 7.240000E+00 -1.340300E-01 +6 2.922000E+00 -9.477420E-02 +7 1.022000E+00 2.622890E-01 +8 3.818000E-01 5.646670E-01 +9 1.301000E-01 3.412500E-01 +P 1 +1 1.022000E+00 1.000000E+00 +P 1 +1 1.301000E-01 1.000000E+00 +P 1 +1 0.0419000 1.0000000 +D 1 +1 1.046000E+00 1.000000E+00 +D 1 +1 3.440000E-01 1.000000E+00 +D 1 +1 0.1350000 1.0000000 +F 1 +1 7.060000E-01 1.0000000 +F 1 +1 0.3120000 1.0000000 + +CHLORINE +S 15 +1 4.561000E+05 4.929700E-05 +2 6.833000E+04 3.830290E-04 +3 1.555000E+04 2.008540E-03 +4 4.405000E+03 8.385580E-03 +5 1.439000E+03 2.947030E-02 +6 5.204000E+02 8.783250E-02 +7 2.031000E+02 2.114730E-01 +8 8.396000E+01 3.653640E-01 +9 3.620000E+01 3.408840E-01 +10 1.583000E+01 1.021330E-01 +11 6.334000E+00 3.116750E-03 +12 2.694000E+00 1.057510E-03 +13 9.768000E-01 -3.780000E-04 +14 4.313000E-01 1.561360E-04 +15 1.625000E-01 -5.141260E-05 +S 15 +1 4.561000E+05 -1.383040E-05 +2 6.833000E+04 -1.072790E-04 +3 1.555000E+04 -5.650830E-04 +4 4.405000E+03 -2.361350E-03 +5 1.439000E+03 -8.458860E-03 +6 5.204000E+02 -2.596380E-02 +7 2.031000E+02 -6.863620E-02 +8 8.396000E+01 -1.418740E-01 +9 3.620000E+01 -1.993190E-01 +10 1.583000E+01 -1.956620E-02 +11 6.334000E+00 4.997410E-01 +12 2.694000E+00 5.637360E-01 +13 9.768000E-01 7.903250E-02 +14 4.313000E-01 -8.350910E-03 +15 1.625000E-01 2.324560E-03 +S 15 +1 4.561000E+05 4.185460E-06 +2 6.833000E+04 3.243950E-05 +3 1.555000E+04 1.711050E-04 +4 4.405000E+03 7.141760E-04 +5 1.439000E+03 2.567050E-03 +6 5.204000E+02 7.885520E-03 +7 2.031000E+02 2.108670E-02 +8 8.396000E+01 4.422640E-02 +9 3.620000E+01 6.516700E-02 +10 1.583000E+01 6.030120E-03 +11 6.334000E+00 -2.064950E-01 +12 2.694000E+00 -4.058710E-01 +13 9.768000E-01 7.595580E-02 +14 4.313000E-01 7.256610E-01 +15 1.625000E-01 3.944230E-01 +S 1 +1 9.768000E-01 1.000000E+00 +S 1 +1 1.625000E-01 1.000000E+00 +S 1 +1 0.0591000 1.0000000 +P 9 +1 6.633000E+02 2.404480E-03 +2 1.568000E+02 1.921480E-02 +3 4.998000E+01 8.850970E-02 +4 1.842000E+01 2.560200E-01 +5 7.240000E+00 4.369270E-01 +6 2.922000E+00 3.503340E-01 +7 1.022000E+00 5.854950E-02 +8 3.818000E-01 -4.584230E-03 +9 1.301000E-01 2.269700E-03 +P 9 +1 6.633000E+02 -6.521450E-04 +2 1.568000E+02 -5.194450E-03 +3 4.998000E+01 -2.469380E-02 +4 1.842000E+01 -7.281670E-02 +5 7.240000E+00 -1.340300E-01 +6 2.922000E+00 -9.477420E-02 +7 1.022000E+00 2.622890E-01 +8 3.818000E-01 5.646670E-01 +9 1.301000E-01 3.412500E-01 +P 1 +1 1.022000E+00 1.000000E+00 +P 1 +1 1.301000E-01 1.000000E+00 +P 1 +1 0.0419000 1.0000000 +D 1 +1 1.046000E+00 1.000000E+00 +D 1 +1 3.440000E-01 1.000000E+00 +D 1 +1 0.1350000 1.0000000 +F 1 +1 7.060000E-01 1.0000000 +F 1 +1 0.3120000 1.0000000 + +BROMINE +S 20 +1 1.063900E+07 5.900000E-06 +2 1.593400E+06 4.610000E-05 +3 3.626100E+05 2.422000E-04 +4 1.027000E+05 1.022600E-03 +5 3.350100E+04 3.711300E-03 +6 1.209300E+04 1.197850E-02 +7 4.715900E+03 3.469270E-02 +8 1.955600E+03 8.912390E-02 +9 8.526100E+02 1.934557E-01 +10 3.876700E+02 3.209019E-01 +11 1.826800E+02 3.299233E-01 +12 8.824500E+01 1.494121E-01 +13 3.926300E+01 1.499380E-02 +14 1.923400E+01 -9.165000E-04 +15 9.405700E+00 4.380000E-04 +16 4.160100E+00 -2.398000E-04 +17 1.899500E+00 7.360000E-05 +18 6.047200E-01 -3.670000E-05 +19 3.011400E-01 2.390000E-05 +20 1.251500E-01 -5.600000E-06 +S 20 +1 1.063900E+07 -1.900000E-06 +2 1.593400E+06 -1.450000E-05 +3 3.626100E+05 -7.610000E-05 +4 1.027000E+05 -3.210000E-04 +5 3.350100E+04 -1.170900E-03 +6 1.209300E+04 -3.796800E-03 +7 4.715900E+03 -1.123070E-02 +8 1.955600E+03 -2.992770E-02 +9 8.526100E+02 -7.127060E-02 +10 3.876700E+02 -1.403136E-01 +11 1.826800E+02 -2.030763E-01 +12 8.824500E+01 -9.609850E-02 +13 3.926300E+01 3.558086E-01 +14 1.923400E+01 5.921792E-01 +15 9.405700E+00 2.215977E-01 +16 4.160100E+00 1.376480E-02 +17 1.899500E+00 8.395000E-04 +18 6.047200E-01 -4.510000E-05 +19 3.011400E-01 -8.500000E-06 +20 1.251500E-01 -1.240000E-05 +S 20 +1 1.063900E+07 7.000000E-07 +2 1.593400E+06 5.700000E-06 +3 3.626100E+05 3.030000E-05 +4 1.027000E+05 1.275000E-04 +5 3.350100E+04 4.659000E-04 +6 1.209300E+04 1.509600E-03 +7 4.715900E+03 4.485200E-03 +8 1.955600E+03 1.198350E-02 +9 8.526100E+02 2.895710E-02 +10 3.876700E+02 5.815660E-02 +11 1.826800E+02 8.881330E-02 +12 8.824500E+01 4.452440E-02 +13 3.926300E+01 -2.060387E-01 +14 1.923400E+01 -5.127017E-01 +15 9.405700E+00 -1.509349E-01 +16 4.160100E+00 6.789203E-01 +17 1.899500E+00 5.817697E-01 +18 6.047200E-01 4.675550E-02 +19 3.011400E-01 -1.118250E-02 +20 1.251500E-01 2.440200E-03 +S 20 +1 1.063900E+07 -2.000000E-07 +2 1.593400E+06 -1.800000E-06 +3 3.626100E+05 -9.300000E-06 +4 1.027000E+05 -3.910000E-05 +5 3.350100E+04 -1.428000E-04 +6 1.209300E+04 -4.628000E-04 +7 4.715900E+03 -1.375000E-03 +8 1.955600E+03 -3.678400E-03 +9 8.526100E+02 -8.898100E-03 +10 3.876700E+02 -1.795290E-02 +11 1.826800E+02 -2.757320E-02 +12 8.824500E+01 -1.409530E-02 +13 3.926300E+01 6.725610E-02 +14 1.923400E+01 1.766928E-01 +15 9.405700E+00 5.288610E-02 +16 4.160100E+00 -3.075955E-01 +17 1.899500E+00 -4.700658E-01 +18 6.047200E-01 2.558761E-01 +19 3.011400E-01 6.980341E-01 +20 1.251500E-01 2.967256E-01 +S 1 +1 6.047200E-01 1.000000E+00 +S 1 +1 1.251500E-01 1.000000E+00 +S 1 +1 0.0455930 1.0000000 +P 13 +1 8.676500E+03 4.357000E-04 +2 2.055900E+03 3.781500E-03 +3 6.662300E+02 2.047820E-02 +4 2.531000E+02 7.928340E-02 +5 1.061200E+02 2.178473E-01 +6 4.724200E+01 3.878585E-01 +7 2.182500E+01 3.594350E-01 +8 9.968400E+00 1.121995E-01 +9 4.517100E+00 4.387400E-03 +10 1.998200E+00 1.780900E-03 +11 7.098800E-01 -4.576000E-04 +12 2.814500E-01 2.122000E-04 +13 1.020400E-01 -7.340000E-05 +P 13 +1 8.676500E+03 -1.748000E-04 +2 2.055900E+03 -1.526300E-03 +3 6.662300E+02 -8.339900E-03 +4 2.531000E+02 -3.322030E-02 +5 1.061200E+02 -9.541800E-02 +6 4.724200E+01 -1.824026E-01 +7 2.182500E+01 -1.558308E-01 +8 9.968400E+00 1.867899E-01 +9 4.517100E+00 5.427733E-01 +10 1.998200E+00 3.873309E-01 +11 7.098800E-01 4.530690E-02 +12 2.814500E-01 -4.378400E-03 +13 1.020400E-01 1.811100E-03 +P 13 +1 8.676500E+03 4.510000E-05 +2 2.055900E+03 3.964000E-04 +3 6.662300E+02 2.155500E-03 +4 2.531000E+02 8.672000E-03 +5 1.061200E+02 2.486800E-02 +6 4.724200E+01 4.854720E-02 +7 2.182500E+01 3.961560E-02 +8 9.968400E+00 -6.057490E-02 +9 4.517100E+00 -1.871699E-01 +10 1.998200E+00 -1.377757E-01 +11 7.098800E-01 2.928021E-01 +12 2.814500E-01 5.760896E-01 +13 1.020400E-01 3.078617E-01 +P 1 +1 7.098800E-01 1.000000E+00 +P 1 +1 1.020400E-01 1.000000E+00 +P 1 +1 0.0351420 1.0000000 +D 8 +1 4.038300E+02 1.473200E-03 +2 1.211700E+02 1.267250E-02 +3 4.634500E+01 5.804510E-02 +4 1.972100E+01 1.705103E-01 +5 8.862400E+00 3.185958E-01 +6 3.996200E+00 3.845023E-01 +7 1.763600E+00 2.737737E-01 +8 7.061900E-01 7.439670E-02 +D 1 +1 7.061900E-01 1.000000E+00 +D 1 +1 2.639000E-01 1.000000E+00 +D 1 +1 0.1047000 1.0000000 +F 1 +1 5.515000E-01 1.0000000 +F 1 +1 0.2580000 1.0000000 + + #+END_example + +** Electron coordinates + + Electron coordinates are stored in atomic units in normal format. + +#+begin_src c :tangle ../tests/chbrclf.h +#define chbrclf_elec_up_num ((int64_t) 34) +#define chbrclf_elec_dn_num ((int64_t) 34) +#define chbrclf_elec_num ((int64_t) 68) +#define chbrclf_walk_num ((int64_t) 2) + +const double chbrclf_elec_coord[chbrclf_walk_num][chbrclf_elec_num][3] = { { + {-2.26995253563, -5.15737533569, -2.22940072417}, + { 3.51983380318, -1.08717381954, -1.19617708027}, + {-1.66791832447, -3.11651110649, 2.11557179689}, + {-2.54040765762, -6.29868507385, 1.97103276849}, + {-2.29463744164, -3.35111081600, -5.44719845057}, + {-2.78860569000, -3.85001629591, 1.48611024022}, + { 1.26378631592, 3.41475939751, -2.98826307058}, + { 1.09431362152, 8.47581565380, 7.57644295692}, + { 3.76009845734, -1.30891036987, -1.30899637938}, + {-2.40264558792, -4.04087215662, 9.50866565108}, + { 3.04867124557, -6.51501715183, -4.97306495905}, + { 3.84830522537, -1.05451405048, -2.95348644257}, + { 3.50539922714, -1.34033131599, -4.16487485170}, + {-2.73639702797, -4.54458445311, 4.83948200941}, + {-2.10262560844, 4.50256705284, 8.65258097649}, + {-2.21880722046, -1.73338234425, -9.46770235896}, + {-1.88443505764, -3.78501087427, -4.88811969757}, + {-2.49273109436, -8.57867524028, -3.68066996336}, + {-3.13859176636, 1.89580932260, -7.63508498668}, + {-2.14591693878, -6.56111717224, -6.69820383191}, + {-1.92061448097, -1.09247815609, 6.60725891589}, + { 6.78668081760, 1.96723997593, 4.59519505501}, + { 3.13553071022, -1.15522086620, 5.73987923563}, + {-2.29674005508, -3.97602945566, -8.58206078410}, + { 1.61597287655, 7.94150531292, 1.39395284653}, + { 9.63889718056, 3.76062178612, -2.30398878455}, + { 1.49050402641, 2.90106987953, -1.05920815468}, + { 8.01355421543, 2.98550319672, -1.37276327610}, + { 4.67240428925, -1.42258465290, -7.31541633606}, + { 4.78209877014, -1.97110056877, -6.36375367641}, + { 3.47065544128, -1.58680915833, 8.09270441532}, + { 2.78402256966, -1.61627101898, -1.14950299263}, + {-2.43154764175, -4.92580950260, -5.94577729702}, + {-2.07331848145, -8.07791411877, -5.79017937183}, + {-2.20136833191, -2.79306620359, 1.49220023304}, + { 3.50646018982, -1.30311572552, -1.54289022088}, + {-2.57634282112, -2.89503604174, -1.62051007152}, + {-2.28945779800, -3.16228151321, 1.58045440912}, + {-1.96759450436, -1.22897170484, 3.13766419888}, + {-2.32720947266, -4.58756178617, -1.04388400912}, + { 1.34714412689, 3.28201150894, -3.74540209770}, + { 1.02136373520, 8.49682748317, 8.75190198421}, + { 3.80308532715, -9.79767143726, -7.25016415118}, + {-1.72578215599, 5.16327172518, -1.34024426341}, + { 2.54812169075, -1.19696271420, -4.35636699200}, + { 3.44056987762, -1.42631483078, -1.80410727859}, + { 3.56901502609, -1.25196957588, 2.14892253280}, + {-2.25152993202, -3.58026176691, -2.36085981131}, + {-1.81981575489, -1.61404407024, 6.01518213749}, + {-2.34611868858, 5.45890212059, 1.05074942112}, + {-2.48747754097, 3.01646441221, -2.20733918250}, + {-2.08608031273, -4.99503910542, 1.22879549861}, + {-2.62009620667, -4.38899755478, -2.94447898865}, + {-2.46968364716, -2.14957594872, -4.44929867983}, + {-2.27588725090, -4.03139829636, -1.34325772524}, + { 2.07814240456, 3.53174304962, 4.32420790195}, + { 3.19689464569, -1.74846553802, -9.51488316059}, + {-1.77437961102, 8.71710777282, 7.98717916012}, + {-2.43852794170, 1.01102793217, 1.67076694965}, + { 1.68815839291, 2.98291635513, -3.98489713669}, + { 6.72981083393, 3.35125422478, -8.33267033100}, + { 1.64096879959, 3.26126050949, -6.12493693829}, + { 3.77453780174, 4.46122527122, 6.66481316090}, + { 3.43319153786, -1.30005681515, 3.92319053411}, + { 2.63329458237, -1.30157423019, -8.17687213421}, + { 3.57572197914, -1.07295131683, -4.24419552088}, + {-2.67209243774, -1.46084114909, -1.16768456995}, + {-2.09756040573, -9.31840538979, -3.85717511177}}, + { + {-2.34410619736, -3.20016115904, -1.53496759012}, + { 3.17996025085, -1.40260577202, 1.49473607540}, + {-2.23076605797, -2.83948600292, 9.49927791953}, + {-2.43097519875, -8.68766680360, 1.60800144076}, + {-2.30478429794, -3.56340646744, -4.09480594099}, + {-2.14133548737, -1.02651178837, 4.94684696198}, + { 1.62508022785, 2.60330677032, -8.47915709019}, + { 1.27408051491, 3.01226794720, 4.51113164425}, + { 3.35605812073, -1.12264251709, -3.33058685064}, + {-2.37143301964, -5.74941754341, 8.54486040771}, + { 3.18820738792, -1.45322322845, -2.11916580796}, + { 3.41001844406, -1.34255969524, -1.54219895601}, + { 4.52576208115, -6.47054672241, -2.16511666775}, + {-2.40094542503, -7.25721180439, -1.55527725816}, + {-2.77491641045, -1.10882985592, 5.76599717140}, + {-2.20180344582, -1.91131502390, 2.21937447786}, + {-2.13283038139, -2.67622411251, -3.17741572857}, + {-2.18208360672, 5.69592237473, -2.07313925028}, + {-2.77465915680, -5.78670740128, 4.42580580711}, + {-1.85710799694, -7.07677602768, 1.04370221496}, + {-2.38139748573, -4.66007351875, -9.08390283585}, + { 2.70240306854, 4.33306598663, -4.81943219900}, + { 2.12172913551, -1.01243197918, 1.90536692739}, + {-2.59672832489, 1.63385756314, -4.87916678190}, + { 9.92364227772, 1.40893876553, 1.16456234455}, + { 1.39175999165, 3.11557602882, -4.44381356239}, + { 2.11633038521, 2.02847170830, -1.00864779949}, + { 1.14409208298, 3.74614620209, -7.69796907901}, + { 3.99155473709, -1.15835893154, -5.75888492167}, + { 3.81746459007, -1.76095283031, 3.65874171257}, + { 2.39833283424, -1.97481775284, 1.68805599213}, + { 3.50797653198, -9.54507589340, -7.73615688086}, + {-2.22397685051, -2.59196788073, -5.47018386424}, + {-2.05891585350, 5.35349249840, 8.92746448517}, + {-2.42279815674, -4.47994381189, 4.74890284240}, + { 3.47718238831, -1.31481623650, -1.13119445741}, + {-2.13573265076, -3.77991527319, 9.89178344607}, + {-2.39205574989, -4.24590885639, -2.14120149612}, + {-2.32959675789, -1.04270493612, -2.64487534761}, + {-2.28894376755, -3.51045638323, -4.60519827902}, + { 1.60694050789, 3.09509325027, -3.17743927240}, + { 8.79046201706, 1.23586606979, 1.10633921623}, + { 3.66632819176, -7.73513436317, -2.82783180475}, + {-1.56432127953, -8.28551828861, -1.27556353807}, + { 3.64514565468, -8.48878860474, 1.50680422783}, + { 3.56896424294, -1.43446743488, 2.74687930942}, + { 3.87763309479, -1.23341560364, -8.10135483742}, + {-2.39496254921, -3.45572710037, -4.26582060754}, + {-2.46606898308, -7.99975514412, 2.00696870685}, + {-2.78703904152, -5.71972310543, -1.65262192488}, + {-2.10356879234, -5.14238119124, -1.54197901487}, + {-1.46284854412, 6.09897315502, -8.87724041939}, + {-2.40337014198, 4.84354734421, 3.36634337902}, + {-2.31666541100, -3.93751084805, -5.00837624073}, + {-2.69825482368, 1.31541609764, -2.08565697074}, + { 9.76799368858, 2.24494481087, 6.91881835461}, + { 2.17129302025, -1.59818923473, 2.69582271576}, + {-1.90924882889, 1.96396946907, 1.97196662426}, + { 1.54570734501, 9.02010202408, 8.17995429039}, + { 1.24686288834, 3.31178450584, 1.26904413104}, + { 2.53851819038, 3.38208723068, -4.56276416779}, + { 9.43495273590, 3.29948759079, -1.81205761433}, + { 3.28666305542, -1.16521859169, 6.84504806995}, + { 4.27903270721, 7.15266764164, 1.18705637753}, + { 3.30623006821, -1.17509567738, -2.75256365538}, + { 4.33063077927, -6.61120176315, 1.08258962631}, + {-3.12304520607, 4.37339305878, 1.31159663200}, + {-2.16836428642, -6.58241450787, -1.20764113963}} +}; + + +#+END_src diff --git a/org/qmckl_utils.org b/org/qmckl_utils.org new file mode 100644 index 0000000..2a31e22 --- /dev/null +++ b/org/qmckl_utils.org @@ -0,0 +1,229 @@ +#+TITLE: Utility functions +#+SETUPFILE: ../tools/theme.setup +#+INCLUDE: ../tools/lib.org + +* Headers :noexport: + #+begin_src elisp :noexport :results none +(org-babel-lob-ingest "../tools/lib.org") +#+end_src + + #+begin_src c :comments link :tangle (eval c_test) :noweb yes +#include "qmckl.h" +#include "assert.h" +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif +int main() { + qmckl_context context; + context = qmckl_context_create(); + + #+end_src + +* Matrix operations + +** ~qmckl_transpose~ + + Transposes a matrix: $B_{ji} = A_{ij}$ + + #+NAME: qmckl_transpose_args + | qmckl_context | context | in | Global state | + | int64_t | m | in | Number of rows of the input matrix | + | int64_t | n | in | Number of columns of the input matrix | + | double | A[][lda] | in | Array containing the $m \times n$ matrix $A$ | + | int64_t | lda | in | Leading dimension of array ~A~ | + | double | B[][ldb] | out | Array containing the $n \times m$ matrix $B$ | + | int64_t | ldb | in | Leading dimension of array ~B~ | + +*** Requirements + + - ~context~ is not ~QMCKL_NULL_CONTEXT~ + - ~m > 0~ + - ~n > 0~ + - ~lda >= m~ + - ~ldb >= n~ + - ~A~ is allocated with at least $m \times n \times 8$ bytes + - ~B~ is allocated with at least $n \times m \times 8$ bytes + +*** C header + + #+CALL: generate_c_header(table=qmckl_transpose_args,rettyp="qmckl_exit_code",fname="qmckl_transpose") + + #+RESULTS: + #+begin_src c :tangle (eval h_func) :comments org + qmckl_exit_code qmckl_transpose ( + const qmckl_context context, + const int64_t m, + const int64_t n, + const double* A, + const int64_t lda, + double* const B, + const int64_t ldb ); + #+end_src + +*** Source + #+begin_src f90 :tangle (eval f) +integer function qmckl_transpose_f(context, m, n, A, LDA, B, LDB) & + result(info) + use qmckl + implicit none + integer(qmckl_context) , intent(in) :: context + integer*8 , intent(in) :: m, n + integer*8 , intent(in) :: lda + real*8 , intent(in) :: A(lda,*) + integer*8 , intent(in) :: ldb + real*8 , intent(out) :: B(ldb,*) + + integer*8 :: i,j + + info = QMCKL_SUCCESS + + if (context == QMCKL_NULL_CONTEXT) then + info = QMCKL_INVALID_CONTEXT + return + endif + + if (m <= 0_8) then + info = QMCKL_INVALID_ARG_2 + return + endif + + if (n <= 0_8) then + info = QMCKL_INVALID_ARG_3 + return + endif + + if (LDA < m) then + info = QMCKL_INVALID_ARG_5 + return + endif + + if (LDB < n) then + info = QMCKL_INVALID_ARG_7 + return + endif + + do j=1,m + do i=1,n + B(i,j) = A(j,i) + end do + end do + +end function qmckl_transpose_f + #+end_src + +*** C interface :noexport: + + #+CALL: generate_c_interface(table=qmckl_transpose_args,rettyp="qmckl_exit_code",fname="qmckl_transpose") + + #+RESULTS: + #+begin_src f90 :tangle (eval f) :comments org :exports none + integer(c_int32_t) function qmckl_transpose & + (context, m, n, A, lda, B, ldb) & + bind(C) result(info) + + use, intrinsic :: iso_c_binding + implicit none + + integer (c_int64_t) , intent(in) , value :: context + integer (c_int64_t) , intent(in) , value :: m + integer (c_int64_t) , intent(in) , value :: n + real (c_double ) , intent(in) :: A(lda,*) + integer (c_int64_t) , intent(in) , value :: lda + real (c_double ) , intent(out) :: B(ldb,*) + integer (c_int64_t) , intent(in) , value :: ldb + + integer(c_int32_t), external :: qmckl_transpose_f + info = qmckl_transpose_f & + (context, m, n, A, lda, B, ldb) + + end function qmckl_transpose + #+end_src + + #+CALL: generate_f_interface(table=qmckl_transpose_args,rettyp="qmckl_exit_code",fname="qmckl_transpose") + + #+RESULTS: + #+begin_src f90 :tangle (eval fh_func) :comments org :exports none + interface + integer(c_int32_t) function qmckl_transpose & + (context, m, n, A, lda, B, ldb) & + bind(C) + use, intrinsic :: iso_c_binding + import + implicit none + + integer (c_int64_t) , intent(in) , value :: context + integer (c_int64_t) , intent(in) , value :: m + integer (c_int64_t) , intent(in) , value :: n + real (c_double ) , intent(in) :: A(lda,*) + integer (c_int64_t) , intent(in) , value :: lda + real (c_double ) , intent(out) :: B(ldb,*) + integer (c_int64_t) , intent(in) , value :: ldb + + end function qmckl_transpose + end interface + #+end_src + +*** Test :noexport: + #+begin_src f90 :tangle (eval f_test) +integer(qmckl_exit_code) function test_qmckl_transpose(context) bind(C) + use qmckl + implicit none + integer(qmckl_context), intent(in), value :: context + + double precision, allocatable :: A(:,:), B(:,:) + integer*8 :: m, n, LDA, LDB + integer*8 :: i,j + double precision :: x + + m = 5 + n = 6 + LDA = m+3 + LDB = n+1 + + allocate( A(LDA,n), B(LDB,m) ) + + A = 0.d0 + B = 0.d0 + do j=1,n + do i=1,m + A(i,j) = -10.d0 + dble(i+j) + end do + end do + + test_qmckl_transpose = qmckl_transpose(context, m, n, A, LDA, B, LDB) + + if (test_qmckl_transpose /= QMCKL_SUCCESS) return + + test_qmckl_transpose = QMCKL_FAILURE + + x = 0.d0 + do j=1,n + do i=1,m + x = x + (A(i,j)-B(j,i))**2 + end do + end do + if (dabs(x) <= 1.d-15) then + test_qmckl_transpose = QMCKL_SUCCESS + endif + + deallocate(A,B) +end function test_qmckl_transpose + #+end_src + + #+begin_src c :comments link :tangle (eval c_test) +qmckl_exit_code test_qmckl_transpose(qmckl_context context); +assert(QMCKL_SUCCESS == test_qmckl_transpose(context)); + #+end_src + +* End of files :noexport: + + #+begin_src c :comments link :tangle (eval c_test) + assert (qmckl_context_destroy(context) == QMCKL_SUCCESS); + return 0; +} + + #+end_src + + +# -*- mode: org -*- +# vim: syntax=c diff --git a/org/table_of_contents b/org/table_of_contents index ad39ea1..e9e2190 100644 --- a/org/table_of_contents +++ b/org/table_of_contents @@ -7,3 +7,5 @@ qmckl_nucleus.org qmckl_electron.org qmckl_ao.org qmckl_distance.org +qmckl_utils.org +qmckl_tests.org