Added 'N' and 'T' in coord

This commit is contained in:
Anthony Scemama 2021-05-19 00:28:56 +02:00
parent d78890e31b
commit 1504330500
8 changed files with 1164 additions and 158 deletions

View File

@ -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 \

BIN
org/chbrclf.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 28 KiB

View File

@ -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:

View File

@ -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 ; i<walk_num ; ++i) {
for (int k=0 ; k<3 ; ++k) {
for (int j=0 ; j<num ; ++j) {
,*x = coord_in[i][j][k];
x += 1;
}
}
}
/* --- */
@ -533,7 +545,7 @@ for (size_t i=0 ; i<3*num ; ++i) {
#+end_src
* Computation
The computed data is stored in the context so that it can be reused
by different kernels. To ensure that the data is valid, for each
computed data the date of the context is stored when it is computed.
@ -553,12 +565,12 @@ qmckl_exit_code qmckl_get_electron_ee_distance(qmckl_context context, double* di
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_electron_ee_distance(qmckl_context context, double* distance)
{
qmckl_exit_code rc;
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT;
}
qmckl_exit_code rc;
rc = qmckl_provide_ee_distance(context);
if (rc != QMCKL_SUCCESS) return rc;
@ -581,7 +593,6 @@ qmckl_exit_code qmckl_provide_ee_distance(qmckl_context context);
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_provide_ee_distance(qmckl_context context)
{
qmckl_exit_code rc;
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT;
@ -672,16 +683,15 @@ integer function qmckl_compute_ee_distance_f(context, elec_num, walk_num, coord,
return
endif
!$OMP PARALLEL DO DEFAULT(NONE) &
!$OMP SHARED(elec_num, walk_num, coord, ee_distance)
!$OMP PRIVATE(k)
do k=1,walk_num
info = qmckl_distance(context, 'T', 'T', elec_num, elec_num, &
coord(1,1,k), elec_num, &
coord(1,1,k), elec_num, &
ee_distance(1,1,k), elec_num)
if (info /= QMCKL_SUCCESS) then
exit
endif
end do
!$OMP END PARALLEL DO
end function qmckl_compute_ee_distance_f
#+end_src
@ -722,19 +732,33 @@ qmckl_exit_code qmckl_compute_ee_distance (
*** Test
#+begin_src c :tangle (eval c_test)
/* Reference input data */
assert(qmckl_electron_provided(context));
double ee_distance[walk_num*num*num];
rc = qmckl_get_electron_ee_distance(context, ee_distance);
// (e1,e2,w)
// (0,0,0) == 0.
assert(ee_distance[0] == 0.);
// (1,0,0) == (0,1,0)
assert(ee_distance[1] == ee_distance[num]);
assert(fabs(ee_distance[1]-8.6114953086801) < 1.e-12);
// value of (1,0,0)
assert(fabs(ee_distance[1]-7.152322512964209) < 1.e-12);
// (0,0,1) == 0.
assert(ee_distance[num*num] == 0.);
// (1,0,1) == (0,1,1)
assert(ee_distance[num*num+1] == ee_distance[num*num+num]);
// value of (1,0,1)
assert(fabs(ee_distance[num*num+1]-6.5517646321055665) < 1.e-12);
#+end_src
** Electron-nucleus distances
** Electron-nucleus distances
*** Get
@ -745,12 +769,13 @@ qmckl_exit_code qmckl_get_electron_en_distance(qmckl_context context, double* di
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_electron_en_distance(qmckl_context context, double* distance)
{
qmckl_exit_code rc;
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT;
}
qmckl_exit_code rc;
rc = qmckl_provide_en_distance(context);
if (rc != QMCKL_SUCCESS) return rc;
@ -773,7 +798,6 @@ qmckl_exit_code qmckl_provide_en_distance(qmckl_context context);
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_provide_en_distance(qmckl_context context)
{
qmckl_exit_code rc;
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT;
@ -838,7 +862,7 @@ qmckl_exit_code qmckl_provide_en_distance(qmckl_context context)
| int64_t | nucl_num | in | Number of nuclei |
| int64_t | walk_num | in | Number of walkers |
| double | elec_coord[walk_num][3][elec_num] | in | Electron coordinates |
| double | nucl_coord[walk_num][3][elec_num] | in | Nuclear coordinates |
| double | nucl_coord[3][elec_num] | in | Nuclear coordinates |
| double | en_distance[walk_num][nucl_num][elec_num] | out | Electron-nucleus distances |
#+begin_src f90 :comments org :tangle (eval f) :noweb yes
@ -851,7 +875,7 @@ integer function qmckl_compute_en_distance_f(context, elec_num, nucl_num, walk_n
integer*8 , intent(in) :: nucl_num
integer*8 , intent(in) :: walk_num
double precision , intent(in) :: elec_coord(elec_num,3,walk_num)
double precision , intent(in) :: nucl_coord(nucl_num,3,walk_num)
double precision , intent(in) :: nucl_coord(nucl_num,3)
double precision , intent(out) :: en_distance(elec_num,nucl_num,walk_num)
integer*8 :: k
@ -878,16 +902,15 @@ integer function qmckl_compute_en_distance_f(context, elec_num, nucl_num, walk_n
return
endif
!$OMP PARALLEL DO DEFAULT(NONE) &
!$OMP SHARED(elec_num, nucl_num, walk_num, elec_coord, nucl_coord, en_distance)
!$OMP PRIVATE(k)
do k=1,walk_num
info = qmckl_distance(context, 'T', 'T', elec_num, nucl_num, &
elec_coord(1,1,k), elec_num, &
nucl_coord, nucl_num, &
en_distance(1,1,k), elec_num)
if (info /= QMCKL_SUCCESS) then
exit
endif
end do
!$OMP END PARALLEL DO
end function qmckl_compute_en_distance_f
#+end_src
@ -919,7 +942,7 @@ qmckl_exit_code qmckl_compute_en_distance (
integer (c_int64_t) , intent(in) , value :: nucl_num
integer (c_int64_t) , intent(in) , value :: walk_num
real (c_double ) , intent(in) :: elec_coord(elec_num,3,walk_num)
real (c_double ) , intent(in) :: nucl_coord(elec_num,3,walk_num)
real (c_double ) , intent(in) :: nucl_coord(elec_num,3)
real (c_double ) , intent(out) :: en_distance(elec_num,nucl_num,walk_num)
integer(c_int32_t), external :: qmckl_compute_en_distance_f
@ -932,13 +955,40 @@ qmckl_exit_code qmckl_compute_en_distance (
*** Test
#+begin_src c :tangle (eval c_test)
/* Reference input data */
assert(!qmckl_nucleus_provided(context));
assert(qmckl_electron_provided(context));
double en_distance[walk_num*num*num];
rc = qmckl_get_electron_en_distance(context, en_distance);
assert(rc == QMCKL_NOT_PROVIDED);
rc = qmckl_set_nucleus_num (context, nucl_num);
assert(rc == QMCKL_SUCCESS);
rc = qmckl_set_nucleus_charge (context, charge);
assert (rc == QMCKL_SUCCESS);
rc = qmckl_set_nucleus_coord (context, &(nucl_coord[0]));
assert (rc == QMCKL_SUCCESS);
assert(qmckl_nucleus_provided(context));
double en_distance[walk_num][nucl_num][num];
rc = qmckl_get_electron_en_distance(context, &(en_distance[0][0][0]));
assert (rc == QMCKL_SUCCESS);
// (e,n,w) in Fortran notation
// (1,1,1)
assert(fabs(en_distance[0][0][0] - 6.855624268793153) < 1.e-12);
// (2,1,1)
assert(fabs(en_distance[0][0][1] - 3.698922010513608) < 1.e-12);
// (1,2,1)
assert(fabs(en_distance[0][1][0] - 8.143800105434433) < 1.e-12);
// (2,2,1)
assert(fabs(en_distance[0][1][1] - 5.16360835635664) < 1.e-12);
// (4,3,2)
assert(fabs(en_distance[1][2][3] - 12.599138999960012) < 1.e-12);
#+end_src

View File

@ -1,10 +1,10 @@
#i+TITLE: Nucleus
# +SETUPFILE: ../tools/theme.setup
# +INCLUDE: ../tools/lib.org
#+TITLE: Nucleus
#+SETUPFILE: ../tools/theme.setup
#+INCLUDE: ../tools/lib.org
A ll the data relative to the molecular geometry is described here.
All the data relative to the molecular geometry is described here.
* Headers :noexport:
* Headers :noexport:
#+begin_src elisp :noexport :results none
( org-babel-lob-ingest "../tools/lib.org")
#+end_src
@ -24,6 +24,9 @@ A ll the data relative to the molecular geometry is described here.
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include "chbrclf.h"
int main() {
qmckl_context context;
context = qmckl_context_create();
@ -63,7 +66,7 @@ int main() {
| ~num~ | int64_t | Total number of nuclei |
| ~provided~ | bool | If true, ~nucleus~ is valid |
| ~charge~ | double[num] | Nuclear charges |
| ~coord~ | double[3][num] | Nuclear coordinates |
| ~coord~ | double[3][num] | Nuclear coordinates, in transposed format |
| ~nn_distance~ | double[num][num] | Nucleus-nucleus distances |
| ~nn_distance_date~ | int64_t | Date when Nucleus-nucleus distances were computed |
| ~repulsion~ | double | Nuclear repulsion energy |
@ -93,9 +96,9 @@ typedef struct qmckl_nucleus_struct {
** Access functions
#+begin_src c :comments org :tangle (eval h_func) :exports none
qmckl_exit_code 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_exit_code qmckl_get_nucleus_coord (const qmckl_context context, double* coord);
qmckl_exit_code qmckl_get_nucleus_num (const qmckl_context context, int64_t* const num);
qmckl_exit_code qmckl_get_nucleus_charge (const qmckl_context context, double* const charge);
qmckl_exit_code qmckl_get_nucleus_coord (const qmckl_context context, const char transp, double* const coord);
#+end_src
#+NAME:post
@ -107,19 +110,29 @@ if ( (ctx->nucleus.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) {
<<pre2>>
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);
<<post2>>
}
#+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) {
<<pre2>>
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));
}
<<post2>>
}
#+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 ; i<nucl_num ; ++i) {
assert( nucl_coord[nucl_num*k+i] == nucl_coord2[3*i+k] );
}
}
rc = qmckl_get_nucleus_charge(context, charge2);
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2);
assert(rc == QMCKL_SUCCESS);
for (size_t i=0 ; i<3*nucl_num ; ++i) {
assert( nucl_coord[i] == nucl_coord2[i] );
}
double nucl_charge2[nucl_num];
rc = qmckl_get_nucleus_charge(context, nucl_charge2);
assert(rc == QMCKL_NOT_PROVIDED);
rc = qmckl_set_nucleus_charge(context, charge);
rc = qmckl_set_nucleus_charge(context, nucl_charge);
assert(rc == QMCKL_SUCCESS);
rc = qmckl_get_nucleus_charge(context, charge2);
rc = qmckl_get_nucleus_charge(context, nucl_charge2);
assert(rc == QMCKL_SUCCESS);
for (size_t i=0 ; i<num ; ++i) {
assert( charge[i] == charge2[i] );
for (size_t i=0 ; i<nucl_num ; ++i) {
assert( nucl_charge[i] == nucl_charge2[i] );
}
assert(qmckl_nucleus_provided(context));
#+end_src
@ -462,19 +512,19 @@ qmckl_exit_code qmckl_provide_nn_distance(qmckl_context context)
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return (char) 0;
}
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL);
if (!ctx->nucleus.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 -*-

677
org/qmckl_tests.org Normal file
View File

@ -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

229
org/qmckl_utils.org Normal file
View File

@ -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

View File

@ -7,3 +7,5 @@ qmckl_nucleus.org
qmckl_electron.org
qmckl_ao.org
qmckl_distance.org
qmckl_utils.org
qmckl_tests.org