TREX/qmckl
1
0
Fork 0
mirror of https://github.com/TREX-CoE/qmckl.git synced 2025-04-30 04:15:00 +02:00
Code Issues Projects Releases Wiki Activity

Merge branch 'jastrow_c' into qmckl_compute_een_rescaled_n

Browse Source
This commit is contained in:
Gianfranco Abrusci 2022-02-22 14:55:55 +01:00
commit 3bd249176f
16 changed files with 1515 additions and 423 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

28
.github/workflows/test-build.yml vendored
View File

@ -40,7 +40,9 @@ jobs:
- name: Build QMCkl - name: Build QMCkl
run: | run: |
./autogen.sh ./autogen.sh
./configure --enable-silent-rules --enable-debug mkdir _build
cd _build
../configure --enable-silent-rules --enable-debug
make -j 4 make -j 4
- name: Run test - name: Run test
@ -72,7 +74,7 @@ jobs:
- uses: actions/checkout@v2 - uses: actions/checkout@v2
- name: install dependencies - name: install dependencies
run: brew install emacs hdf5 automake pkg-config run: brew install emacs hdf5 automake pkg-config
- name: Symlink gfortran (macOS) - name: Symlink gfortran (macOS)
if: runner.os == 'macOS' if: runner.os == 'macOS'
run: | run: |
@ -89,6 +91,28 @@ jobs:
git clone https://github.com/TREX-CoE/trexio.git git clone https://github.com/TREX-CoE/trexio.git
cd trexio cd trexio
./autogen.sh ./autogen.sh
mkdir _build
cd _build
../configure --enable-silent-rules
make -j 4
- name: Run test
run: make -j 4 check
- name: Archive test log file
if: failure()
uses: actions/upload-artifact@v2
with:
name: test-report-ubuntu
path: test-suite.log
- name: Dist test
run: make distcheck
- name: Archive test log file
if: failure()
uses: actions/upload-artifact@v2
with:
./configure --prefix=${PWD}/_install --enable-silent-rules ./configure --prefix=${PWD}/_install --enable-silent-rules
make -j 4 make -j 4
make install make install

16
Makefile.am
View File

@ -82,18 +82,17 @@ htmldir_loc = share/doc/qmckl/html
textdir_loc = share/doc/qmckl/text textdir_loc = share/doc/qmckl/text
htmlize_el = $(htmldir_loc)/htmlize.el htmlize_el = $(htmldir_loc)/htmlize.el
dist_html_DATA = $(HTML_FILES) \ dist_html_DATA = $(HTML_FILES) \
share/doc/qmckl/html/qmckl.css \ $(htmldir_loc)/qmckl.css \
share/doc/qmckl/html/index.html $(htmldir_loc)/index.html
dist_text_DATA = $(TEXT_FILES) dist_text_DATA = $(TEXT_FILES)
html-local: $(dist_html_DATA) html-local: $(htmlize_el) $(dist_html_DATA)
text: $(dist_text_DATA) text: $(htmlize_el) $(dist_text_DATA)
doc: html text doc: html text
if QMCKL_DEVEL if QMCKL_DEVEL
if VFC_CI if VFC_CI
@ -102,9 +101,9 @@ endif
EXTRA_DIST += $(ORG_FILES) $(TANGLED_FILES) $(EXPORTED_FILES) $(test_qmckl_f) EXTRA_DIST += $(ORG_FILES) $(TANGLED_FILES) $(EXPORTED_FILES) $(test_qmckl_f)
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) $(src_qmckl_f) $(test_qmckl_f) $(qmckl_h) $(header_tests) 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) $(src_qmckl_f) $(test_qmckl_f) $(qmckl_h) $(header_tests) $(htmlize_el)
CLEANFILES += $(BUILT_SOURCES) $(C_TEST_FILES) $(F_TEST_FILES) $(TANGLED_FILES) $(C_TEST_FILES) $(F_TEST_FILES) $(src_qmckl_f) $(test_qmckl_f) $(qmckl_h) $(HTML_FILES) $(TEXT_FILES) $(EXPORTED_FILES) $(header_tests) CLEANFILES += $(BUILT_SOURCES) $(C_TEST_FILES) $(F_TEST_FILES) $(TANGLED_FILES) $(C_TEST_FILES) $(F_TEST_FILES) $(src_qmckl_f) $(test_qmckl_f) $(qmckl_h) $(HTML_FILES) $(TEXT_FILES) $(EXPORTED_FILES) $(header_tests) $(htmlize_el)
EXTRA_DIST += \ EXTRA_DIST += \
tools/build_doc.sh \ tools/build_doc.sh \
@ -152,8 +151,7 @@ $(src_qmckl_f): $(FH_FUNC_FILES) $(FH_TYPE_FILES)
$(cat_h_verbose)top_builddir=$(abs_top_builddir) srcdir=$(abs_srcdir) src_qmckl_f=$(src_qmckl_f) $(srcdir)/tools/build_qmckl_f.sh $(cat_h_verbose)top_builddir=$(abs_top_builddir) srcdir=$(abs_srcdir) src_qmckl_f=$(src_qmckl_f) $(srcdir)/tools/build_qmckl_f.sh
$(htmlize_el): $(htmlize_el):
$(MKDIR_P) $(htmldir_loc) $(MKDIR_P) $(htmldir_loc) $(textdir_loc)
$(MKDIR_P) $(textdir_loc)
abs_srcdir=$(abs_srcdir) $(srcdir)/tools/install_htmlize.sh $(htmlize_el) abs_srcdir=$(abs_srcdir) $(srcdir)/tools/install_htmlize.sh $(htmlize_el)
tests/chbrclf.h: $(qmckl_h) tests/chbrclf.h: $(qmckl_h)

42
README.md
View File

@ -36,12 +36,39 @@ make check
## For users ## For users
Obtain a source distribution and run Obtain a source distribution.
To build the documentation version:
``` ```
./configure ./configure
make ```
make check
To build an optimized version with Intel compilers:
```
./configure \
--with-icc \
--with-ifort \
--enable-hpc \
--with-openmp
```
To build an optimized version with GCC:
```
./configure \
CC=gcc \
CFLAGS="-g -O2 -march=native -flto -fno-trapping-math -fno-math-errno -ftree-vectorize" \
FC=gfortran \
FCFLAGS="-g -O2 -march=native -flto -ftree-vectorize" \
--enable-hpc \
--with-openmp
```
Then, compile with:
```
make -j
make -j check
sudo make install sudo make install
sudo make installcheck sudo make installcheck
``` ```
@ -54,7 +81,12 @@ by the preprocessor otherwise). To enable vfc_ci support, the library should be
configured with the following command : configured with the following command :
``` ```
./configure --prefix=$PWD/_install \ --enable-vfc_ci --host=x86_64 \ CC="verificarlo-f" FC="verificarlo-f" ./configure \
CC="verificarlo-f" \
FC="verificarlo-f" \
--prefix=$PWD/_install \
--enable-vfc_ci \
--host=x86_64 \
``` ```
where CC and FC are set to verificarlo-f, and support is explicitely enabled where CC and FC are set to verificarlo-f, and support is explicitely enabled

49
configure.ac
View File

@ -46,6 +46,22 @@ AS_IF([test -d ${srcdir}/.git],
[enable_maintainer_mode="no"] [enable_maintainer_mode="no"]
) )
AC_ARG_WITH(ifort, [AS_HELP_STRING([--with-ifort],[Use Intel Fortran compiler])], with_ifort=$withval, with_ifort=no)
AS_IF([test "$with_ifort" == "yes"], [
FC=ifort
FCFLAGS="-xHost -ip -O2 -ftz -finline -g -mkl=sequential" ])
AC_ARG_WITH(icc, [AS_HELP_STRING([--with-icc],[Use Intel C compiler])], with_icc=$withval, with_icc=no)
AS_IF([test "$with_icc" == "yes"], [
CC=icc
CFLAGS="-xHost -ip -O2 -ftz -finline -g -mkl=sequential" ])
AS_IF([test "$with_icc"."$with_ifort" == "yes.yes"], [
ax_blas_ok="yes"
ax_lapack_ok="yes"
BLAS_LIBS=""
LAPACK_LIBS=""])
AM_PROG_AR AM_PROG_AR
AM_MAINTAINER_MODE() AM_MAINTAINER_MODE()
LT_INIT LT_INIT
@ -100,13 +116,14 @@ AC_CHECK_HEADERS([assert.h errno.h math.h pthread.h stdbool.h stdint.h stdio.h s
AC_CHECK_LIB([pthread], [pthread_create]) AC_CHECK_LIB([pthread], [pthread_create])
# OpenMP # OpenMP
#AC_ARG_WITH(openmp, [AS_HELP_STRING([--with-openmp],[enable OpenMP])], with_omp=$withval, with_omp=no) AC_ARG_WITH(openmp, [AS_HELP_STRING([--with-openmp],[activate OpenMP])], with_omp=$withval, with_omp=no)
#if test "x$with_omp" = xyes; then if test "x$with_omp" = xyes; then
# AC_DEFINE([HAVE_OPENMP], [1], [Define to use OpenMP threading.]) AC_DEFINE([HAVE_OPENMP], [1], [Define to use OpenMP threading.])
# AX_OPENMP([], AX_OPENMP([],
# [AC_MSG_ERROR([Could not find OpenMP flags; configure with --without-openmp])]) [AC_MSG_ERROR([Could not find OpenMP flags; configure with --without-openmp])])
# CFLAGS="${CFLAGS} ${OPENMP_CFLAGS}" CFLAGS="${CFLAGS} ${OPENMP_CFLAGS}"
#fi FCFLAGS="${CFLAGS} ${OPENMP_FCFLAGS}"
fi
# CHAMELEON # CHAMELEON
AC_ARG_WITH(chameleon, AC_ARG_WITH(chameleon,
@ -177,6 +194,7 @@ AX_BLAS([], [AC_MSG_ERROR([BLAS was not found.])])
## LAPACK ## LAPACK
AX_LAPACK([], [AC_MSG_ERROR([LAPACK was not found.])]) AX_LAPACK([], [AC_MSG_ERROR([LAPACK was not found.])])
AS_IF([test "$BLAS_LIBS" == "$LAPACK_LIBS"], [BLAS_LIBS=""])
# Specific options required with some compilers # Specific options required with some compilers
@ -186,10 +204,10 @@ case $FC in
FCFLAGS="$FCFLAGS -nofor-main" FCFLAGS="$FCFLAGS -nofor-main"
;; ;;
gfortran*) #gfortran*)
# Order is important here # # Order is important here
FCFLAGS="-cpp $FCFLAGS" # FCFLAGS="-cpp $FCFLAGS"
;; # ;;
esac esac
@ -204,12 +222,13 @@ AC_ARG_ENABLE(debug, [AS_HELP_STRING([--enable-debug],[compile for debugging])],
if test "$ok" = "yes"; then if test "$ok" = "yes"; then
if test "$GCC" = "yes"; then if test "$GCC" = "yes"; then
CFLAGS="$CFLAGS \ CFLAGS="$CFLAGS \
-Wall -W -Wbad-function-cast -Wcast-qual \ -g -Wall -W -Wbad-function-cast -Wcast-qual -Warray-bounds -Wdisabled-optimization \
-Wpointer-arith -Wcast-align -Wpedantic -Wextra -fmax-errors=3" -Wpointer-arith -Wcast-align -Wpedantic -Wextra -Walloc-zero -Werror \
"
fi fi
if test "$GFC" = "yes"; then if test "$GFC" = "yes"; then
FCFLAGS="$FCFLAGS \ FCFLAGS="$FCFLAGS \
-fcheck=all -Waliasing -Wampersand -Wconversion \ -g -fcheck=all -Waliasing -Wampersand -Wconversion \
-Wsurprising -ffpe-trap=zero,overflow,underflow \ -Wsurprising -ffpe-trap=zero,overflow,underflow \
-Wintrinsics-std -Wno-tabs -Wintrinsic-shadow -Wline-truncation \ -Wintrinsics-std -Wno-tabs -Wintrinsic-shadow -Wline-truncation \
-Wreal-q-constant -Wuninitialized -fbacktrace -finit-real=nan" -Wreal-q-constant -Wuninitialized -fbacktrace -finit-real=nan"
@ -314,7 +333,7 @@ fi
#mkl-dynamic-lp64-seq #mkl-dynamic-lp64-seq
PKG_LIBS="$PKG_LIBS $LIBS" PKG_LIBS="$PKG_LIBS $LIBS"
LIBS="$BLAS_LIBS $LAPACK_LIBS $BLAS_LIBS $PKG_LIBS" LIBS="$BLAS_LIBS $LAPACK_LIBS $PKG_LIBS"
CFLAGS="$CFLAGS $PKG_CFLAGS" CFLAGS="$CFLAGS $PKG_CFLAGS"
AC_SUBST([PKG_LIBS]) AC_SUBST([PKG_LIBS])
AC_SUBST([PKG_CFLAGS]) AC_SUBST([PKG_CFLAGS])

1218
org/qmckl_ao.org
View File

File diff suppressed because it is too large Load Diff

42
org/qmckl_blas.org
View File

@ -111,26 +111,26 @@ qmckl_vector_alloc( qmckl_context context,
#+begin_src c :comments org :tangle (eval h_private_func) #+begin_src c :comments org :tangle (eval h_private_func)
qmckl_exit_code qmckl_exit_code
qmckl_vector_free( qmckl_context context, qmckl_vector_free( qmckl_context context,
qmckl_vector vector); qmckl_vector* vector);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) #+begin_src c :comments org :tangle (eval c)
qmckl_exit_code qmckl_exit_code
qmckl_vector_free( qmckl_context context, qmckl_vector_free( qmckl_context context,
qmckl_vector vector) qmckl_vector* vector)
{ {
/* Always true */ /* Always true */
assert (vector.data != NULL); assert (vector->data != NULL);
qmckl_exit_code rc; qmckl_exit_code rc;
rc = qmckl_free(context, vector.data); rc = qmckl_free(context, vector->data);
if (rc != QMCKL_SUCCESS) { if (rc != QMCKL_SUCCESS) {
return rc; return rc;
} }
vector.size = (int64_t) 0; vector->size = (int64_t) 0;
vector.data = NULL; vector->data = NULL;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
#+end_src #+end_src
@ -192,26 +192,26 @@ qmckl_matrix_alloc( qmckl_context context,
#+begin_src c :comments org :tangle (eval h_private_func) #+begin_src c :comments org :tangle (eval h_private_func)
qmckl_exit_code qmckl_exit_code
qmckl_matrix_free( qmckl_context context, qmckl_matrix_free( qmckl_context context,
qmckl_matrix matrix); qmckl_matrix* matrix);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) #+begin_src c :comments org :tangle (eval c)
qmckl_exit_code qmckl_exit_code
qmckl_matrix_free( qmckl_context context, qmckl_matrix_free( qmckl_context context,
qmckl_matrix matrix) qmckl_matrix* matrix)
{ {
/* Always true */ /* Always true */
assert (matrix.data != NULL); assert (matrix->data != NULL);
qmckl_exit_code rc; qmckl_exit_code rc;
rc = qmckl_free(context, matrix.data); rc = qmckl_free(context, matrix->data);
if (rc != QMCKL_SUCCESS) { if (rc != QMCKL_SUCCESS) {
return rc; return rc;
} }
matrix.data = NULL; matrix->data = NULL;
matrix.size[0] = (int64_t) 0; matrix->size[0] = (int64_t) 0;
matrix.size[1] = (int64_t) 0; matrix->size[1] = (int64_t) 0;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -286,25 +286,25 @@ qmckl_tensor_alloc( qmckl_context context,
#+begin_src c :comments org :tangle (eval h_private_func) #+begin_src c :comments org :tangle (eval h_private_func)
qmckl_exit_code qmckl_exit_code
qmckl_tensor_free( qmckl_context context, qmckl_tensor_free( qmckl_context context,
qmckl_tensor tensor); qmckl_tensor* tensor);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) #+begin_src c :comments org :tangle (eval c)
qmckl_exit_code qmckl_exit_code
qmckl_tensor_free( qmckl_context context, qmckl_tensor_free( qmckl_context context,
qmckl_tensor tensor) qmckl_tensor* tensor)
{ {
/* Always true */ /* Always true */
assert (tensor.data != NULL); assert (tensor->data != NULL);
qmckl_exit_code rc; qmckl_exit_code rc;
rc = qmckl_free(context, tensor.data); rc = qmckl_free(context, tensor->data);
if (rc != QMCKL_SUCCESS) { if (rc != QMCKL_SUCCESS) {
return rc; return rc;
} }
memset(&tensor, 0, sizeof(qmckl_tensor)); memset(tensor, 0, sizeof(qmckl_tensor));
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -712,7 +712,7 @@ qmckl_tensor_of_double(const qmckl_context context,
for (int64_t i=0 ; i<p ; ++i) for (int64_t i=0 ; i<p ; ++i)
assert ( qmckl_vec(vec2, i) == qmckl_vec(vec, i) ) ; assert ( qmckl_vec(vec2, i) == qmckl_vec(vec, i) ) ;
qmckl_vector_free(context, vec); qmckl_vector_free(context, &vec);
} }
#+end_src #+end_src
@ -2121,8 +2121,8 @@ qmckl_transpose (qmckl_context context,
for (int i=0 ; i<2 ; ++i) for (int i=0 ; i<2 ; ++i)
assert (qmckl_mat(A, i, j) == qmckl_mat(At, j, i)); assert (qmckl_mat(A, i, j) == qmckl_mat(At, j, i));
qmckl_matrix_free(context, A); qmckl_matrix_free(context, &A);
qmckl_matrix_free(context, At); qmckl_matrix_free(context, &At);
} }
#+end_src #+end_src

38
org/qmckl_context.org
View File

@ -141,15 +141,6 @@ typedef struct qmckl_context_struct {
} qmckl_context_struct; } qmckl_context_struct;
#+end_src #+end_src
The context keeps a ``date'' that allows to check which data needs
to be recomputed. The date is incremented when the electron
coordinates are updated.
When a new element is added to the context, the functions
[[Creation][qmckl_context_create]], [[Destroy][qmckl_context_destroy]] and [[Copy][qmckl_context_copy]]
should be updated inorder to make deep copies.
A tag is used internally to check if the memory domain pointed A tag is used internally to check if the memory domain pointed
by a pointer is a valid context. This allows to check that even if by a pointer is a valid context. This allows to check that even if
the pointer associated with a context is non-null, we can still the pointer associated with a context is non-null, we can still
@ -185,6 +176,35 @@ qmckl_context qmckl_context_check(const qmckl_context context) {
} }
#+end_src #+end_src
The context keeps a ``date'' that allows to check which data needs
to be recomputed. The date is incremented when the context is touched.
When a new element is added to the context, the functions
[[Creation][qmckl_context_create]], [[Destroy][qmckl_context_destroy]] and [[Copy][qmckl_context_copy]]
should be updated in order to make deep copies.
#+begin_src c :comments org :tangle (eval h_func) :noexport
qmckl_exit_code qmckl_context_touch(const qmckl_context context) ;
#+end_src
#+begin_src c :tangle (eval c)
qmckl_exit_code qmckl_context_touch(const qmckl_context context) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return qmckl_failwith( context,
QMCKL_INVALID_CONTEXT,
"qmckl_context_touch",
NULL);
}
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
ctx->date += 1UL;
ctx->point.date += 1UL;
return QMCKL_SUCCESS;
}
#+end_src
** Creation ** Creation
To create a new context, ~qmckl_context_create()~ should be used. To create a new context, ~qmckl_context_create()~ should be used.

2
org/qmckl_determinant.org
View File

@ -316,7 +316,7 @@ int64_t* qmckl_get_determinant_mo_index_alpha (const qmckl_context context) {
int32_t mask = 1 << 4; int32_t mask = 1 << 4;
if ( (ctx->det.uninitialized & mask) != 0) { if ( (ctx->det.uninitialized & mask) != 0) {
return (int64_t) 0; return NULL;
} }
assert (ctx->det.mo_index_alpha != NULL); assert (ctx->det.mo_index_alpha != NULL);

15
org/qmckl_electron.org
View File

@ -499,18 +499,15 @@ ctx->electron.provided = (ctx->electron.uninitialized == 0);
if (ctx->electron.provided) { if (ctx->electron.provided) {
if (ctx->electron.coord_new.data != NULL) { if (ctx->electron.coord_new.data != NULL) {
const qmckl_exit_code rc = qmckl_matrix_free(context, ctx->electron.coord_new); const qmckl_exit_code rc = qmckl_matrix_free(context, &(ctx->electron.coord_new));
assert (rc == QMCKL_SUCCESS); assert (rc == QMCKL_SUCCESS);
} }
if (ctx->electron.coord_old.data != NULL) { if (ctx->electron.coord_old.data != NULL) {
const qmckl_exit_code rc = qmckl_matrix_free(context, ctx->electron.coord_old); const qmckl_exit_code rc = qmckl_matrix_free(context, &(ctx->electron.coord_old));
assert (rc == QMCKL_SUCCESS); assert (rc == QMCKL_SUCCESS);
} }
const int64_t walk_num = ctx->electron.walk_num; const int64_t point_num = ctx->electron.walk_num * ctx->electron.num;
const int64_t elec_num = ctx->electron.num;
const int64_t point_num = walk_num * elec_num;
qmckl_matrix coord_new = qmckl_matrix_alloc(context, point_num, 3); qmckl_matrix coord_new = qmckl_matrix_alloc(context, point_num, 3);
@ -1667,8 +1664,7 @@ qmckl_exit_code qmckl_provide_ee_potential(qmckl_context context)
ctx->electron.ee_pot = ee_pot; ctx->electron.ee_pot = ee_pot;
} }
qmckl_exit_code rc = rc = qmckl_compute_ee_potential(context,
qmckl_compute_ee_potential(context,
ctx->electron.num, ctx->electron.num,
ctx->electron.walk_num, ctx->electron.walk_num,
ctx->electron.ee_distance, ctx->electron.ee_distance,
@ -2712,8 +2708,7 @@ qmckl_exit_code qmckl_provide_en_potential(qmckl_context context)
ctx->electron.en_pot = en_pot; ctx->electron.en_pot = en_pot;
} }
qmckl_exit_code rc = rc = qmckl_compute_en_potential(context,
qmckl_compute_en_potential(context,
ctx->electron.num, ctx->electron.num,
ctx->nucleus.num, ctx->nucleus.num,
ctx->electron.walk_num, ctx->electron.walk_num,

400
org/qmckl_jastrow.org
View File

@ -9,10 +9,10 @@
nuclear ($\mathbf{R}$) coordinates. Its defined as $\exp(J(\mathbf{r},\mathbf{R}))$, where nuclear ($\mathbf{R}$) coordinates. Its defined as $\exp(J(\mathbf{r},\mathbf{R}))$, where
\[ \[
J(\mathbf{r},\mathbf{R}) = J_{\text{eN}}(\mathbf{r},\mathbf{R}) + J_{\text{ee}}(\mathbf{r}) + J_{\text{eeN}}(\mathbf{r},\mathbf{R}) J(\mathbf{r},\mathbf{R}) = J_{\text{eN}}(\mathbf{r},\mathbf{R}) + J_{\text{ee}}(\mathbf{r}) + J_{\text{eeN}}(\mathbf{r},\mathbf{R})
\] \]
In the following, we us the notations $r_{ij} = |\mathbf{r}_i - \mathbf{r}_j|$ and In the following, we us the notations $r_{ij} = |\mathbf{r}_i - \mathbf{r}_j|$ and
$R_{i\alpha} = |\mathbf{r}_i - \mathbf{R}_\alpha|$. $R_{i\alpha} = |\mathbf{r}_i - \mathbf{R}_\alpha|$.
$J_{\text{eN}}$ contains electron-nucleus terms: $J_{\text{eN}}$ contains electron-nucleus terms:
@ -23,9 +23,9 @@
\sum_{p=2}^{N_\text{ord}^a} a_{p+1}\, [f(R_{i\alpha})]^p - J_{eN}^\infty \sum_{p=2}^{N_\text{ord}^a} a_{p+1}\, [f(R_{i\alpha})]^p - J_{eN}^\infty
\] \]
$J_{\text{ee}}$ contains electron-electron terms: $J_{\text{ee}}$ contains electron-electron terms:
\[ \[
J_{\text{ee}}(\mathbf{r}) = J_{\text{ee}}(\mathbf{r}) =
\sum_{i=1}^{N_\text{elec}} \sum_{j=1}^{i-1} \sum_{i=1}^{N_\text{elec}} \sum_{j=1}^{i-1}
\frac{b_1\, f(r_{ij})}{1+b_2\, f(r_{ij})} + \frac{b_1\, f(r_{ij})}{1+b_2\, f(r_{ij})} +
\sum_{p=2}^{N_\text{ord}^b} a_{p+1}\, [f(r_{ij})]^p - J_{ee}^\infty \sum_{p=2}^{N_\text{ord}^b} a_{p+1}\, [f(r_{ij})]^p - J_{ee}^\infty
@ -43,7 +43,7 @@
\sum_{l=0}^{p-k-2\delta_{k,0}} \sum_{l=0}^{p-k-2\delta_{k,0}}
c_{lkp\alpha} \left[ g({r}_{ij}) \right]^k c_{lkp\alpha} \left[ g({r}_{ij}) \right]^k
\left[ \left[ g({R}_{i\alpha}) \right]^l + \left[ g({R}_{j\alpha}) \right]^l \right] \left[ \left[ g({R}_{i\alpha}) \right]^l + \left[ g({R}_{j\alpha}) \right]^l \right]
\left[ g({R}_{i\,\alpha}) \, g({R}_{j\alpha}) \right]^{(p-k-l)/2} \left[ g({R}_{i\,\alpha}) \, g({R}_{j\alpha}) \right]^{(p-k-l)/2}
\] \]
$c_{lkp\alpha}$ are non-zero only when $p-k-l$ is even. $c_{lkp\alpha}$ are non-zero only when $p-k-l$ is even.
@ -55,10 +55,10 @@
g(r) = e^{-\kappa\, r}. g(r) = e^{-\kappa\, r}.
\] \]
The terms $J_{\text{ee}}^\infty$ and $J_{\text{eN}}^\infty$ are shifts to ensure that The terms $J_{\text{ee}}^\infty$ and $J_{\text{eN}}^\infty$ are shifts to ensure that
$J_{\text{ee}}$ and $J_{\text{eN}}$ have an asymptotic value of zero. $J_{\text{ee}}$ and $J_{\text{eN}}$ have an asymptotic value of zero.
* Headers :noexport: * Headers :noexport:
#+begin_src elisp :noexport :results none #+begin_src elisp :noexport :results none
(org-babel-lob-ingest "../tools/lib.org") (org-babel-lob-ingest "../tools/lib.org")
@ -412,10 +412,10 @@ qmckl_exit_code qmckl_get_jastrow_aord_num (qmckl_context context, int
qmckl_exit_code qmckl_get_jastrow_bord_num (qmckl_context context, int64_t* const bord_num); qmckl_exit_code qmckl_get_jastrow_bord_num (qmckl_context context, int64_t* const bord_num);
qmckl_exit_code qmckl_get_jastrow_cord_num (qmckl_context context, int64_t* const bord_num); qmckl_exit_code qmckl_get_jastrow_cord_num (qmckl_context context, int64_t* const bord_num);
qmckl_exit_code qmckl_get_jastrow_type_nucl_num (qmckl_context context, int64_t* const type_nucl_num); qmckl_exit_code qmckl_get_jastrow_type_nucl_num (qmckl_context context, int64_t* const type_nucl_num);
qmckl_exit_code qmckl_get_jastrow_type_nucl_vector (qmckl_context context, int64_t* const type_nucl_num, int64_t* size_max); qmckl_exit_code qmckl_get_jastrow_type_nucl_vector (qmckl_context context, int64_t* const type_nucl_num, const int64_t size_max);
qmckl_exit_code qmckl_get_jastrow_aord_vector (qmckl_context context, double * const aord_vector, int64_t* size_max); qmckl_exit_code qmckl_get_jastrow_aord_vector (qmckl_context context, double * const aord_vector, const int64_t size_max);
qmckl_exit_code qmckl_get_jastrow_bord_vector (qmckl_context context, double * const bord_vector, int64_t* size_max); qmckl_exit_code qmckl_get_jastrow_bord_vector (qmckl_context context, double * const bord_vector, const int64_t size_max);
qmckl_exit_code qmckl_get_jastrow_cord_vector (qmckl_context context, double * const cord_vector, int64_t* size_max); qmckl_exit_code qmckl_get_jastrow_cord_vector (qmckl_context context, double * const cord_vector, const int64_t size_max);
#+end_src #+end_src
Along with these core functions, calculation of the jastrow factor Along with these core functions, calculation of the jastrow factor
@ -474,7 +474,7 @@ qmckl_exit_code qmckl_get_jastrow_aord_num (const qmckl_context context, int64_t
} }
assert (ctx->jastrow.aord_num > 0); assert (ctx->jastrow.aord_num > 0);
*aord_num = ctx->jastrow.aord_num; ,*aord_num = ctx->jastrow.aord_num;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -501,7 +501,7 @@ qmckl_exit_code qmckl_get_jastrow_bord_num (const qmckl_context context, int64_t
} }
assert (ctx->jastrow.bord_num > 0); assert (ctx->jastrow.bord_num > 0);
*bord_num = ctx->jastrow.bord_num; ,*bord_num = ctx->jastrow.bord_num;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -528,7 +528,7 @@ qmckl_exit_code qmckl_get_jastrow_cord_num (const qmckl_context context, int64_t
} }
assert (ctx->jastrow.cord_num > 0); assert (ctx->jastrow.cord_num > 0);
*cord_num = ctx->jastrow.cord_num; ,*cord_num = ctx->jastrow.cord_num;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -555,11 +555,15 @@ qmckl_exit_code qmckl_get_jastrow_type_nucl_num (const qmckl_context context, in
} }
assert (ctx->jastrow.type_nucl_num > 0); assert (ctx->jastrow.type_nucl_num > 0);
*type_nucl_num = ctx->jastrow.type_nucl_num; ,*type_nucl_num = ctx->jastrow.type_nucl_num;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
qmckl_exit_code qmckl_get_jastrow_type_nucl_vector (const qmckl_context context, int64_t * const type_nucl_vector, int64_t* size_max) { qmckl_exit_code
qmckl_get_jastrow_type_nucl_vector (const qmckl_context context,
int64_t* const type_nucl_vector,
const int64_t size_max)
{
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return (char) 0; return (char) 0;
@ -582,12 +586,21 @@ qmckl_exit_code qmckl_get_jastrow_type_nucl_vector (const qmckl_context context,
} }
assert (ctx->jastrow.type_nucl_vector != NULL); assert (ctx->jastrow.type_nucl_vector != NULL);
if (size_max < ctx->jastrow.type_nucl_num) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_type_nucl_vector",
"Array too small. Expected jastrow.type_nucl_num");
}
memcpy(type_nucl_vector, ctx->jastrow.type_nucl_vector, ctx->jastrow.type_nucl_num*sizeof(int64_t)); memcpy(type_nucl_vector, ctx->jastrow.type_nucl_vector, ctx->jastrow.type_nucl_num*sizeof(int64_t));
(*size_max) = ctx->jastrow.type_nucl_num;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
qmckl_exit_code qmckl_get_jastrow_aord_vector (const qmckl_context context, double * const aord_vector, int64_t* size_max) { qmckl_exit_code
qmckl_get_jastrow_aord_vector (const qmckl_context context,
double * const aord_vector,
const int64_t size_max) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return (char) 0; return (char) 0;
@ -611,12 +624,20 @@ qmckl_exit_code qmckl_get_jastrow_aord_vector (const qmckl_context context, doub
assert (ctx->jastrow.aord_vector != NULL); assert (ctx->jastrow.aord_vector != NULL);
int64_t sze = (ctx->jastrow.aord_num + 1)*ctx->jastrow.type_nucl_num; int64_t sze = (ctx->jastrow.aord_num + 1)*ctx->jastrow.type_nucl_num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_aord_vector",
"Array too small. Expected (ctx->jastrow.aord_num + 1)*ctx->jastrow.type_nucl_num");
}
memcpy(aord_vector, ctx->jastrow.aord_vector, sze*sizeof(double)); memcpy(aord_vector, ctx->jastrow.aord_vector, sze*sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
qmckl_exit_code qmckl_get_jastrow_bord_vector (const qmckl_context context, double * const bord_vector, int64_t* size_max) { qmckl_exit_code
qmckl_get_jastrow_bord_vector (const qmckl_context context,
double * const bord_vector,
const int64_t size_max) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return (char) 0; return (char) 0;
@ -640,12 +661,20 @@ qmckl_exit_code qmckl_get_jastrow_bord_vector (const qmckl_context context, doub
assert (ctx->jastrow.bord_vector != NULL); assert (ctx->jastrow.bord_vector != NULL);
int64_t sze=ctx->jastrow.bord_num +1; int64_t sze=ctx->jastrow.bord_num +1;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_bord_vector",
"Array too small. Expected (ctx->jastrow.bord_num + 1)");
}
memcpy(bord_vector, ctx->jastrow.bord_vector, sze*sizeof(double)); memcpy(bord_vector, ctx->jastrow.bord_vector, sze*sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
qmckl_exit_code qmckl_get_jastrow_cord_vector (const qmckl_context context, double * const cord_vector, int64_t* size_max) { qmckl_exit_code
qmckl_get_jastrow_cord_vector (const qmckl_context context,
double * const cord_vector,
const int64_t size_max) {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return (char) 0; return (char) 0;
@ -668,14 +697,19 @@ qmckl_exit_code qmckl_get_jastrow_cord_vector (const qmckl_context context, doub
} }
assert (ctx->jastrow.cord_vector != NULL); assert (ctx->jastrow.cord_vector != NULL);
int64_t dim_cord_vect; int64_t dim_cord_vect;
qmckl_exit_code rc = qmckl_get_jastrow_dim_cord_vect(context, &dim_cord_vect); qmckl_exit_code rc = qmckl_get_jastrow_dim_cord_vect(context, &dim_cord_vect);
if (rc != QMCKL_SUCCESS) return rc; if (rc != QMCKL_SUCCESS) return rc;
int64_t sze=dim_cord_vect * ctx->jastrow.type_nucl_num; int64_t sze=dim_cord_vect * ctx->jastrow.type_nucl_num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_cord_vector",
"Array too small. Expected dim_cord_vect * jastrow.type_nucl_num");
}
memcpy(cord_vector, ctx->jastrow.cord_vector, sze*sizeof(double)); memcpy(cord_vector, ctx->jastrow.cord_vector, sze*sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -690,9 +724,9 @@ qmckl_exit_code qmckl_get_jastrow_cord_vector (const qmckl_context context, doub
qmckl_exit_code qmckl_set_jastrow_ord_num (qmckl_context context, const int64_t aord_num, const int64_t bord_num, const int64_t cord_num); qmckl_exit_code qmckl_set_jastrow_ord_num (qmckl_context context, const int64_t aord_num, const int64_t bord_num, const int64_t cord_num);
qmckl_exit_code qmckl_set_jastrow_type_nucl_num (qmckl_context context, const int64_t type_nucl_num); qmckl_exit_code qmckl_set_jastrow_type_nucl_num (qmckl_context context, const int64_t type_nucl_num);
qmckl_exit_code qmckl_set_jastrow_type_nucl_vector (qmckl_context context, const int64_t* type_nucl_vector, const int64_t nucl_num); qmckl_exit_code qmckl_set_jastrow_type_nucl_vector (qmckl_context context, const int64_t* type_nucl_vector, const int64_t nucl_num);
qmckl_exit_code qmckl_set_jastrow_aord_vector (qmckl_context context, const double * aord_vector, int64_t size_max); qmckl_exit_code qmckl_set_jastrow_aord_vector (qmckl_context context, const double * aord_vector, const int64_t size_max);
qmckl_exit_code qmckl_set_jastrow_bord_vector (qmckl_context context, const double * bord_vector, int64_t size_max); qmckl_exit_code qmckl_set_jastrow_bord_vector (qmckl_context context, const double * bord_vector, const int64_t size_max);
qmckl_exit_code qmckl_set_jastrow_cord_vector (qmckl_context context, const double * cord_vector, int64_t size_max); qmckl_exit_code qmckl_set_jastrow_cord_vector (qmckl_context context, const double * cord_vector, const int64_t size_max);
#+end_src #+end_src
#+NAME:pre2 #+NAME:pre2
@ -718,7 +752,12 @@ return QMCKL_SUCCESS;
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_set_jastrow_ord_num(qmckl_context context, const int64_t aord_num, const int64_t bord_num, const int64_t cord_num) { qmckl_exit_code
qmckl_set_jastrow_ord_num(qmckl_context context,
const int64_t aord_num,
const int64_t bord_num,
const int64_t cord_num)
{
<<pre2>> <<pre2>>
if (aord_num <= 0) { if (aord_num <= 0) {
@ -750,7 +789,10 @@ qmckl_exit_code qmckl_set_jastrow_ord_num(qmckl_context context, const int64_t a
<<post2>> <<post2>>
} }
qmckl_exit_code qmckl_set_jastrow_type_nucl_num(qmckl_context context, const int64_t type_nucl_num) {
qmckl_exit_code
qmckl_set_jastrow_type_nucl_num(qmckl_context context, const int64_t type_nucl_num)
{
<<pre2>> <<pre2>>
if (type_nucl_num <= 0) { if (type_nucl_num <= 0) {
@ -766,7 +808,12 @@ qmckl_exit_code qmckl_set_jastrow_type_nucl_num(qmckl_context context, const int
<<post2>> <<post2>>
} }
qmckl_exit_code qmckl_set_jastrow_type_nucl_vector(qmckl_context context, int64_t const * type_nucl_vector, const int64_t nucl_num) {
qmckl_exit_code
qmckl_set_jastrow_type_nucl_vector(qmckl_context context,
int64_t const * type_nucl_vector,
const int64_t nucl_num)
{
<<pre2>> <<pre2>>
int32_t mask = 1 << 2; int32_t mask = 1 << 2;
@ -816,7 +863,12 @@ qmckl_exit_code qmckl_set_jastrow_type_nucl_vector(qmckl_context context, int64_
<<post2>> <<post2>>
} }
qmckl_exit_code qmckl_set_jastrow_aord_vector(qmckl_context context, double const * aord_vector, int64_t size_max) {
qmckl_exit_code
qmckl_set_jastrow_aord_vector(qmckl_context context,
double const * aord_vector,
const int64_t size_max)
{
<<pre2>> <<pre2>>
int32_t mask = 1 << 3; int32_t mask = 1 << 3;
@ -849,13 +901,13 @@ qmckl_exit_code qmckl_set_jastrow_aord_vector(qmckl_context context, double cons
return qmckl_failwith( context, rc, return qmckl_failwith( context, rc,
"qmckl_set_ord_vector", "qmckl_set_ord_vector",
NULL); NULL);
} }
} }
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
mem_info.size = (aord_num + 1) * type_nucl_num * sizeof(double); mem_info.size = (aord_num + 1) * type_nucl_num * sizeof(double);
if (size_max < mem_info.size/sizeof(double)) { if ((size_t) size_max < mem_info.size/sizeof(double)) {
return qmckl_failwith( context, return qmckl_failwith( context,
QMCKL_INVALID_ARG_3, QMCKL_INVALID_ARG_3,
"qmckl_set_jastrow_aord_vector", "qmckl_set_jastrow_aord_vector",
@ -878,7 +930,12 @@ qmckl_exit_code qmckl_set_jastrow_aord_vector(qmckl_context context, double cons
<<post2>> <<post2>>
} }
qmckl_exit_code qmckl_set_jastrow_bord_vector(qmckl_context context, double const * bord_vector, int64_t size_max) {
qmckl_exit_code
qmckl_set_jastrow_bord_vector(qmckl_context context,
double const * bord_vector,
const int64_t size_max)
{
<<pre2>> <<pre2>>
int32_t mask = 1 << 4; int32_t mask = 1 << 4;
@ -913,7 +970,7 @@ qmckl_exit_code qmckl_set_jastrow_bord_vector(qmckl_context context, double cons
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
mem_info.size = (bord_num + 1) * sizeof(double); mem_info.size = (bord_num + 1) * sizeof(double);
if (size_max < mem_info.size/sizeof(double)) { if ((size_t) size_max < mem_info.size/sizeof(double)) {
return qmckl_failwith( context, return qmckl_failwith( context,
QMCKL_INVALID_ARG_3, QMCKL_INVALID_ARG_3,
"qmckl_set_jastrow_bord_vector", "qmckl_set_jastrow_bord_vector",
@ -936,7 +993,12 @@ qmckl_exit_code qmckl_set_jastrow_bord_vector(qmckl_context context, double cons
<<post2>> <<post2>>
} }
qmckl_exit_code qmckl_set_jastrow_cord_vector(qmckl_context context, double const * cord_vector, int64_t size_max) {
qmckl_exit_code
qmckl_set_jastrow_cord_vector(qmckl_context context,
double const * cord_vector,
const int64_t size_max)
{
<<pre2>> <<pre2>>
int32_t mask = 1 << 5; int32_t mask = 1 << 5;
@ -978,7 +1040,7 @@ qmckl_exit_code qmckl_set_jastrow_cord_vector(qmckl_context context, double cons
qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero;
mem_info.size = dim_cord_vect * type_nucl_num * sizeof(double); mem_info.size = dim_cord_vect * type_nucl_num * sizeof(double);
if (size_max < mem_info.size/sizeof(double)) { if ((size_t) size_max < mem_info.size/sizeof(double)) {
return qmckl_failwith( context, return qmckl_failwith( context,
QMCKL_INVALID_ARG_3, QMCKL_INVALID_ARG_3,
"qmckl_set_jastrow_cord_vector", "qmckl_set_jastrow_cord_vector",
@ -1069,6 +1131,7 @@ double* elec_coord = &(n2_elec_coord[0][0][0]);
const double* nucl_charge = n2_charge; const double* nucl_charge = n2_charge;
int64_t nucl_num = n2_nucl_num; int64_t nucl_num = n2_nucl_num;
double* nucl_coord = &(n2_nucl_coord[0][0]); double* nucl_coord = &(n2_nucl_coord[0][0]);
int64_t size_max;
/* Provide Electron data */ /* Provide Electron data */
@ -1246,11 +1309,17 @@ assert(qmckl_nucleus_provided(context));
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_asymp_jasb(qmckl_context context, double* const asymp_jasb, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_asymp_jasb(qmckl_context context,
double* const asymp_jasb,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_asymp_jasb(qmckl_context context, double* const asymp_jasb, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_asymp_jasb(qmckl_context context,
double* const asymp_jasb,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -1264,9 +1333,14 @@ qmckl_exit_code qmckl_get_jastrow_asymp_jasb(qmckl_context context, double* cons
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
size_t sze = 2; int64_t sze = 2;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_asymp_jasb",
"Array too small. Expected 2");
}
memcpy(asymp_jasb, ctx->jastrow.asymp_jasb, sze * sizeof(double)); memcpy(asymp_jasb, ctx->jastrow.asymp_jasb, sze * sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -1409,7 +1483,7 @@ qmckl_exit_code qmckl_compute_asymp_jasb (
} }
asym_one = bord_vector[0] * kappa_inv / (1.0 + bord_vector[1] * kappa_inv); asym_one = bord_vector[0] * kappa_inv / (1.0 + bord_vector[1] * kappa_inv);
asymp_jasb[0] = asym_one; asymp_jasb[0] = asym_one;
asymp_jasb[1] = 0.5 * asym_one; asymp_jasb[1] = 0.5 * asym_one;
for (int i = 0 ; i <= 1; ++i) { for (int i = 0 ; i <= 1; ++i) {
@ -1418,7 +1492,7 @@ qmckl_exit_code qmckl_compute_asymp_jasb (
x = x * kappa_inv; x = x * kappa_inv;
asymp_jasb[i] = asymp_jasb[i] + bord_vector[p + 1] * x; asymp_jasb[i] = asymp_jasb[i] + bord_vector[p + 1] * x;
} }
} }
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -1433,7 +1507,7 @@ qmckl_exit_code qmckl_compute_asymp_jasb (
const int64_t bord_num, const int64_t bord_num,
const double* bord_vector, const double* bord_vector,
const double rescale_factor_kappa_ee, const double rescale_factor_kappa_ee,
double* const asymp_jasb ); double* const asymp_jasb );
#+end_src #+end_src
@ -1457,7 +1531,7 @@ print("asym_one : ", asym_one)
print("asymp_jasb[0] : ", asymp_jasb[0]) print("asymp_jasb[0] : ", asymp_jasb[0])
print("asymp_jasb[1] : ", asymp_jasb[1]) print("asymp_jasb[1] : ", asymp_jasb[1])
#+end_src #+end_src
#+RESULTS: asymp_jasb #+RESULTS: asymp_jasb
: asym_one : 0.43340325572525706 : asym_one : 0.43340325572525706
: asymp_jasb[0] : 0.5323750557252571 : asymp_jasb[0] : 0.5323750557252571
@ -1500,8 +1574,7 @@ assert(rc == QMCKL_SUCCESS);
assert(qmckl_jastrow_provided(context)); assert(qmckl_jastrow_provided(context));
double asymp_jasb[2]; double asymp_jasb[2];
int64_t size_max=0; rc = qmckl_get_jastrow_asymp_jasb(context, asymp_jasb,2);
rc = qmckl_get_jastrow_asymp_jasb(context, asymp_jasb,&size_max);
// calculate asymp_jasb // calculate asymp_jasb
assert(fabs(asymp_jasb[0]-0.5323750557252571) < 1.e-12); assert(fabs(asymp_jasb[0]-0.5323750557252571) < 1.e-12);
@ -1521,11 +1594,17 @@ f_{ee} = \sum_{i,j<i} \left\{ \frac{ \eta B_0 C_{ij}}{1 - B_1 C_{ij}} - J_{asym
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_factor_ee(qmckl_context context, double* const factor_ee, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_factor_ee(qmckl_context context,
double* const factor_ee,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_factor_ee(qmckl_context context, double* const factor_ee, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_factor_ee(qmckl_context context,
double* const factor_ee,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -1540,8 +1619,13 @@ qmckl_exit_code qmckl_get_jastrow_factor_ee(qmckl_context context, double* const
assert (ctx != NULL); assert (ctx != NULL);
int64_t sze=ctx->electron.walk_num; int64_t sze=ctx->electron.walk_num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_ee",
"Array too small. Expected walk_num");
}
memcpy(factor_ee, ctx->jastrow.factor_ee, sze*sizeof(double)); memcpy(factor_ee, ctx->jastrow.factor_ee, sze*sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -1828,8 +1912,7 @@ print("factor_ee :",factor_ee)
assert(qmckl_jastrow_provided(context)); assert(qmckl_jastrow_provided(context));
double factor_ee[walk_num]; double factor_ee[walk_num];
size_max=0; rc = qmckl_get_jastrow_factor_ee(context, factor_ee, walk_num);
rc = qmckl_get_jastrow_factor_ee(context, factor_ee, &size_max);
// calculate factor_ee // calculate factor_ee
assert(fabs(factor_ee[0]+4.282760865958113) < 1.e-12); assert(fabs(factor_ee[0]+4.282760865958113) < 1.e-12);
@ -1848,11 +1931,17 @@ assert(fabs(factor_ee[0]+4.282760865958113) < 1.e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_factor_ee_deriv_e(qmckl_context context, double* const factor_ee_deriv_e, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_factor_ee_deriv_e(qmckl_context context,
double* const factor_ee_deriv_e,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_factor_ee_deriv_e(qmckl_context context, double* const factor_ee_deriv_e, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_factor_ee_deriv_e(qmckl_context context,
double* const factor_ee_deriv_e,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -1867,8 +1956,14 @@ qmckl_exit_code qmckl_get_jastrow_factor_ee_deriv_e(qmckl_context context, doubl
assert (ctx != NULL); assert (ctx != NULL);
int64_t sze = ctx->electron.walk_num * 4 * ctx->electron.num; int64_t sze = ctx->electron.walk_num * 4 * ctx->electron.num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_ee_deriv_e",
"Array too small. Expected 4*walk_num*elec_num");
}
memcpy(factor_ee_deriv_e, ctx->jastrow.factor_ee_deriv_e, sze * sizeof(double)); memcpy(factor_ee_deriv_e, ctx->jastrow.factor_ee_deriv_e, sze * sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -2026,7 +2121,7 @@ integer function qmckl_compute_factor_ee_deriv_e_f(context, walk_num, elec_num,
dx(1) = ee_distance_rescaled_deriv_e(1, i, j, nw) dx(1) = ee_distance_rescaled_deriv_e(1, i, j, nw)
dx(2) = ee_distance_rescaled_deriv_e(2, i, j, nw) dx(2) = ee_distance_rescaled_deriv_e(2, i, j, nw)
dx(3) = ee_distance_rescaled_deriv_e(3, i, j, nw) dx(3) = ee_distance_rescaled_deriv_e(3, i, j, nw)
dx(4) = ee_distance_rescaled_deriv_e(4, i, j, nw) dx(4) = ee_distance_rescaled_deriv_e(4, i, j, nw)
if((i .LE. up_num .AND. j .LE. up_num ) .OR. & if((i .LE. up_num .AND. j .LE. up_num ) .OR. &
(i .GT. up_num .AND. j .GT. up_num)) then (i .GT. up_num .AND. j .GT. up_num)) then
@ -2241,8 +2336,7 @@ assert(qmckl_jastrow_provided(context));
// calculate factor_ee_deriv_e // calculate factor_ee_deriv_e
double factor_ee_deriv_e[walk_num][4][elec_num]; double factor_ee_deriv_e[walk_num][4][elec_num];
size_max=0; rc = qmckl_get_jastrow_factor_ee_deriv_e(context, &(factor_ee_deriv_e[0][0][0]),walk_num*4*elec_num);
rc = qmckl_get_jastrow_factor_ee_deriv_e(context, &(factor_ee_deriv_e[0][0][0]),&size_max);
// check factor_ee_deriv_e // check factor_ee_deriv_e
assert(fabs(factor_ee_deriv_e[0][0][0]-0.16364894652107934) < 1.e-12); assert(fabs(factor_ee_deriv_e[0][0][0]-0.16364894652107934) < 1.e-12);
@ -2264,11 +2358,17 @@ f_{en} = \sum_{i,j<i} \left\{ \frac{ A_0 C_{ij}}{1 - A_1 C_{ij}} + \sum^{nord}_{
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_factor_en(qmckl_context context, double* const factor_en, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_factor_en(qmckl_context context,
double* const factor_en,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_factor_en(qmckl_context context, double* const factor_en, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_factor_en(qmckl_context context,
double* const factor_en,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -2283,8 +2383,13 @@ qmckl_exit_code qmckl_get_jastrow_factor_en(qmckl_context context, double* const
assert (ctx != NULL); assert (ctx != NULL);
int64_t sze=ctx->electron.walk_num; int64_t sze=ctx->electron.walk_num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_en",
"Array too small. Expected walk_num");
}
memcpy(factor_en, ctx->jastrow.factor_en, sze*sizeof(double)); memcpy(factor_en, ctx->jastrow.factor_en, sze*sizeof(double));
(*size_max)=sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -2562,8 +2667,7 @@ print("factor_en :",factor_en)
assert(qmckl_jastrow_provided(context)); assert(qmckl_jastrow_provided(context));
double factor_en[walk_num]; double factor_en[walk_num];
size_max=0; rc = qmckl_get_jastrow_factor_en(context, factor_en,walk_num);
rc = qmckl_get_jastrow_factor_en(context, factor_en,&size_max);
// calculate factor_en // calculate factor_en
assert(fabs(factor_en[0]+5.865822569188727) < 1.e-12); assert(fabs(factor_en[0]+5.865822569188727) < 1.e-12);
@ -2579,11 +2683,17 @@ assert(fabs(factor_en[0]+5.865822569188727) < 1.e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context, double* const factor_en_deriv_e, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context,
double* const factor_en_deriv_e,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context, double* const factor_en_deriv_e, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context,
double* const factor_en_deriv_e,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -2598,8 +2708,13 @@ qmckl_exit_code qmckl_get_jastrow_factor_en_deriv_e(qmckl_context context, doubl
assert (ctx != NULL); assert (ctx != NULL);
int64_t sze = ctx->electron.walk_num * 4 * ctx->electron.num; int64_t sze = ctx->electron.walk_num * 4 * ctx->electron.num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_en_deriv_e",
"Array too small. Expected 4*walk_num*elec_num");
}
memcpy(factor_en_deriv_e, ctx->jastrow.factor_en_deriv_e, sze*sizeof(double)); memcpy(factor_en_deriv_e, ctx->jastrow.factor_en_deriv_e, sze*sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -2966,8 +3081,7 @@ assert(qmckl_jastrow_provided(context));
// calculate factor_en_deriv_e // calculate factor_en_deriv_e
double factor_en_deriv_e[walk_num][4][elec_num]; double factor_en_deriv_e[walk_num][4][elec_num];
size_max=0; rc = qmckl_get_jastrow_factor_en_deriv_e(context, &(factor_en_deriv_e[0][0][0]),walk_num*4*elec_num);
rc = qmckl_get_jastrow_factor_en_deriv_e(context, &(factor_en_deriv_e[0][0][0]),&size_max);
// check factor_en_deriv_e // check factor_en_deriv_e
assert(fabs(factor_en_deriv_e[0][0][0]-0.11609919541763383) < 1.e-12); assert(fabs(factor_en_deriv_e[0][0][0]-0.11609919541763383) < 1.e-12);
@ -2991,11 +3105,17 @@ assert(fabs(factor_en_deriv_e[0][3][0]+0.9667363412285741 ) < 1.e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_een_rescaled_e(qmckl_context context, double* const distance_rescaled, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_een_rescaled_e(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_een_rescaled_e(qmckl_context context, double* const distance_rescaled, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_een_rescaled_e(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -3009,9 +3129,14 @@ qmckl_exit_code qmckl_get_jastrow_een_rescaled_e(qmckl_context context, double*
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
size_t sze = ctx->electron.num * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1); int64_t sze = ctx->electron.num * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1);
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_een_rescaled_e",
"Array too small. Expected ctx->electron.num * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1)");
}
memcpy(distance_rescaled, ctx->jastrow.een_rescaled_e, sze * sizeof(double)); memcpy(distance_rescaled, ctx->jastrow.een_rescaled_e, sze * sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -3223,7 +3348,7 @@ end function qmckl_compute_een_rescaled_e_f
#+end_src #+end_src
*** Test *** Test
#+begin_src python :results output :exports none :noweb yes #+begin_src python :results output :exports none :noweb yes
import numpy as np import numpy as np
@ -3287,8 +3412,7 @@ assert(qmckl_electron_provided(context));
double een_rescaled_e[walk_num][(cord_num + 1)][elec_num][elec_num]; double een_rescaled_e[walk_num][(cord_num + 1)][elec_num][elec_num];
size_max=0; rc = qmckl_get_jastrow_een_rescaled_e(context, &(een_rescaled_e[0][0][0][0]),elec_num*elec_num*(cord_num+1)*walk_num);
rc = qmckl_get_jastrow_een_rescaled_e(context, &(een_rescaled_e[0][0][0][0]),&size_max);
// value of (0,2,1) // value of (0,2,1)
assert(fabs(een_rescaled_e[0][1][0][2]-0.08084493981483197) < 1.e-12); assert(fabs(een_rescaled_e[0][1][0][2]-0.08084493981483197) < 1.e-12);
@ -3313,11 +3437,17 @@ assert(fabs(een_rescaled_e[0][2][1][5]-0.3424402276009091) < 1.e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_een_rescaled_e_deriv_e(qmckl_context context, double* const distance_rescaled, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_een_rescaled_e_deriv_e(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_een_rescaled_e_deriv_e(qmckl_context context, double* const distance_rescaled, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_een_rescaled_e_deriv_e(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -3331,9 +3461,14 @@ qmckl_exit_code qmckl_get_jastrow_een_rescaled_e_deriv_e(qmckl_context context,
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
size_t sze = ctx->electron.num * 4 * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1); int64_t sze = ctx->electron.num * 4 * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1);
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_een_deriv_e",
"Array too small. Expected ctx->electron.num * 4 * ctx->electron.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1)");
}
memcpy(distance_rescaled, ctx->jastrow.een_rescaled_e_deriv_e, sze * sizeof(double)); memcpy(distance_rescaled, ctx->jastrow.een_rescaled_e_deriv_e, sze * sizeof(double));
(*size_max) = sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -3654,8 +3789,9 @@ for l in range(0,cord_num+1):
#+begin_src c :tangle (eval c_test) #+begin_src c :tangle (eval c_test)
//assert(qmckl_electron_provided(context)); //assert(qmckl_electron_provided(context));
double een_rescaled_e_deriv_e[walk_num][(cord_num + 1)][elec_num][4][elec_num]; double een_rescaled_e_deriv_e[walk_num][(cord_num + 1)][elec_num][4][elec_num];
size_max=0; size_max=walk_num*(cord_num + 1)*elec_num*4*elec_num;
rc = qmckl_get_jastrow_een_rescaled_e_deriv_e(context, &(een_rescaled_e_deriv_e[0][0][0][0][0]),&size_max); rc = qmckl_get_jastrow_een_rescaled_e_deriv_e(context,
&(een_rescaled_e_deriv_e[0][0][0][0][0]),size_max);
// value of (0,0,0,2,1) // value of (0,0,0,2,1)
assert(fabs(een_rescaled_e_deriv_e[0][1][0][0][2] + 0.05991352796887283 ) < 1.e-12); assert(fabs(een_rescaled_e_deriv_e[0][1][0][0][2] + 0.05991352796887283 ) < 1.e-12);
@ -3680,11 +3816,17 @@ assert(fabs(een_rescaled_e_deriv_e[0][2][1][0][5] + 0.5880599146214673 ) < 1.
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_een_rescaled_n(qmckl_context context, double* const distance_rescaled, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_een_rescaled_n(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_een_rescaled_n(qmckl_context context, double* const distance_rescaled, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_een_rescaled_n(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -3698,9 +3840,14 @@ qmckl_exit_code qmckl_get_jastrow_een_rescaled_n(qmckl_context context, double*
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
size_t sze = ctx->electron.num * ctx->nucleus.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1); int64_t sze = ctx->electron.num * ctx->nucleus.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1);
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_een_deriv_e",
"Array too small. Expected ctx->electron.num * ctx->nucleus.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1)");
}
memcpy(distance_rescaled, ctx->jastrow.een_rescaled_n, sze * sizeof(double)); memcpy(distance_rescaled, ctx->jastrow.een_rescaled_n, sze * sizeof(double));
(*size_max)=sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -3982,8 +4129,8 @@ print(" een_rescaled_n[1, 5, 2] = ",een_rescaled_n[1, 5, 2])
assert(qmckl_electron_provided(context)); assert(qmckl_electron_provided(context));
double een_rescaled_n[walk_num][(cord_num + 1)][nucl_num][elec_num]; double een_rescaled_n[walk_num][(cord_num + 1)][nucl_num][elec_num];
size_max=0; size_max=walk_num*(cord_num + 1)*nucl_num*elec_num;
rc = qmckl_get_jastrow_een_rescaled_n(context, &(een_rescaled_n[0][0][0][0]),&size_max); rc = qmckl_get_jastrow_een_rescaled_n(context, &(een_rescaled_n[0][0][0][0]),size_max);
// value of (0,2,1) // value of (0,2,1)
assert(fabs(een_rescaled_n[0][1][0][2]-0.10612983920006765) < 1.e-12); assert(fabs(een_rescaled_n[0][1][0][2]-0.10612983920006765) < 1.e-12);
@ -4003,11 +4150,17 @@ assert(fabs(een_rescaled_n[0][2][1][5]-0.01343938025140174) < 1.e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_een_rescaled_n_deriv_e(qmckl_context context, double* const distance_rescaled, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_een_rescaled_n_deriv_e(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_een_rescaled_n_deriv_e(qmckl_context context, double* const distance_rescaled, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_een_rescaled_n_deriv_e(qmckl_context context,
double* const distance_rescaled,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -4021,9 +4174,14 @@ qmckl_exit_code qmckl_get_jastrow_een_rescaled_n_deriv_e(qmckl_context context,
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
size_t sze = ctx->electron.num * 4 * ctx->nucleus.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1); int64_t sze = ctx->electron.num * 4 * ctx->nucleus.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1);
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_een_deriv_e",
"Array too small. Expected ctx->electron.num * 4 * ctx->nucleus.num * ctx->electron.walk_num * (ctx->jastrow.cord_num + 1)");
}
memcpy(distance_rescaled, ctx->jastrow.een_rescaled_n_deriv_e, sze * sizeof(double)); memcpy(distance_rescaled, ctx->jastrow.een_rescaled_n_deriv_e, sze * sizeof(double));
(*size_max)=sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -4282,9 +4440,9 @@ end function qmckl_compute_factor_een_rescaled_n_deriv_e_f
end function qmckl_compute_factor_een_rescaled_n_deriv_e end function qmckl_compute_factor_een_rescaled_n_deriv_e
#+end_src #+end_src
*** Test *** Test
#+begin_src python :results output :exports none :noweb yes #+begin_src python :results output :exports none :noweb yes
import numpy as np import numpy as np
@ -4354,8 +4512,8 @@ print(" een_rescaled_n_deriv_e[2, 1, 6, 2] = ",een_rescaled_n_deriv_e[5, 0, 1, 2
assert(qmckl_electron_provided(context)); assert(qmckl_electron_provided(context));
double een_rescaled_n_deriv_e[walk_num][(cord_num + 1)][nucl_num][4][elec_num]; double een_rescaled_n_deriv_e[walk_num][(cord_num + 1)][nucl_num][4][elec_num];
size_max=0; size_max=walk_num*(cord_num + 1)*nucl_num*4*elec_num;
rc = qmckl_get_jastrow_een_rescaled_n_deriv_e(context, &(een_rescaled_n_deriv_e[0][0][0][0][0]),&size_max); rc = qmckl_get_jastrow_een_rescaled_n_deriv_e(context, &(een_rescaled_n_deriv_e[0][0][0][0][0]),size_max);
// value of (0,2,1) // value of (0,2,1)
assert(fabs(een_rescaled_n_deriv_e[0][1][0][0][2]+0.07633444246999128 ) < 1.e-12); assert(fabs(een_rescaled_n_deriv_e[0][1][0][0][2]+0.07633444246999128 ) < 1.e-12);
@ -5397,11 +5555,17 @@ assert(fabs(dtmp_c[0][1][0][0][0][0] - 0.237440520852232) < 1e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_factor_een(qmckl_context context, double* const factor_een, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_factor_een(qmckl_context context,
double* const factor_een,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_factor_een(qmckl_context context, double* const factor_een, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_factor_een(qmckl_context context,
double* const factor_een,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -5415,9 +5579,14 @@ qmckl_exit_code qmckl_get_jastrow_factor_een(qmckl_context context, double* cons
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
int64_t sze = ctx->electron.walk_num * ctx->electron.num; int64_t sze = ctx->electron.walk_num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_een",
"Array too small. Expected walk_num");
}
memcpy(factor_een, ctx->jastrow.factor_een, sze*sizeof(double)); memcpy(factor_een, ctx->jastrow.factor_een, sze*sizeof(double));
(*size_max)=sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -5871,8 +6040,7 @@ print("factor_een:",factor_een)
assert(qmckl_jastrow_provided(context)); assert(qmckl_jastrow_provided(context));
double factor_een[walk_num]; double factor_een[walk_num];
size_max=0; rc = qmckl_get_jastrow_factor_een(context, &(factor_een[0]),walk_num);
rc = qmckl_get_jastrow_factor_een(context, &(factor_een[0]),&size_max);
assert(fabs(factor_een[0] + 0.37407972141304213) < 1e-12); assert(fabs(factor_een[0] + 0.37407972141304213) < 1e-12);
#+end_src #+end_src
@ -5886,11 +6054,17 @@ assert(fabs(factor_een[0] + 0.37407972141304213) < 1e-12);
*** Get *** Get
#+begin_src c :comments org :tangle (eval h_func) :noweb yes #+begin_src c :comments org :tangle (eval h_func) :noweb yes
qmckl_exit_code qmckl_get_jastrow_factor_een_deriv_e(qmckl_context context, double* const factor_een_deriv_e, int64_t* size_max); qmckl_exit_code
qmckl_get_jastrow_factor_een_deriv_e(qmckl_context context,
double* const factor_een_deriv_e,
const int64_t size_max);
#+end_src #+end_src
#+begin_src c :comments org :tangle (eval c) :noweb yes :exports none #+begin_src c :comments org :tangle (eval c) :noweb yes :exports none
qmckl_exit_code qmckl_get_jastrow_factor_een_deriv_e(qmckl_context context, double* const factor_een_deriv_e, int64_t* size_max) qmckl_exit_code
qmckl_get_jastrow_factor_een_deriv_e(qmckl_context context,
double* const factor_een_deriv_e,
const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return QMCKL_NULL_CONTEXT; return QMCKL_NULL_CONTEXT;
@ -5904,9 +6078,14 @@ qmckl_exit_code qmckl_get_jastrow_factor_een_deriv_e(qmckl_context context, doub
qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; qmckl_context_struct* const ctx = (qmckl_context_struct* const) context;
assert (ctx != NULL); assert (ctx != NULL);
int64_t sze = ctx->electron.walk_num * ctx->electron.num; int64_t sze = ctx->electron.walk_num * 4 * ctx->electron.num;
if (size_max < sze) {
return qmckl_failwith( context,
QMCKL_INVALID_ARG_3,
"qmckl_get_jastrow_factor_een_deriv_e",
"Array too small. Expected 4*walk_num*elec_num");
}
memcpy(factor_een_deriv_e, ctx->jastrow.factor_een_deriv_e, sze*sizeof(double)); memcpy(factor_een_deriv_e, ctx->jastrow.factor_een_deriv_e, sze*sizeof(double));
(*size_max)=sze;
return QMCKL_SUCCESS; return QMCKL_SUCCESS;
} }
@ -6394,7 +6573,7 @@ end function qmckl_compute_factor_een_deriv_e_f
end function qmckl_compute_factor_een_deriv_e end function qmckl_compute_factor_een_deriv_e
#+end_src #+end_src
*** Test *** Test
#+begin_src python :results output :exports none :noweb yes #+begin_src python :results output :exports none :noweb yes
import numpy as np import numpy as np
@ -6451,11 +6630,10 @@ print("factor_een:",factor_een)
/* Check if Jastrow is properly initialized */ /* Check if Jastrow is properly initialized */
assert(qmckl_jastrow_provided(context)); assert(qmckl_jastrow_provided(context));
double factor_een_deriv_e[walk_num][elec_num]; double factor_een_deriv_e[4][walk_num][elec_num];
size_max=0; rc = qmckl_get_jastrow_factor_een_deriv_e(context, &(factor_een_deriv_e[0][0][0]),4*walk_num*elec_num);
rc = qmckl_get_jastrow_factor_een_deriv_e(context, &(factor_een_deriv_e[0][0]),&size_max);
assert(fabs(factor_een_deriv_e[0][0] + 0.0005481671107226865) < 1e-12); assert(fabs(factor_een_deriv_e[0][0][0] + 0.0005481671107226865) < 1e-12);
#+end_src #+end_src
* End of files :noexport: * End of files :noexport:

2
org/qmckl_mo.org
View File

@ -104,7 +104,7 @@ typedef struct qmckl_mo_basis_struct {
double * coefficient; double * coefficient;
double * mo_vgl; double * mo_vgl;
int64_t mo_vgl_date; uint64_t mo_vgl_date;
int32_t uninitialized; int32_t uninitialized;
bool provided; bool provided;

16
org/qmckl_nucleus.org
View File

@ -333,7 +333,7 @@ qmckl_get_nucleus_coord (const qmckl_context context,
rc = qmckl_transpose(context, ctx->nucleus.coord, At); rc = qmckl_transpose(context, ctx->nucleus.coord, At);
if (rc != QMCKL_SUCCESS) return rc; if (rc != QMCKL_SUCCESS) return rc;
rc = qmckl_double_of_matrix(context, At, coord, size_max); rc = qmckl_double_of_matrix(context, At, coord, size_max);
qmckl_matrix_free(context, At); qmckl_matrix_free(context, &At);
} else { } else {
rc = qmckl_double_of_matrix(context, ctx->nucleus.coord, coord, size_max); rc = qmckl_double_of_matrix(context, ctx->nucleus.coord, coord, size_max);
} }
@ -470,6 +470,12 @@ qmckl_set_nucleus_charge(qmckl_context context,
ctx->nucleus.charge = qmckl_vector_alloc(context, num); ctx->nucleus.charge = qmckl_vector_alloc(context, num);
rc = qmckl_vector_of_double(context, charge, num, &(ctx->nucleus.charge)); rc = qmckl_vector_of_double(context, charge, num, &(ctx->nucleus.charge));
if (rc != QMCKL_SUCCESS) {
return qmckl_failwith( context,
QMCKL_FAILURE,
"qmckl_set_nucleus_charge",
"Error in vector->double* conversion");
}
<<post2>> <<post2>>
} }
@ -518,7 +524,7 @@ qmckl_set_nucleus_coord(qmckl_context context,
const int64_t nucl_num = (int64_t) ctx->nucleus.num; const int64_t nucl_num = (int64_t) ctx->nucleus.num;
if (ctx->nucleus.coord.data != NULL) { if (ctx->nucleus.coord.data != NULL) {
rc = qmckl_matrix_free(context, ctx->nucleus.coord); rc = qmckl_matrix_free(context, &(ctx->nucleus.coord));
if (rc != QMCKL_SUCCESS) return rc; if (rc != QMCKL_SUCCESS) return rc;
} }
@ -662,14 +668,14 @@ assert(!qmckl_nucleus_provided(context));
rc = qmckl_get_nucleus_coord (context, 'N', nucl_coord2, 3*nucl_num); rc = qmckl_get_nucleus_coord (context, 'N', nucl_coord2, 3*nucl_num);
assert(rc == QMCKL_SUCCESS); assert(rc == QMCKL_SUCCESS);
for (size_t k=0 ; k<3 ; ++k) { for (size_t k=0 ; k<3 ; ++k) {
for (size_t i=0 ; i<nucl_num ; ++i) { for (int64_t i=0 ; i<nucl_num ; ++i) {
assert( nucl_coord[nucl_num*k+i] == nucl_coord2[3*i+k] ); assert( nucl_coord[nucl_num*k+i] == nucl_coord2[3*i+k] );
} }
} }
rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*nucl_num); rc = qmckl_get_nucleus_coord (context, 'T', nucl_coord2, 3*nucl_num);
assert(rc == QMCKL_SUCCESS); assert(rc == QMCKL_SUCCESS);
for (size_t i=0 ; i<3*nucl_num ; ++i) { for (int64_t i=0 ; i<3*nucl_num ; ++i) {
assert( nucl_coord[i] == nucl_coord2[i] ); assert( nucl_coord[i] == nucl_coord2[i] );
} }
@ -683,7 +689,7 @@ assert(rc == QMCKL_SUCCESS);
rc = qmckl_get_nucleus_charge(context, nucl_charge2, nucl_num); rc = qmckl_get_nucleus_charge(context, nucl_charge2, nucl_num);
assert(rc == QMCKL_SUCCESS); assert(rc == QMCKL_SUCCESS);
for (size_t i=0 ; i<nucl_num ; ++i) { for (int64_t i=0 ; i<nucl_num ; ++i) {
assert( nucl_charge[i] == nucl_charge2[i] ); assert( nucl_charge[i] == nucl_charge2[i] );
} }
assert(qmckl_nucleus_provided(context)); assert(qmckl_nucleus_provided(context));

19
org/qmckl_point.org
View File

@ -224,7 +224,7 @@ qmckl_get_point(const qmckl_context context,
if (rc != QMCKL_SUCCESS) return rc; if (rc != QMCKL_SUCCESS) return rc;
rc = qmckl_double_of_matrix( context, At, coord, size_max); rc = qmckl_double_of_matrix( context, At, coord, size_max);
if (rc != QMCKL_SUCCESS) return rc; if (rc != QMCKL_SUCCESS) return rc;
rc = qmckl_matrix_free(context, At); rc = qmckl_matrix_free(context, &At);
} else { } else {
rc = qmckl_double_of_matrix( context, ctx->point.coord, coord, size_max); rc = qmckl_double_of_matrix( context, ctx->point.coord, coord, size_max);
} }
@ -301,7 +301,7 @@ qmckl_set_point (qmckl_context context,
if (ctx->point.num < num) { if (ctx->point.num < num) {
if (ctx->point.coord.data != NULL) { if (ctx->point.coord.data != NULL) {
rc = qmckl_matrix_free(context, ctx->point.coord); rc = qmckl_matrix_free(context, &(ctx->point.coord));
assert (rc == QMCKL_SUCCESS); assert (rc == QMCKL_SUCCESS);
} }
@ -318,8 +318,17 @@ qmckl_set_point (qmckl_context context,
ctx->point.num = num; ctx->point.num = num;
if (transp == 'T') { if (transp == 'T') {
memcpy(ctx->point.coord.data, coord, 3*num*sizeof(double)); double *a = ctx->point.coord.data;
#ifdef HAVE_OPENMP
#pragma omp for
#endif
for (int64_t i=0 ; i<3*num ; ++i) {
a[i] = coord[i];
}
} else { } else {
#ifdef HAVE_OPENMP
#pragma omp for
#endif
for (int64_t i=0 ; i<num ; ++i) { for (int64_t i=0 ; i<num ; ++i) {
qmckl_mat(ctx->point.coord, i, 0) = coord[3*i ]; qmckl_mat(ctx->point.coord, i, 0) = coord[3*i ];
qmckl_mat(ctx->point.coord, i, 1) = coord[3*i+1]; qmckl_mat(ctx->point.coord, i, 1) = coord[3*i+1];
@ -328,8 +337,8 @@ qmckl_set_point (qmckl_context context,
} }
/* Increment the date of the context */ /* Increment the date of the context */
ctx->date += 1UL; rc = qmckl_context_touch(context);
ctx->point.date = ctx->date; assert (rc == QMCKL_SUCCESS);
return QMCKL_SUCCESS; return QMCKL_SUCCESS;

36
org/qmckl_trexio.org
View File

@ -429,6 +429,11 @@ qmckl_trexio_read_ao_X(qmckl_context context, trexio_t* const file)
/* Reformat data */ /* Reformat data */
rc = qmckl_set_ao_basis_nucleus_index(context, nucleus_index, nucleus_num); rc = qmckl_set_ao_basis_nucleus_index(context, nucleus_index, nucleus_num);
if (rc != QMCKL_SUCCESS) {
qmckl_free(context, nucleus_index);
nucleus_index = NULL;
return rc;
}
for (int i=shell_num-1 ; i>=0 ; --i) { for (int i=shell_num-1 ; i>=0 ; --i) {
const int k = tmp_array[i]; const int k = tmp_array[i];
@ -1048,21 +1053,44 @@ qmckl_trexio_read_mo_X(qmckl_context context, trexio_t* const file)
* Read everything * Read everything
#+begin_src c :tangle (eval h_func) #+begin_src c :tangle (eval h_func)
qmckl_exit_code qmckl_trexio_read(const qmckl_context context, const char* file_name); qmckl_exit_code
qmckl_trexio_read(const qmckl_context context,
const char* file_name,
const int64_t size_max);
#+end_src
#+begin_src f90 :tangle (eval fh_func) :comments org :exports none
interface
integer(c_int32_t) function qmckl_trexio_read &
(context, file_name, size_max) &
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 :: size_max
character(c_char ) , intent(in) :: file_name(size_max)
end function qmckl_trexio_read
end interface
#+end_src #+end_src
#+begin_src c :tangle (eval c) #+begin_src c :tangle (eval c)
qmckl_exit_code qmckl_exit_code
qmckl_trexio_read(const qmckl_context context, const char* file_name) qmckl_trexio_read(const qmckl_context context, const char* file_name, const int64_t size_max)
{ {
if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) { if (qmckl_context_check(context) == QMCKL_NULL_CONTEXT) {
return false; return false;
} }
qmckl_exit_code rc; qmckl_exit_code rc;
char file_name_new[size_max+1];
strncpy(file_name_new, file_name, size_max+1);
file_name_new[size_max] = '\0';
#ifdef HAVE_TREXIO #ifdef HAVE_TREXIO
trexio_t* file = qmckl_trexio_open_X(file_name, &rc); trexio_t* file = qmckl_trexio_open_X(file_name_new, &rc);
if (file == NULL) { if (file == NULL) {
trexio_close(file); trexio_close(file);
return qmckl_failwith( context, return qmckl_failwith( context,
@ -1141,7 +1169,7 @@ strncat(fname, "/chbrclf", 255);
printf("Test file: %s\n", fname); printf("Test file: %s\n", fname);
rc = qmckl_set_electron_walk_num(context, walk_num); rc = qmckl_set_electron_walk_num(context, walk_num);
rc = qmckl_trexio_read(context, fname); rc = qmckl_trexio_read(context, fname, 255);
if (rc != QMCKL_SUCCESS) { if (rc != QMCKL_SUCCESS) {
printf("%s\n", qmckl_string_of_error(rc)); printf("%s\n", qmckl_string_of_error(rc));

9
tools/build_doc.sh
View File

@ -9,12 +9,7 @@ if [[ -z ${srcdir} ]] ; then
exit 1 exit 1
fi fi
if [[ -z ${top_builddir} ]] ; then readonly DOCS=share/doc/qmckl/
echo "Error: srcdir environment variable is not defined"
exit 1
fi
readonly DOCS=${top_builddir}/share/doc/qmckl/
readonly ORG=${srcdir}/org/ readonly ORG=${srcdir}/org/
readonly HTMLIZE=${DOCS}/html/htmlize.el readonly HTMLIZE=${DOCS}/html/htmlize.el
readonly CONFIG_DOC=${srcdir}/tools/config_doc.el readonly CONFIG_DOC=${srcdir}/tools/config_doc.el
@ -47,7 +42,7 @@ function extract_doc()
for i in $@ for i in $@
do do
exported=${i%.org}.exported exported=${i%.org}.exported
exported=${top_builddir}/src/$(basename $exported) exported=${srcdir}/src/$(basename $exported)
NOW=$(date +"%m%d%H%M.%S") NOW=$(date +"%m%d%H%M.%S")
extract_doc ${i} &> $exported extract_doc ${i} &> $exported
rc=$? rc=$?

6
tools/build_makefile.py
View File

@ -49,8 +49,8 @@ def main():
f_test_o = "tests/test_"+i+"_f.$(OBJEXT)" f_test_o = "tests/test_"+i+"_f.$(OBJEXT)"
c_test = "tests/test_"+i+".c" c_test = "tests/test_"+i+".c"
f_test = "tests/test_"+i+"_f.F90" f_test = "tests/test_"+i+"_f.F90"
html = "share/doc/qmckl/html/"+i+".html" html = "$(htmldir_loc)/"+i+".html"
text = "share/doc/qmckl/text/"+i+".txt" text = "$(textdir_loc)/"+i+".txt"
i="src/"+i i="src/"+i
@ -267,7 +267,7 @@ def main():
for f in sorted(DEPS_DOC.keys()): for f in sorted(DEPS_DOC.keys()):
output += [ DEPS_DOC[f] + ": $(srcdir)/" + f + " $(htmlize_el)", output += [ DEPS_DOC[f] + ": $(srcdir)/" + f + " $(htmlize_el)",
"\t$(export_verbose)top_builddir=$(abs_top_builddir) srcdir=$(abs_srcdir) $(srcdir)/tools/missing bash $(srcdir)/tools/build_doc.sh $(srcdir)/"+f, "\tsrcdir=$(abs_srcdir) $(srcdir)/tools/missing bash $(srcdir)/tools/build_doc.sh $(srcdir)/"+f,
"" ] "" ]
output += ["endif"] output += ["endif"]